MATLAB多智能体编队控制：从一致性协议到仿真实现

简介：面向多智能体编队控制与一致性研究的初学者，这份Matlab程序包以一篇IEEE控制领域权威期刊论文为蓝本，由作者独立编程实现并完成运行验证。程序围绕时变编队控制问题展开，涵盖控制律设计、队形生成、仿真模型搭建和结果绘图等环节，能够直观展示多智能体在二维或三维空间内的编队演化过程。压缩包内共收录8个文件，包括4个m格式的源码主程序与辅助函数、1个slx格式的Simulink仿真模型、1个txt格式的详细使用说明、1篇论文原文以及1个资料快捷链接，整体体积约980KB，轻量易部署。配合文档中的操作指引和原文献，读者可逐段对照算法理解关键代码，快速上手复现实验，也可在此基础上修改参数探究不同编队策略。该资源已有6572人学习使用，适合作为课程设计、毕业设计或科研入门的仿真参考。

1. 从两台 AGV 的走廊相遇说起：MATLAB 里跑通的多智能体编队控制

先设想这样一个现场：三台搬运机器人同时从三个库位出发，目标是要排成“一字型”穿过一条 2 米宽的通道。如果每台机器人各自用路径规划冲到目的地，走廊里必然出现死锁。更常见的做法是给机群一个“队形状态量”，让每台机器人与相邻机器人保持固定间距，同时整支队伍沿给定轨迹移动。机器人之间的间距由控制协议去纠正，不需要中央调度的指令队列。这正是多智能体编队控制要解决的核心问题：只依赖局部通信，让所有智能体在期望队形上收敛，并保持队形整体运动。

标题里的“多智能体”“编队控制”“MATLAB”三个关键词，组合起来就是一套可以在实现前先用仿真验证的控制闭环。对研究机器人、无人机或集群算法的工程师来说，MATLAB 的优势在于线性代数、微分方程求解和绘图都是原生操作，编队控制的矩阵协议可以直接写成一页脚本。新手能照的最小例子理解一致性协议，熟练的人则能拿同一套框架快速验证不同拓扑下的收敛速度、切换队形和避让逻辑。下面顺着“模型 → 最小实现 → 参数整定 → 验证进阶”这条路线，把一套能复现、能改着用的多智能体编队控制 MATLAB 程序拆开讲清楚。

2. 编队控制的一致性数学模型：图拉普拉斯矩阵如何约束队形

多智能体编队控制的前提是“个体动力学 + 通信拓扑 + 一致性协议”三件事的准确建模。通信拓扑告诉每个智能体它能感知谁；一致性协议把拓扑信息变成位置修正量；个体动力学则决定这个修正量能否在真实约束下生效。MATLAB 程序中最先要写的，不是控制律，而是描述通信关系的邻接矩阵和拉普拉斯矩阵。

2.1 从单车模型到多智能体系统的状态表达

一维情况下，每个智能体有自己的位置 (x_i) 和速度 (v_i)。如果是小扰动场景，可以近似成双积分器模型：

[ \dot{x}_i = v_i, \quad \dot{v}_i = u_i ]

其中 (u_i) 是控制输入，也就是我们要设计的编队控制律。这里的“位置”不一定是物理坐标，也可以解释为队形参数，比如无人机编队的偏航角、无人物流车的纵向距离。MATLAB 编程里，所有智能体的位置和速度通常组合成一个大向量 state，方便用 ode45 或自写欧拉积分统一更新。

编队控制的目标不是让所有 (x_i) 相等，而是让它们之间的差收敛到一组期望值。因此定义队形偏置 (d_{ij} = d_i - d_j)，其中 (d_i) 是智能体 i 在目标队形中的坐标。以一字队形为例，两个相邻智能体如果期望间距是 1 米，那么 (d_1=0, d_2=1, d_3=2)。控制输入要作用于 (e_i = x_i - d_i)，而不是 (x_i)。这样一来，队形跟随就变成让所有 (e_i) 趋向相等，仍然是一个一致性收敛问题。

2.2 拉普拉斯矩阵与一致性协议的关系

假设通信拓扑是无向图，智能体 i 能收到邻居 (j \in N_i) 的相对位置信息。经典一致性协议为：

[ \dot{x}i = u_i = -\sum{j \in N_i} a_{ij}(e_i - e_j) ]

写成向量形式就是 (\dot{e} = -L e)，其中 L 是通信图的拉普拉斯矩阵：

[ L = D - A ]

A 是邻接矩阵，D 是度矩阵。这一关系是编队控制 MATLAB 程序里最常被直接矩阵化的部分。使用谱分析能得到非常直观的结论：如果通信图连通，拉普拉斯矩阵只有一个零特征值，其余特征值都大于零。因此 (e) 会渐近收敛到零特征值对应的特征向量空间，也就是所有分量相等的状态。

MATLAB 里验证上述性质通常用 eig 函数。构造一个有向环的邻接矩阵，再调用eig查看特征值分布，可以立刻判断出通信图是否连通。这个检查步骤在正式仿真之前就应该完成，否则后面跑出的发散曲线会让人误以为是控制参数问题。

2.3 为什么“通信拓扑 + 队形偏置”就能完成编队

有人会问：既然每对智能体都只控制相对误差，那整个队伍整体的位置由谁保证？答案是不需要全局位置。只要拓扑连通，所有 (e_i) 收敛到共同值，队形就固定了。队伍整体的平移则由队形中心的期望轨迹决定。

在 MATLAB 程序里，队形与轨迹常常解耦：控制输入里加上一项 (\dot{r}(t) - \sum a_{ij}(d_i - d_j))，其中 (r(t)) 是队形中心的期望轨迹。这时每个智能体控制律变成：

[ u_i = \dot{r} - \sum_{j \in N_i} a_{ij}[(x_i - d_i) - (x_j - d_j)] ]

这里 (\dot{r}) 作为前馈速度加入，保证队伍整体跟随轨迹。这个结构在后续第四章参数调优和第五章动态编队变换中都很关键。实际写代码时，可以先把 r(t) 设为常速直线，等这个框架跑通后再升级成曲线轨迹或多段路径。

3. 用 MATLAB 编写多智能体编队控制程序：最小可运行框架

现在开始实现一套最小可运行的多智能体编队控制 MATLAB 程序。这里先把智能体数设为 4，通信拓扑定为环形，每个智能体只与左右邻居通信。目标是让 4 个智能体在二维平面内形成一个正方形编队，并按固定速度向 x 轴正方向移动。

3.1 初始化通信拓扑与目标队形

程序的第一部分处理拓扑和队形参数。打开 MATLAB，新建脚本 formation_control_demo.m，按下面的代码开始。

clear; clc; close all; % 智能体数量 N = 4; % 环形通信拓扑：i 与 i+1、i-1 通信 A = zeros(N,N); for i = 1:N j = mod(i, N) + 1; % 右邻居 k = mod(i-2, N) + 1; % 左邻居 A(i,j) = 1; A(i,k) = 1; end % 度矩阵 D 和拉普拉斯矩阵 L D = diag(sum(A,2)); L = D - A; % 期望队形：正方形四个顶点，中心在原点 Dx = [1; -1; -1; 1]; Dy = [1; 1; -1; -1];

A 是邻接矩阵，A(i,j)=1 表示智能体 i 能接收 j 的状态信息。D 是用 sum(A,2) 计算每个节点的度数，再放到对角线上。拉普拉斯矩阵 L 可以直接用 D - A 得到，后面控制律计算全部基于 L，不会再显式循环邻居求和。Dx 和 Dy 分别是四个智能体在 x 轴和 y 轴上的队形偏置，它们共同组成边长为 2 的正方形。

3.2 主循环：欧拉积分与编队控制律

编队控制很多场景下并不需要 ode45 那种高精度求解器，固定步长的欧拉积分配合合适的 dt 已经足够。下面的代码完成初始化与动态演化。

% 仿真参数 dt = 0.02; % 步长 T = 6; % 总仿真时长 t = 0:dt:T; nStep = length(t); % 状态初始化：位置和速度 x = zeros(N,1); y = zeros(N,1); vx = zeros(N,1); vy = zeros(N,1); % 初始位置轻轻偏离目标队形，便于观察收敛过程 x0 = Dx + [0.3; -0.2; 0.1; -0.4]; y0 = Dy + [0.2; 0.4; -0.3; -0.1]; x = x0; y = y0; % 期望编队中心轨迹：匀速直线 r(t) = t r_dot = [1; 0]; % 中心速度 % 记录轨迹供绘图 x_hist = zeros(N,nStep); y_hist = zeros(N,nStep); x_hist(:,1) = x; y_hist(:,1) = y; % 主循环 for k = 1:nStep-1 % 误差状态 = 当前位置 - 目标队形位置 ex = x - Dx; ey = y - Dy; % 一致性控制律：u = r_dot - L * e ux = r_dot(1) - L * ex; uy = r_dot(2) - L * ey; % 双积分器模型，这里近似成 vx, vy 一阶模型 % 为了简化，速度直接设为控制输入 vx = ux; vy = uy; % 欧拉积分更新位置 x = x + vx * dt; y = y + vy * dt; % 保存历史数据 x_hist(:,k+1) = x; y_hist(:,k+1) = y; end

这里使用的一阶模型 (\dot{x}=u) 而不是双积分器模型，因为这样能让最小框架只保留控制协议本身，避免速度状态带来的额外参数。如果你的无人机模型需要速度约束，可以在后面加一个v = v + u * dt的积分环节，但注意稳定性边界会发生变化。代码里L * ex就是上一章一致性协议的矩阵形式，它一次性把当前误差向邻居均化的需求计算完毕。

3.3 绘图与结果分析

仿真结束后，用两幅图来观察收敛情况：一幅是 xy 平面的运动轨迹，另一幅是队形误差随时间下降的曲线。

% 轨迹图 figure; plot(x_hist', y_hist', 'LineWidth', 1.5); hold on; plot(x0, y0, 'ko', 'MarkerSize', 8, 'MarkerFaceColor', 'k'); plot(x_hist(:,end), y_hist(:,end), 'rs', 'MarkerSize', 10, 'LineWidth', 1.5); title('多智能体编队控制轨迹'); legend('智能体1','智能体2','智能体3','智能体4','初始位置','终止位置'); xlabel('x'); ylabel('y'); axis equal; grid on; % 队形误差图 error = sqrt((x_hist - repmat(Dx,1,nStep)).^2 + ... (y_hist - repmat(Dy,1,nStep)).^2); figure; plot(t, error', 'LineWidth', 1.2); title('各智能体与目标队形位置的误差'); xlabel('t/s'); ylabel('误差/m'); grid on; legend('智能体1','智能体2','智能体3','智能体4');

绘图代码里使用了 repmat 将队形偏置广播到每个时间步，从而按元素计算误差模长。第一条轨迹图如果看到一圈线从黑点开始最终压到红色方块上，说明编队形成。误差图中四条曲线应当单调下降并趋于零。如果出现振荡或误差不收敛，需要回到通信拓扑或时间步长这两处检查。

提示：这段程序里 diag(sum(A,2)) 假设无向图。如果使用有向图，度和入度的区别会直接改变拉普拉斯矩阵的定义，并影响收敛行为。在修改代码前先确认你的拓扑是有向还是无向。

4. 参数调优与 MATLAB 常见坑位：拉普拉斯特征值、步长与通信时延

能跑通最小框架之后，下一步是围绕程序中的关键参数做调优。很多人从网上拿到的“完美运行”程序，并不是因为用了多高深算法，而是把参数卡在了稳定范围内。这一章把最重要的几个参数和窗口期特征讲透。

4.1 拉普拉斯矩阵特征值与收敛速度的关系

从 (\dot{e} = -L e) 可以推出收敛速度由 L 的最小非零特征值 (\lambda_2) 决定。(\lambda_2) 越大，队形误差收敛越快。MATLAB 中直接用 eig(L) 就能看到特征值。对上面的环状拓扑，N=4 时 (\lambda_2=2)。如果把通信图改成全连接，(\lambda_2) 变成 4，收敛速度会变快，但通信代价也更大。

下表总结了不同通信拓扑对实现的影响，适合在实际编码前做选择：

如果你需要在 MATLAB 仿真中让编队快速成形，可以把拓扑从中等连通度开始，先观察误差曲线，再通过稀疏化逐步减少链路。判断链路是否可删除，就看生成图的第二小特征值是否仍大于零。这就是图论中连通性的数值判据。

4.2 欧拉积分的步长 dt 与稳定性上限

欧拉积分是显式方法，稳定性受步长限制。对一致性协议，误差状态更新为 (e_{k+1} = (I - dt \cdot L)e_k)。离散系统矩阵 (I - dt \cdot L) 的特征值必须落在单位圆内。因此 dt 需要满足：

[ dt < \frac{2}{\lambda_{\max}(L)} ]

对环状拓扑 N=4，拉普拉斯最大特征值为 4，所以 dt 必须小于 0.5。代码里的 dt=0.02，远低于上限，欧拉积分非常稳定。但如果你把智能体数增加到 100，链式/环形的最大特征值会变大，dt 要相应调小。这也能解释为什么同一套代码把 N 从 4 改成 20 后开始发散，往往不是控制律的问题，而是积分步长越界。

在 MATLAB 中可以用下面的代码验证当前参数是否满足稳定条件：

dt = 0.02; maxEig = max(eig(L)); disp(['拉普拉斯最大特征值: ', num2str(maxEig)]); disp(['理论步长上限: ', num2str(2/maxEig)]); if dt < 2/maxEig disp('步长满足稳定性约束'); else disp('步长过大，需要减小 dt'); end

4.3 通信时延的建模：从“理想程序”到“真实系统”

网络存在时延时的多智能体编队控制是另一个分开讨论的话题，但 MATLAB 最小框架里可以先用缓存队列模拟固定时延。以一步固定时延为例，每个智能体控制输入只能使用上一时刻收到的邻居状态。这相当于把控制律变成：

[ u_i(k) = -\sum a_{ij} [e_i(k) - e_j(k-1)] ]

在编队控制 MATLAB 程序里，可以简单地把ex(j)在计算前替换成ex_prev(j)。这样修改后很容易观察到误差曲线出现超调甚至振荡。通常减小增益或让每个智能体对历史状态做加权平均，能缓解时延影响。

下面是一段带时延的代码片段，可以直接在之前的框架上替换核心控制律部分：

% 用上一时刻误差构造时延影响 ex_prev = ex; ey_prev = ey; % 模拟一跳时延：只对邻居状态用历史值 ex_delay = ex; ey_delay = ey; for i = 1:N % 找出邻居索引 neighbors = find(A(i,:)); ex_delay(neighbors) = ex_prev(neighbors); ey_delay(neighbors) = ey_prev(neighbors); end ux = r_dot(1) - L * ex_delay; uy = r_dot(2) - L * ey_delay;

这里ex_delay中当前智能体自身的误差仍用当前值，邻居状态用上一帧值。这在通信链路局部时延模拟中是常见处理。对比有无时延的误差图，可以看到收敛时间明显拉长。此时为了避免发散，可以把dt减小到原来的 1/2，或者把控制律中乘一个比例增益 (k_p < 1)，等效于减弱控制作用。

4.4 多智能体博弈与编队控制中的角色互换

搜索多智能体相关内容时经常看到“多智能体博弈”这个热词，它和编队控制并不是同一件事，但两者可以结合。编队控制中的领航者-跟随者结构本身就有一种博弈关系：领航者决定轨迹，跟随者决定队形收敛。在 MATLAB 程序里，你可以把一个智能体设定为虚拟领航者，它不参与队形误差计算，只输出期望位置；其余跟随者通过跟踪虚拟中心来保持队形。这样在参数调优时，领航者的轨迹增益和跟随者的一致性增益就可以分开整定，避免所有智能体相互纠偏导致整体轨迹漂移。

5. 验证与进阶：动态队形变换、避免碰撞和 Simulink 联合仿真

最后一章把最小框架延伸到能放进论文或项目里的状态。重点讲两个方向：一是怎么验证“完美运行”不是恰好巧合，二是怎么把程序升级成更接近实际系统的动态编队控制。

5.1 用队形误差指标量化“完美运行”

与其“看起来收敛了”，不如定义一个可计算的编队误差指标 (E_{form})。最简单有效的定义是：

[ E_{form}(t) = \max_i \sqrt{(x_i(t) - Dx_i)^2 + (y_i(t) - Dy_i)^2} ]

当 (E_{form}(t)) 小于一个预设阈值（比如 0.05 米），就认为编队形成。在 MATLAB 里，这个指标可以放到每个时间步循环里，用来测量首次进入阈值的时间，也就是“收敛时间”。这个收敛时间可以作为后续调参的量化目标。

% 计算队形误差指标 formation_error = zeros(1,nStep); for k = 1:nStep ex = x_hist(:,k) - Dx; ey = y_hist(:,k) - Dy; formation_error(k) = max(sqrt(ex.^2 + ey.^2)); end % 找到首次小于0.05的时间 idx = find(formation_error < 0.05, 1); if ~isempty(idx) disp(['编队收敛时间为: ', num2str(t(idx)), ' 秒']); else disp('仿真结束前未达到期望编队误差阈值'); end

判定“完美运行”的标准不应只有一条轨迹图，而应包括：队形误差单调下降、稳态误差小于阈值、不同初始位置下均能重收敛。第三个条件需要把初始化放进一个 for 循环里跑多次，这是验证鲁棒性最便宜的自动化方式。

5.2 动态队形变换：把 (d_i) 从常量变成时变量

当编队需要在不同队形间切换（从一字队形切换到菱形队形），队形偏置不再是固定的 Dx、Dy，而是随时间变化的 (d_i(t))。常见做法是设计一个“队形插值器”，让队形偏置从一个矩阵平滑过渡到另一个矩阵。这样比直接切换更不容易激发大超调。

% 队形变换：一字形到 V 字形，插值时间 2 秒 switch_start = 2; switch_end = 4; D1x = [0; 1; 2; 3]; D1y = [0; 0; 0; 0]; D2x = [0; 1; 2; 3]; D2y = [0; 1; 0; -1]; for k = 1:nStep-1 tc = t(k); if tc < switch_start Dx = D1x; Dy = D1y; elseif tc > switch_end Dx = D2x; Dy = D2y; else alpha = (tc - switch_start) / (switch_end - switch_start); Dx = (1 - alpha) * D1x + alpha * D2x; Dy = (1 - alpha) * D1y + alpha * D2y; end % 后续控制律与主循环一致 ex = x - Dx; ey = y - Dy; ux = r_dot(1) - L * ex; uy = r_dot(2) - L * ey; x = x + ux * dt; y = y + uy * dt; end

关键参数是切换时间窗口switch_start到switch_end。窗口越长，队形变换越平滑，但占用路径长度更长。在 MATLAB 里观察队形变换是否过冲，可以画瞬时四边形的“形状面积”随时间变化。

5.3 在 Simulink 中复用 MATLAB 编队控制算法

如果项目后续要接电机模型或传感器模型，把纯 MATLAB 脚本转成 Simulink 模型会方便很多。最直接的办法是在 Simulink 里用 MATLAB Function 模块，把上面主循环里的控制律封装成函数：

function u = formation_controller(x, y, Dx, Dy, L, r_dot) ex = x - Dx; ey = y - Dy; ux = r_dot(1) - L * ex; uy = r_dot(2) - L * ey; u = [ux; uy]; end

在 Simulink 中，x 和 y 是由被控对象模块反馈回来的信号，Dx 和 Dy 是外部输入的队形指令。这样原本的 MATLAB 脚本就变成了控制器设计文档，而 Simulink 模型可以继续添加执行器饱和、噪声和时延。多智能体博弈的“技术选型”在这里表现为中心化与分布式控制结构的选择：Simulink 中如果所有智能体共用一个控制器模块，那就是中心化；如果每个智能体单独建一个子系统且只接收邻居信号，那就是分布式。后者的 MATLAB Function 输入引脚要去掉全局量，只用本地邻居信号。

5.4 常见验证技巧：分组测试拉普拉斯矩阵拓扑

最后给一个很容易被忽略的验证技巧：把拓扑和目标队形一起测试。在 MATLAB 中写一个循环，依次测试环状、链状、全连接三种拓扑，记录每种拓扑下的收敛时间。你会发现环状拓扑收敛速度较慢，但链路少，适合真实机器人；全连接收敛最快，但通信负担大。用语言表述就是：实际项目中先画通信图，再决定控制参数，顺序不能反。所有仿真脚本都可以通过主循环前的A矩阵换行来切换到不同拓扑，而不用改控制器代码。这个习惯能让你的编队控制 MATLAB 程序从“跑通一次”变成“改哪里都出结果”。

本文还有配套的精品资源，点击获取