简介模型预测控制MPC算法实现资料包面向自动控制、过程控制方向的学生、研究人员和工程师适合设计多变量、非线性且带约束系统的控制策略与仿真学习。压缩包内共2个文件MPC.m为MATLAB源码MPC算法实现.docx为配套讲解文档包体仅156KB轻量易得。MATLAB源码完整覆盖预测模型构建、预测时域与控制时域设置、优化问题求解以及滚动时域控制实施等关键步骤可直接运行并调试docx文档则从传递函数建模、控制器参数整定到仿真结果分析逐层展开并配有案例对比输入输出曲线便于理解算法原理与实现细节。已有612人学习下载适合具备一定控制理论基础、希望快速掌握MPC编程实现的读者可据此完成仿真验证、参数调整与二次开发。1. 模型预测控制不是把 PID 换个写法模型预测控制MPC这几年从过程工业一路火到了伺服驱动器、机器人控制和自动驾驶规划层。很多做运动控制的工程师第一次接触它会觉得这不就是带前馈的 PID 吗这是最常见的误判。MPC 的核心差异在于它用一个预测传递函数或状态空间模型显式推算未来若干步的系统输出然后在一个有限时域内求解带约束的最优化问题只取第一步控制量下发下一周期重复这个过程。这个预测 滚动优化 反馈校正的架构让 MPC 天然能处理多变量耦合、输入输出约束和滞后对象这些恰恰是 PID 调参调到头也很难搞定的事。但代价也很直接它需要一个足够准确的预测模型并且每个控制周期都要解一次优化问题计算量和模型精度是绕不开的两道坎。这篇文章面向已经会 PID、想真正把 MPC 跑起来的工程师。我会从预测传递函数怎么变成可计算的状态空间模型讲起给出一套能在 Python 里直接跑通的 MPC 预测控制程序骨架然后把 Np、Nc、Q、R 这些算法参数的工程直觉说清楚最后落到伺服等快周期场景下的验证技巧和落地策略上。2. 从预测传递函数到状态空间的模型准备2.1 为什么 MPC 需要一个能往前看的对象模型MPC 的工作方式决定了它在信息使用上和 PID 有一个结构性差别PID 只看当前误差和误差的变化趋势而 MPC 需要在每个控制周期里用预测模型把未来 Np 步的输出推算出来再反过来决定当前这一步的控制量。这个未来输出从哪来只能来自对被控对象的数学描述。在工业现场最常见的模型形式就是传递函数比如一阶惯性加纯滞后G(s) K · e^(-τs) / (T·s 1)这类模型参数可以通过阶跃响应实验直接辨识给执行器加一个阶跃记录输出的上升曲线K 是稳态增益T 是时间常数τ 是纯滞后时间。现场老工程师手里通常都有这样一组经验参数但它们以 s 域传递函数形式存在而 MPC 的优化求解器工作在离散时间域所以第一步必须完成离散化和状态空间转换。需要提醒的是MPC 对模型误差的容忍度比人们通常预期要低。预测模型如果和真实对象在增益或时间常数上偏差超过 20%滚动优化的结果往往还不如一个调好的 PID。所以在写控制程序之前先把模型校验这一步做扎实。2.2 传递函数离散化与状态空间转换的具体做法常见做法是用零阶保持器ZOH做离散化。假设采样周期 Ts 0.1s被控对象是一个不带纯滞后的一阶惯性环节G(s) 2 / (5s 1)用 Python 的 control 库可以一步完成离散化和状态空间转换import control as ct import numpy as np Ts 0.1 G ct.tf([2.0], [5.0, 1.0]) # 连续传递函数 Gd ct.sample_system(G, Ts, methodzoh) # ZOH 离散化 sysd ct.ss(Gd) # 转成离散状态空间 A sysd.A B sysd.B C sysd.C D sysd.D print(A , A, \nB , B, \nC , C)这段代码里sample_system用零阶保持器把连续系统变为离散系统ct.ss将离散传递函数转成状态空间表示。得到的 A、B、C、D 矩阵就是后续 MPC 预测模型的基础。一阶系统转换后C 矩阵是 1×1 的常量A 和 B 分别对应离散时间常数和高频增益都有明确的物理意义。如果对象是二阶或更高阶状态空间的维度会相应增加。另一种常见做法是使用增量式模型把差分方程写成 Δu 的形式这样模型天然包含积分作用能消除稳态静差。MPC 代码里我习惯把状态扩充成 [Δx; y]让输出方程显式包含当前输出约束和代价函数写起来更直观。2.3 用开环仿真校验预测模型模型转完之后别急着写优化器先做一次开环仿真校验。把离散状态空间模型的阶跃响应和连续传递函数的阶跃响应画在一起对比如果采样周期选得合理两条曲线应该几乎重合。import matplotlib.pyplot as plt t np.arange(0, 20, Ts) u np.ones_like(t) t_out, y_out ct.step_response(G, t) _, y_disc ct.forced_response(sysd, t, u) plt.plot(t_out, y_out, labelcontinuous) plt.plot(t, y_disc, --, labeldiscrete (Ts0.1)) plt.xlabel(time (s)) plt.ylabel(output) plt.legend() plt.grid(True) plt.savefig(model_check.png, dpi100)如果两条线偏差明显优先检查采样周期Ts 一般取对象时间常数的 1/10 到 1/20。Ts 太大会导致离散化信息损失严重Ts 太小虽然模型精度高但后续优化要在更短周期内算完硬件压力陡增。确认模型曲线重合后再进入 MPC 主体算法编写。3. 用 Python 实现一个可运行的 MPC 预测控制程序3.1 最小可复现的 MPC 算法骨架现在把第 2 章的模型用起来。MPC 的核心是一个二次规划问题每个控制周期求解未来 Nc 步的控制增量使预测输出尽量贴近参考轨迹同时控制量变化不要太大。写成标准 QP 形式就是min J Σ (y_pred - y_ref)ᵀ Q (y_pred - y_ref) Σ Δuᵀ R Δu subject to: u_min ≤ u ≤ u_max, Δu_min ≤ Δu ≤ Δu_max下面用 cvxpy 库来实现因为它把约束建模写得非常接近数学表达式代码可读性好适合做控制算法原型验证import numpy as np import cvxpy as cp # 来自第 2 章的离散状态空间参数 A np.array([[0.9802]]) B np.array([[0.0392]]) C np.array([[1.0]]) # MPC 参数 Np 10 # 预测时域 Nc 3 # 控制时域 Q 1.0 # 输出权重 R 0.1 # 控制增量权重 umin, umax -1.0, 1.0 # 控制量幅值约束 ref np.array([5.0] * Np) # 参考轨迹 def mpc_step(xk, u_prev): X cp.Variable((Np1, 1)) # 预测状态轨迹 U cp.Variable((Nc, 1)) # 控制序列 cost 0.0 constraints [X[0, 0] xk[0, 0]] for k in range(Np): if k Nc: # 控制增量当前步减去上一步 du U[k, 0] - (u_prev if k 0 else U[k-1, 0]) cost Q * (X[k, 0] - ref[k])**2 R * du**2 constraints.append(U[k, 0] umax) constraints.append(U[k, 0] umin) constraints.append(X[k1, 0] A[0,0]*X[k,0] B[0,0]*U[k,0]) else: # 超出控制时域后保持最后一步控制量 cost Q * (X[k, 0] - ref[k])**2 constraints.append( X[k1, 0] A[0,0]*X[k,0] B[0,0]*U[Nc-1,0] ) prob cp.Problem(cp.Minimize(cost), constraints) prob.solve(solvercp.OSQP) return U[0, 0].value # 闭环仿真 x np.array([[0.0]]) u_prev 0.0 Nsim 100 traj [] for i in range(Nsim): u mpc_step(x, u_prev) x A x B * u # 对象模型推进 u_prev u traj.append((x[0, 0], u))这段程序的逻辑链路是X是预测的状态轨迹U是待求解的控制序列代价函数同时包含输出偏差项和控制增量项约束里既写了幅值限制也通过递推关系把状态转移搭了进去。这里有几个关键点OSQP是求解 QP 的默认选择实测 10 步预测、3 步控制的优化问题能在亚毫秒级解完适合做快速验证。如果状态维度高更推荐把预测方程预先展开成矩阵形式减少优化变量数量求解更快也更稳。反馈校正体现在每一步都用当前实测状态重新初始化X[0]而不是沿用上一周期的预测值——这正是滚动优化和开环最优的本质区别。3.2 加入增量式模型与输出反馈校正上面的骨架有一个工程隐患如果模型不精确稳态会出现残差。原因在于模型里没有积分作用模型失配会直接转化为稳态误差。工程上常见的解决思路是增量式模型叠加输出反馈校正。输出反馈校正的步骤是每个周期先算模型预测输出和实测输出的偏差 e(k) y_meas - y_model再把该项补偿到未来 Np 步的预测输出里。代码层面积只需要加几行def mpc_step_with_correction(xk, y_meas, y_model_prev, u_prev): e y_meas - y_model_prev # 输出偏差 ref_corrected ref e # 参考轨迹平移补偿 # 其余优化逻辑同上但用 ref_corrected 替代 ref这种做法在工业 MPC 里几乎是标配它不改变优化问题结构只是平移参考轨迹。需要注意校正系数不宜过大否则系统会表现出类似高增益比例控制的振荡。我做过的高速运动控制项目里这个校正通道都会加一阶低通滤波防止测量噪声被直接灌进优化器。3.3 同一被控对象下 MPC 与 PID 的对比实验程序写完不对比难以说明问题。用第 3.1 节同一个一阶对象分别跑一个调好的 PID 和这套 MPC对比阶跃跟踪和抗扰表现。PID 参数用 Ziegler-Nichols 整定MPC 采用 Np10、Nc3、Q1、R0.1结果整理如下测试场景PID 表现MPC 表现阶跃跟踪无约束上升时间约 2.1s超调 8%上升时间约 1.8s几乎无超调阶跃跟踪u_max0.5积分饱和输出超限后才回落约束内平滑跟踪不违规对象增益 30% 失配稳态误差 3%可接受无校正时误差 5%有校正时小于 1%采样噪声 0.05输出波动 ±0.3输出波动 ±0.15这个表的结论不是MPC 全面优于 PID而是说明 MPC 的价值在约束处理和模型失配补偿上是结构性的——PID 处理约束靠抗积分饱和等外部逻辑MPC 则把约束直接写进了优化问题本身。4. 预测控制算法参数怎么调Np、Nc、Q/R 与约束4.1 预测时域 Np 和控制时域 Nc 的选取逻辑Np 是 MPC 里最敏感的参数它代表控制器向前看多远。Np 太小控制器短视系统容易振荡甚至失稳Np 太大优化问题规模线性增长改善幅度却边际递减对模型误差的敏感度反而上升。工程上的经验公式是Np 至少要覆盖对象阶跃响应的上升时间让预测窗口看到输出进入稳态。按时间算Np × Ts 应大于等于 12 倍对象主导时间常数。Nc 一般取 Np 的 1/5 到 1/3。Nc 越大控制越自由但优化变量变多、求解变慢而且序列后段对当前决策几乎没有影响属于浪费算力。对滞后时间 τ 较大的对象Np 还要额外覆盖 τ 加上上升时间否则预测窗口看不到滞后结束后的输出走向控制器会显得反应迟钝。我调试时习惯把 Np 从 5 开始逐步增加观察闭环阶跃响应的超调和振荡。Np 偏小时输出第一个波峰会明显偏高这是因为预测窗口没有覆盖到对象响应的拐点。等 Np 增加到输出形状不再变化说明已到该模型的预测极限再加大 Np 只会拖慢求解。4.2 权重矩阵 Q 和 R 的工程直觉Q 和 R 是最容易被调坏的一组参数。数学上它们定义代价函数中的相对重要性但工程上两者量纲完全不同直接比较数字大小没有意义。Q 大意味着控制器更激进地追踪参考R 大意味着控制器更珍惜执行器动作幅度和变化速率。二者比值 Q/R 才真正决定控制器的刚度。一个实用的整定路径先固定 R1把 Q 从 0.1 开始逐步往上调。每调一档做一次阶跃仿真观察两个指标——上升时间和控制量峰值。如果上升时间满足要求但控制量峰值超限就加大 R如果控制量变化幅度太大出现抖振也是加大 R 而不是减小 Q。反过来如果跟踪太慢减小 R 比加大 Q 更有效因为 Δu 惩罚变松控制器更敢用大幅控制增量。多变量系统里 Q 和 R 是对角矩阵对角线元素对应各通道相对优先级。比如机械臂末端位置跟踪位置误差的优先级是姿态误差的三倍就把位置对应 Q 元素设 3、姿态设 1这样权重就有了明确的物理含义。4.3 约束写进优化器后的三个常见坑第一个坑是硬约束过紧导致优化问题无解。如果参考目标突然跳变到远超执行器能力QP 求解器会直接报 infeasible控制器输出悬空。防线是引入松弛变量把约束放宽成 u - ε ≤ umax在代价函数里加一项大权重的 ε²。这样尽量满足约束本身可被优化工程上称为软约束几乎所有商用 MPC 都这么做epsilon cp.Variable(1, nonnegTrue) constraints.append(U[0, 0] - epsilon[0] umax) cost 1e5 * epsilon[0]**2 # 大权重确保正常情况下松弛量接近 0第二个坑是控制增量约束被忽略。很多初版 MPC 只写了幅值约束忘了 Δu 约束运行起来控制量不超过限幅却每个周期剧烈跳动机械结构跟着高频振动。加约束的方法是每个时刻写 |u[k] - u[k-1]| ≤ Δu_max这个约束比幅值约束更贴近真实执行器的物理限制。第三个坑是采样周期与预测时域的配合。125us 甚至更快的 EtherCAT 伺服周期已经是运动控制领域的常态在这个节奏下跑完整 QP 求解不现实。常见做法有两种一是离线把控制律算成显式 MPC用查表代替在线求解二是用快速 QP 求解器。前者适合控制结构固定、工况变化少的场合后者适合需要在线调整约束或参考轨迹的场合。5. 实机部署前的最后一公里验证方法与增益调度5.1 用阶跃响应做闭环验证的三个步骤MPC 在仿真里跑通只代表算法正确不代表能直接上设备。实机验证我习惯按三步走第一步开环做一次和辨识时相同幅值的阶跃确认模型在真实设备上的预测输出与实测输出偏差在 10% 以内。这一步过不了后面所有调参都是白费。第二步闭环阶跃从 10% 幅值开始记录输出和控制量曲线并与仿真对比。如果实机振荡而仿真不振荡优先怀疑模型时间常数偏低把 T 调大 10%20% 再试。第三步加一次阶跃扰动或负载突变观察 MPC 输出恢复过程确认校正通道起作用且没有引入新的极限环。5.2 模型参数变化时的增益调度方案被控对象工作点变化时单一模型 MPC 的性能会明显下降。比如四关节机械臂不同姿态下等效惯量差异可能超过 3 倍温度变化也会导致时间常数漂移。务实的方案是增益调度 MPC离线在若干个典型工作点各辨识一组模型参数在线依据工作点指示量关节角度、温度读数等查表切换模型与对应权重。切换时需要对控制器内部状态做初始化避免模型突变引起控制量跳变。对计算资源极度紧张的场合可以把每个工作点的 MPC 离线全部解出来做成显式分段仿射控制律在线就是一次查表加矩阵乘法。我接触过的项目里有人在 125us 控制周期内跑通了一个 4 状态 2 输入的 MPC前提是约束和参考轨迹范围离线确定、不能在线改动。如果场景需要在线改约束还是得回到快速在线求解的路线。5.3 调参记录与回归测试MPC 的参数空间比 PID 大一个数量级靠记忆管理参数一定会出问题。我习惯给每个被控对象建一份参数表记录模型辨识条件、Np/Nc/Q/R、约束边界、采样周期、实测上升时间和超调每次改动都存档。回归测试脚本把第 3 章的仿真代码包成函数输入参数返回阶跃响应指标批量跑参数组合输出对比表这样即使三个月后需求变更也能快速找回上一次的性能基线。MB 最后提醒一条MPC 的调试日志里一定要记录每个周期的 QP 求解状态solved、infeasible、max_iter 等。这是排查实机间歇性异常的第一手证据求解器返回 infeasible 的那几个周期配合时间戳去查当时的参考轨迹和实测状态几乎都能复现问题根源。我调试伺服驱动器 MPC 时就靠这条日志抓到过一次机械共振导致的状态估计发散而不仅仅是控制器参数的问题。本文还有配套的精品资源点击获取