1. 从MPC到二次规划为什么我们绕不开quadprog做模型预测控制MPC的朋友尤其是刚上手实现无人车、机械臂这类系统时大概率会在某个深夜对着一个叫“二次规划”的优化问题挠头。你辛辛苦苦推导了状态空间方程设计了预测时域和控制时域把滚动优化问题转化成了一个标准形式最后发现核心的求解器你还没选。这时候quadprog这个名字就会频繁地出现在各种教程、论文和开源代码里。它不是什么高深莫测的黑科技而是MATLAB和Octave等科学计算环境里一个求解二次规划问题的内置函数。但问题来了MPC里产生的二次规划问题其目标函数的Hessian矩阵就是那个二次型系数矩阵性质可不一样它可能是正定的、半正定的甚至是负定的。不同的性质直接关系到quadprog能不能解、怎么解以及解出来靠不靠谱。今天我们就抛开那些复杂的理论推导直接切入实战聊聊怎么用quadprog这个“老伙计”去搞定MPC中可能遇到的各种性质的二次规划问题特别是那些容易让人栽跟头的半正定和负定情况。2. 二次规划的标准形式与quadprog的输入要求在把MPC问题塞进quadprog之前我们必须先统一语言搞清楚双方各自的标准是什么。一个标准的二次规划问题通常写成这样最小化f(x) 1/2 * x^T * H * x f^T * x满足约束A * x bAeq * x beqlb x ub这里的x就是我们的决策变量在MPC里它通常代表未来一段时域内的控制输入序列。H是一个对称矩阵它就是问题的核心——Hessian矩阵。f是线性项的系数向量。不等式约束A*x b可以表示控制量的幅值限制如方向盘转角、油门开度有上下限等式约束Aeq*x beq可能用于构造动力学方程或其他硬性约束。边界约束lb和ub则是另一种形式的变量上下限。现在看MATLAB的quadprog函数它的基本调用格式是[x, fval, exitflag, output, lambda] quadprog(H, f, A, b, Aeq, beq, lb, ub, x0, options)关键点在于quadprog的官方文档里会强调矩阵H在对应等式约束的子空间上必须是正定的positive definite。这句话是很多问题的根源。为什么因为quadprog内部主要采用有效集法Active-Set Method或内点法Interior-Point Method这些算法在理论上要求问题的Hessian矩阵正定以保证迭代过程中子问题是严格凸的从而有唯一的最小值算法才能稳定收敛。那么MPC问题中的H矩阵是什么它通常来源于对系统性能指标如跟踪误差和控制量加权的二次型求和。一个设计良好的MPC其目标函数通常要求是凸的以保证优化问题的可解性。因此H矩阵理论上应该是**半正定positive semi-definite或正定positive definite**的。负定negative definite的情况则意味着目标函数是凹的求最小化问题会趋向于无穷小这在物理系统中通常对应不稳定的控制律是设计错误。所以我们实战中面对的主要是两类情况理想的正定矩阵和更常见的、带挑战性的半正定矩阵。负定矩阵更多是作为一个需要被识别和修正的错误标志出现。3. 实战场景一H为正定矩阵——最理想的情况当你的MPC问题中控制输入权重矩阵R是严格正定的并且状态权重矩阵Q也是半正定或正定的且系统能观那么最终构造出的整体H矩阵很大概率是正定的。这是quadprog最喜欢的情况。假设我们为一个简单的无人车横向控制设计MPC控制量是前轮转角delta。我们可能这样构造问题% 假设经过推导Hessian矩阵H是正定的 H [2.5, 0.1; 0.1, 1.8]; % 一个小例子实际维度很高 f [-1; 0.5]; % 线性项系数 % 控制量约束转角限制在±30度以内 A [1, 0; -1, 0; 0, 1; 0, -1]; b [0.5236; 0.5236; 0.1; 0.1]; % ±30度±0.1弧度 % 无等式约束 Aeq []; beq []; % 调用quadprog options optimoptions(quadprog, Display, iter, Algorithm, interior-point-convex); [x_opt, fval, exitflag] quadprog(H, f, A, b, Aeq, beq, [], [], [], options);当H正定时quadprog的求解会非常顺畅。exitflag会返回1表示成功收敛。Algorithm选项可以选择‘interior-point-convex’默认或‘trust-region-reflective’需提供梯度等。在正定情况下算法选择相对自由收敛速度快数值稳定性好。注意即使H是正定的也要注意数值精度问题。如果H的条件数condition number非常大即最大特征值和最小特征值比值极大问题就是“病态”的。这会导致求解器数值不稳定可能误报失败或给出精度很差的解。在MPC中这可能源于Q和R的权重比例设置极端不合理。一个实用的技巧是在调用quadprog前检查cond(H)如果过大比如 1e10就需要重新审视你的权重设计了。4. 实战场景二H为半正定矩阵——最常见也最棘手的挑战半正定矩阵意味着H的特征值中至少有一个为零其余非负。在MPC中这太常见了。比如你的状态权重矩阵Q是半正定的可能只对部分状态变量进行惩罚或者控制权重矩阵R在某些控制通道上设为零权重允许该通道控制量自由变化。此时目标函数在某个或某些方向上是“平坦”的存在无穷多组解都能达到相同的最小值。quadprog对于严格的半正定问题其默认算法可能会报错提示 “Hessian is not positive definite”。但这不代表无解而是需要一些技巧来处理。4.1 方法一正则化Regularization——最实用的工程技巧正则化的核心思想是给H矩阵加上一个很小的正定矩阵通常是单位矩阵的倍数从而将其“微调”为正定矩阵。% 原始的、可能半正定的H H_original ...; % 从MPC问题构造而来 % 正则化参数一个很小的正数 epsilon 1e-6; [n, ~] size(H_original); H_regularized H_original epsilon * eye(n); % 使用正则化后的H调用quadprog [x_opt, fval] quadprog(H_regularized, f, A, b, Aeq, beq, lb, ub);为什么这样做是有效的数学上epsilon * eye(n)是一个正定矩阵。任何半正定矩阵加上一个正定矩阵结果一定是正定的。这确保了quadprog的输入满足其严格正定的要求。物理上这相当于在目标函数中增加了一项(epsilon/2) * ||x||^2即对所有控制输入施加了一个极其微小的二次惩罚。它倾向于在所有可能的最优解中选择一个范数最小的解即控制能量最小的解。这在工程上往往是可接受的甚至是有益的因为它能避免控制量无意义的剧烈震荡。参数epsilon的选择至关重要不能太大比如1e-3或更大这会显著改变原问题的解扭曲MPC的性能。你的控制器可能变得过于“懒惰”。不能太小比如1e-12或更小在数值计算中可能无法有效改善H的条件数求解器可能依然报错或数值不稳定。经验范围1e-9到1e-6是一个常见的、安全的起始尝试区间。你需要根据你问题中H矩阵元素的典型大小范数来调整。一个原则是正则化项epsilon * I的Frobenius范数应远小于原始H矩阵的Frobenius范数例如小3到6个数量级。4.2 方法二使用更鲁棒的算法选项MATLAB的quadprog提供了不同的算法。当检测到H非正定时可以尝试指定使用‘active-set’算法。options optimoptions(quadprog, Algorithm, active-set); [x_opt, fval, exitflag] quadprog(H, f, A, b, Aeq, beq, lb, ub, [], options);有效集法Active-Set对于某些半正定问题可能比内点法Interior-Point更鲁棒因为它通过迭代探索约束边界来寻找最优解对Hessian正定性的要求有时可以放宽。但请注意官方文档仍建议Hessian在等式约束零空间正定且有效集法在处理大规模问题时可能较慢。4.3 方法三从问题源头审视——MPC设计修正如果半正定性是由于MPC设计本身导致的比如故意将某个控制输入的权重R(i,i)设为0那么你需要问自己这是否合理如果合理例如某个电机控制通道确实不需要能耗惩罚那么使用正则化方法是合适的。如果不合理可能是建模疏忽那么应该修正你的权重矩阵Q和R确保它们至少是半正定的并且组合后的H在物理意义上应该是正定的。通常给R矩阵的对角线元素赋予一个非常小但非零的正值如1e-4是更规范的做法它既保留了“几乎无惩罚”的意图又从根源上保证了H的正定性。5. 实战场景三H为负定矩阵——识别错误与问题诊断如果你在调用quadprog时遇到了与半正定类似的错误但通过检查发现H的特征值有负数甚至全部为负那么恭喜你你很可能发现了MPC设计中的一个bug。负定的H意味着你的目标函数是凹的求最小值问题在无约束情况下会趋向于负无穷这在实际的物理系统控制中是没有意义的。如何诊断计算特征值eig(H)。如果存在明显的负特征值如 -0.5, -2.3那就是负定或不定。检查权重矩阵回顾你构建H矩阵的代码。H通常由Q状态权重和R控制权重在预测时域内分块对角堆叠而成。确保Q是半正定或正定的。对于状态误差惩罚这通常是成立的。R是正定的。这是最常见的错误来源。R矩阵必须正定即所有对角线元素控制权重必须为正数。如果某个R(i,i) 0就会导致H出现负定或半正定的块。检查矩阵组装过程在复杂MPC中H可能由更复杂的公式生成涉及系统矩阵A,B等。确保矩阵乘法、转置、堆叠的代码正确无误。一个符号错误就可能导致正负号颠倒。如何处理一旦确认是负定不要试图用正则化强行求解。一个巨大的epsilon或许能让H变成正定但得到的“解”对你的控制系统毫无价值甚至危险。正确的做法是修正你的R矩阵将所有控制权重设为正数。重新推导H矩阵的公式检查二次型性能指标的推导过程。验证系统模型确保状态空间模型(A, B)正确特别是在线性化工作点附近。6. quadprog求解失败排查清单与高级选项调优即使H矩阵性质没问题quadprog也可能因为其他原因失败。这里提供一个排查清单检查约束可行性你的约束A*x b,Aeq*x beq,lb x ub是否本身互相矛盾是否存在一个x能同时满足所有约束可以用linprog或简单测试一个初始点x0来验证可行性。提供合理的初始点x0虽然quadprog可以不提供初始点但一个可行的、靠近解空间的初始点x0能显著帮助算法尤其是有效集法和内点法加速收敛并避免陷入局部区域。在MPC中一个很好的初始点猜测是上一时刻求解得到的最优控制序列去掉第一个值再补上一个零或上一时刻的最后一个值热启动。调整优化选项options optimoptions(quadprog, ... OptimalityTolerance, 1e-8, ... % 优化容忍度默认1e-8 ConstraintTolerance, 1e-8, ... % 约束容忍度默认1e-8 StepTolerance, 1e-12, ... % 步长容忍度对高精度问题可调小 MaxIterations, 200, ... % 最大迭代次数复杂问题需增加 Display, final); % 输出信息级别off, final, iter如果求解器在最大迭代次数内未收敛尝试增加‘MaxIterations’。如果解不满足约束可以适当收紧‘ConstraintTolerance’但注意不要小于数值误差水平。对于病态问题可以尝试使用‘interior-point-convex’算法它通常比‘active-set’对病态问题更鲁棒。检查输出信息exitflag和output结构体包含了丰富的诊断信息。exitflag 1成功。exitflag 0超过最大迭代次数。exitflag -2问题不可行。exitflag -3问题无界对于凸QP这通常意味着H半正定且无约束时目标函数可能无下界但结合约束后应可解。output.message会给出更详细的文本描述。7. 在MPC滚动优化中集成quadprog的工程实践在实际的MPC循环中我们不仅仅是调用一次quadprog而是在每个采样周期都调用。这带来了额外的工程考量计算效率MPC要求实时或准实时求解。quadprog对于中小规模问题决策变量维度几百以内速度很快。但对于大规模问题可能需要考虑更专用的QP求解器如OSQP, qpOASES或利用问题结构如稀疏性的求解方法。热启动Warm Start这是提升MPC在线计算速度最关键的技术。将上一个周期求解的最优序列x_opt作为当前周期quadprog的初始点x0。由于相邻时刻的优化问题高度相似热启动能极大减少迭代次数。% 第一个周期没有历史解用零初始化或简单猜测 x0 zeros(n_variables, 1); for k 1 : N_simulation_steps % 更新当前时刻的状态、参考轨迹等构造新的H, f, A, b等 [H_k, f_k, A_k, b_k, Aeq_k, beq_k] build_mpc_problem(current_state, reference); % 使用上一周期的解作为热启动需要做适当的移位 % 假设x_opt_prev是上一周期求得的整个时域的控制序列 if k 1 x0_k x0; else x0_k warm_start_shift(x_opt_prev); % 自定义函数移位并补零 end % 调用quadprog [x_opt_k, ~, exitflag] quadprog(H_k, f_k, A_k, b_k, Aeq_k, beq_k, lb, ub, x0_k, options); if exitflag 0 warning(QP failed at step %d. Exitflag: %d, k, exitflag); % 实施故障安全策略例如使用上一控制量或备用控制器 u_k backup_control; else % 提取当前时刻的控制量解序列的第一个元素 u_k extract_first_control(x_opt_k); % 存储本次解用于下一次热启动 x_opt_prev x_opt_k; end % 将u_k施加给系统并模拟/测量得到下一个状态 current_state simulate_system(current_state, u_k); end求解失败处理在循环中必须对quadprog的失败exitflag 0有预案。简单的策略包括沿用上一时刻的控制量、切换到PD控制器、或者放松约束后重新求解。一个健壮的工业MPC代码必须有这样的容错逻辑。最后我想分享一个在调试MPC与quadprog集成时非常有用的小技巧在开发阶段将每个采样周期构造出的H,f,A,b等矩阵和向量保存下来例如保存到一个.mat文件或结构体数组中。当求解器报错时你可以轻松地复现那个出错时刻的完整QP问题在MATLAB命令行里单独加载和分析它计算H的特征值、检查约束的可行性从而精准定位问题。这比在线调试整个仿真循环要高效得多。记住quadprog只是一个工具而你对MPC问题本身的理解——它的数学性质、物理意义和工程约束——才是成功实现稳定、高效预测控制的关键。