二次规划与积极集法：原理、实现与工程实战

写代码的人应该都有过这种经历：目标函数明明是个漂亮的凸二次函数，求个导、令它等于零，几秒钟就能写完解析解的代码，结果一旦加上约束，解出来的点直接跑到不可行域外面去了。我在做小车MPC轨迹跟踪的时候就被这个问题卡过大半天——二次规划、约束条件、状态限幅都在，算出来的控制量却让执行器饱和，整条轨迹全是毛刺。后来把思路放到积极集法（active set method）上，一步步看它怎么拆解约束、怎么切换工作集，才算真正把二次规划的问题看透。这篇文章把二次规划、积极集法的核心原理、手算推演、可运行代码以及工程落地中常见的问题完整讲一遍，适合刚碰带约束二次规划的学生，也适合已经会调现成求解器、但想搞清楚黑盒内部到底在迭代什么的工程师。

1. 先弄明白二次规划到底难在哪

1.1 标准形式：二次目标加线性约束

二次规划的标准形式可以写成这样：

[ \min_x ; \frac{1}{2} x^T H x + g^T x ] [ s.t. ; A x \le b ]

如果问题里还有等式约束，就再加一项 (A_{eq} x = b_{eq})。这里 (H) 是Hessian矩阵，凸二次规划要求 (H) 至少半正定，也就是说目标函数是个开口向上的碗，碗底是唯一的（严格正定时）。(g) 是线性项系数，约束全部是线性的。

这种形式在工程里出现频率高得惊人。最小二乘回归加约束就是一个典型例子：比如 (y = X\beta) 的回归问题，加上 (\beta \ge 0) 或者 (\sum \beta_i = 1) 之后，目标函数展开就是 (0.5 \beta^T (X^T X) \beta - (X^T y)^T \beta + const)，其中 (H = X^T X) 天然对称半正定，(g = -X^T y)，约束是线性不等式。做投资组合优化时，目标是让组合方差最小，(H) 是协方差矩阵，约束是权重和为1、权重非负，这也是一个典型的凸二次规划。

另一个我接触最多的是模型预测控制（MPC）：每个控制周期都要解一个二次规划，目标函数里含状态偏差和控制增量的惩罚，约束是执行机构饱和幅度、状态量上下限，本质上就是在每个采样时刻求一组最优控制增量。非线性规划领域里的SQP方法，也是把非线性问题在每次迭代点处展开，解一个二次规划子问题来得到搜索方向。可以说，只要优化问题里出现"二次目标 + 线性约束"，二次规划就出现了。

1.2 求导失灵：不等式约束带来的组合难题

无约束二次规划的解法太简单了：目标函数对 (x) 求导并令其等于零，得到 (Hx + g = 0)，一步得到解析解 (x^* = -H^{-1}g)。如果这个解恰好落在可行域内，那问题确实结束了。但绝大多数时候它不会那么听话。

我早年犯过的错误就是拿到带约束的问题后，先算无约束解，再手动把越界的变量"掰"回边界上。简单问题这么糊弄还能勉强对上结果，稍微复杂一点就崩。原因在于不等式约束的最优解并不满足 (\nabla f = 0)，它满足的是KKT条件：

[ Hx + g + A^T \lambda = 0 ] [ \lambda_i \ge 0, \quad \lambda_i (b_i - a_i^T x) = 0 ]

其中 (\lambda) 是拉格朗日乘子，(a_i^T) 是第 (i) 条约束的系数向量。互补条件 (\lambda_i (b_i - a_i^T x) = 0) 的含义是：如果第 (i) 条约束在最优解处没有触到边界（即 (b_i - a_i^T x > 0)），那么它对应的乘子必须为零；反过来，如果约束生效了，目标函数极值点上就一定会有该约束对应的推力，这个推力就是乘子。

难就难在"哪些约束会生效"这件事事先不知道。(m) 条约束里任意组合都可能成为有效约束，这是一个 (2^m) 级的组合选择问题。无约束求导解决不了它，因为它本质上不是连续优化问题里"沿梯度走到极值"的问题，而是要在离散的约束激活模式里找出正确的一套。用个生活类比：出门旅行收拾行李，你不可能把所有季节的衣服都背上，只能先按最可能的天气打包，路上天气变了再往包里加衣服或减衣服。积极集法干的就是这个事。

1.3 积极集法的破局思路：主动"猜"出有效约束

积极集法的核心思路很朴素：既然最优解处只有一部分约束生效，那我就维护一个"我认为生效"的约束集合，叫工作集（working set）。在这个工作集下，把不等式约束当成等式约束来处理，解一个等式约束的二次规划子问题。解完之后做两个检查：检查新点是否满足所有约束，检查工作集约束对应的乘子是否非负。如果出了问题，就调整工作集——要么加一条刚触界的约束，要么踢掉一条乘子为负的约束——然后重新解子问题。

每次迭代只调整工作集里的一个或少数几个约束，逐步逼近真正的最优激活模式。这避免了直接面对 (2^m) 的组合爆炸问题，理论上只要不退化，有限步就能收敛。而且它有一个非常好的工程性质：如果前后两个QP问题的约束激活模式差得不多，上一次解完的工作集可以直接拿来当这次的初始工作集，省掉大量重复迭代。这个性质在MPC里几乎是决定性的优点，后面我会专门展开。

2. 积极集法的两个核心判据：触界与乘子

2.1 工作集(active set)到底是什么

在任意一个可行点 (x) 上，把所有满足 (a_i^T x = b_i) 的约束收集起来，叫这个点的有效约束集。有效约束就是在当前点"触到边界"的约束。积极集法在迭代过程中维护的"工作集"，本质上就是当前迭代点处认为应该生效的约束集合。

为什么这个集合这么重要？因为在最优点处，真正起作用的只有这些约束。假设你已经知道了最终的最优有效约束集合 (W^*)，那么原问题就退化成一个等式约束问题：

[ \min_x ; \frac{1}{2} x^T H x + g^T x, \quad s.t. ; a_i^T x = b_i, ; i \in W^* ]

等式约束问题就好解多了，直接拉格朗日乘子法转线性方程组。所以积极集法整个过程可以概括成：不断猜测 (W)，求解对应的等式约束子问题，再用KKT条件验证猜测是否正确，不正确就修正。

2.2 等式约束子问题是怎么被解出来的

假设当前迭代点是 (x_k)，工作集是 (W)。我不直接解以 (x) 为变量的等式约束问题，而是换成求解一个搜索方向 (p)：

[ \min_p ; \frac{1}{2} p^T H p + (Hx_k + g)^T p ] [ s.t. ; a_i^T p = 0, ; i \in W ]

这个转换的意义在于：工作集内的约束在当前点已经被"钉"住了，所以沿搜索方向走时，这些约束必须保持等式成立，即 (a_i^T (x_k + \alpha p) = b_i)，等价于 (a_i^T p = 0)。目标函数里 (p) 的线性项系数是当前点处的梯度 (Hx_k + g)。

对这个问题写拉格朗日函数，对 (p) 和乘子 (\lambda) 求导置零，会得到一个线性方程组：

[ \begin{bmatrix} H & A_W^T \ A_W & 0 \end{bmatrix} \begin{bmatrix} p \ \lambda \end{bmatrix}

\begin{bmatrix} -(Hx_k + g) \ 0 \end{bmatrix} ]

其中 (A_W) 是工作集约束系数拼成的矩阵，每一行是一条约束。只要 (H) 正定、工作集约束线性无关，这个KKT矩阵就是非奇异的，直接用高斯消去法或Cholesky分解就能解。这也是为什么"约束线性无关"在积极集法里是个前提条件；实际工程中遇到约束线性相关时必须做处理，否则矩阵奇异，一步就报错。每次迭代的核心计算量就集中在这个线性方程组的求解上。

2.3 什么时候加约束，什么时候踢约束

解出 (p) 之后分两种情况。

第一种情况是 (p) 足够接近零向量。这说明在当前工作集下，目标函数已经无法再改进了，当前点是一个"局部最优候选点"。此时需要检查工作集内约束的拉格朗日乘子 (\lambda)。如果所有乘子都大于等于零（数值上允许一定容差），KKT条件全部满足，当前点就是全局最优解。但如果出现了某个 (\lambda_j < 0)，说明第 (j) 条约束在当前点"拉"住了目标函数，就像一只拖后腿的手，把目标困在了一个并不是真正最小的位置。释放这条约束，目标函数还能继续下降。所以要把乘子最负的那条约束从工作集里踢掉，然后重新求解子问题。

第二种情况是 (p) 不为零。(p) 是当前工作集下让目标下降最快的可行方向，沿着它走一步，目标一定下降。能走多远？上限是1，因为子问题沿这个方向的精确最优步长就是1，再远目标就开始反弹。但在到达步长1之前，很可能先碰到某条不属于工作集的约束的边界，不能再往前走了。对所有不在工作集的约束算一下最大可行步长：

[ \alpha_{\max} = \min_{i \notin W, \ a_i^T p > 0} \frac{b_i - a_i^T x_k}{a_i^T p} ]

(a_i^T p > 0) 表示这条约束的边界正在被"逼近"，分母为正；如果 (a_i^T p \le 0)，表示方向在远离这条边界，不需要限制。最终步长取 (\alpha = \min(1, \alpha_{\max}))。如果 (\alpha) 被某条约束卡在小于1的位置，这条约束就是新触界的"拦截者"，要走过去就得把它加进工作集；如果 (\alpha = 1)，那就直接走到子问题最优点，工作集暂时不变，继续下一轮迭代。

3. 一个四步走通的手算例子

3.1 问题、初始点与工作集

光看算法流程还是抽象，我拿一个二维问题完整推演一遍，全部手算都能验算：

[ \min_x ; f(x) = \frac{1}{2}(x_1^2 + x_2^2) - x_1 - x_2 ]

约束： [ x_1 + x_2 \le 1, \quad x_1 \ge 0, \quad x_2 \ge 0 ]

矩阵形式写出来就是 (H = I)（单位阵），(g = (-1, -1))。无约束极小点在 ((1,1))，显然不满足 (x_1 + x_2 \le 1)，所以最优解一定落在约束边界上。几何上看，可行域是三角形区域，目标等值线是以 ((1,1)) 为圆心的同心圆，可行域内离 ((1,1)) 最近的点就是 ((0.5, 0.5))，这个点就是我们最终期望算到的结果。

取初始点 (x_0 = (0, 0.5))。这个点处 (x_1 = 0)，所以约束 (x_1 \ge 0) 是生效的，初始工作集设为 (W = {x_1 \ge 0})。

3.2 完整的迭代过程记录

我按每一轮迭代把关键量列成表，方便对照。

第0轮算出来 (p = (0, 0.5))，意思是保持 (x_1=0) 不动，只把 (x_2) 往上推。它先碰到的约束是 (x_1+x_2 \le 1)，因为从 (x_2=0.5) 出发沿该方向走到边界要走的步长正好是1，所以实际走了单位步长，到达 ((0,1))，并把这条约束加进工作集。这一步直观上也很好理解：既然无约束最优点在 ((1,1))，在当前点最想做的就是增大两个变量，被边界拦住是很自然的事。

第1轮是最关键的一轮。工作集里现在有 (x_1=0) 和 (x_1+x_2=1) 两条约束，当前点 ((0,1)) 是可行域的角点。解子问题得到 (p = (0,0))，说明在被这两条等式约束钉死的状态下，任何方向都无法改善目标。但检查乘子时发现，约束 (x_1 \ge 0) 对应的乘子 (\lambda = -1)，小于零。这个负乘子的几何含义，等下一节细说，操作上是把 (x_1 \ge 0) 这条约束从工作集里移除。

第2轮有意思的地方在于：移除约束后，工作集里只剩 (x_1+x_2 \le 1)，从 ((0,1)) 出发的方向变成 ((0.5, -0.5))。这正好是沿着斜边边界 (x_1+x_2=1) 往下滑的方向。步长计算时，(x_2 \ge 0) 这条约束允许走2步，但单位步长1更小，所以 (\alpha = 1)，新点落在 ((0.5, 0.5))，没有触碰新边界，工作集保持不变。

第3轮验证收尾：在 ((0.5, 0.5)) 处再解一次子问题，得到 (p=(0,0))，乘子 (\lambda = 0.5 > 0)，KKT条件全部满足，计算结束。

3.3 每一步背后在干什么

第1轮的负乘子是最容易困惑的地方。直观解释是这样的：当前点 ((0,1)) 同时被两条约束钉住，其中 (x_1=0) 这条约束其实在"帮倒忙"。如果把这条约束松绑，允许 (x_1) 稍微增大一点，目标函数沿哪个方向能下降？沿着斜边 (x_1+x_2=1) 朝 ((0.5,0.5)) 移动，目标从 (f(0,1) = -0.5) 降到 (f(0.5,0.5) = -0.75)，确实还有改进空间。乘子为负就是算法用数学语言告诉你：这条约束在当前点不是"支撑最优解"，而是"限制你得过头了"。

第2轮的方向 ((0.5,-0.5)) 也有清晰几何意义。((0.5,-0.5)) 正是无约束最优点 ((1,1)) 到斜边垂线的投影方向，所以沿这个方向走单位步长恰好落到最优投影点。而且这轮里虽然 (x_2 \ge 0) 约束还允许继续走，但单位步长已经到极限，因为子问题里沿 (p) 方向的最优步长就是1，再往下走目标函数会开始上升，所以停在这里是合理的。

还有一点值得体会：第1轮和第3轮都得到 (p=0)，但一轮判定"继续迭代"、一轮判定"收敛"，差异完全来自乘子符号。这说明积极集法里"求解子问题"只是半程，"检查乘子"才是判断是否到终点的关键。初学者最容易漏掉乘子检查，结果在 (p=0) 时直接认为已经最优，然后拿着错误答案去调试，怎么都想不通。

4. 手写一个能跑的Python实现

4.1 核心子程序：解KKT系统

前面的推导已经说明，每次迭代的核心是解一个KKT线性方程组。我用NumPy写一个最小可用的实现，不求性能最优，但逻辑清晰、可以直接改来用。

import numpy as np def solve_kkt(H, grad, A_w): """ 解等式约束子问题: min 0.5 * p' H p + grad' p s.t. A_w p = 0 返回 (p, lambda) """ n = H.shape[0] m = A_w.shape[0] KKT = np.block([ [H, A_w.T], [A_w, np.zeros((m, m))] ]) rhs = np.concatenate([-grad, np.zeros(m)]) sol = np.linalg.solve(KKT, rhs) return sol[:n], sol[n:]

这个函数的输入是当前迭代点处的梯度 (grad = Hx + g) 和工作集约束矩阵 (A_w)。返回值中 (p) 是搜索方向，(\lambda) 是对应工作集约束的拉格朗日乘子。注意到我直接用np.linalg.solve，要求KKT矩阵非奇异，这就是前面说的约束线性无关条件。严格来说在约束工作集里可能混入冗余约束导致奇异，完善的实现要先做秩检查，这里先不展开。

4.2 主循环：步长、触界、乘子检查

主循环按之前描述的算法流程组织：解子问题、判断 (p) 是否为零、计算步长、更新工作集。

def active_set_qp(H, g, A, b, x0, W0=None, tol=1e-9, max_iter=200): x = x0.copy() m = A.shape[0] # 初始工作集：默认取所有在当前点触界的约束 if W0 is None: W0 = [i for i in range(m) if abs(A[i] @ x - b[i]) < tol] W = list(W0) for _ in range(max_iter): grad = H @ x + g A_w = A[W, :] p, lam = solve_kkt(H, grad, A_w) # 情况1: 方向接近零，检查乘子 if np.linalg.norm(p) < tol: if len(W) > 0 and lam.min() < -tol: j_local = int(np.argmin(lam)) j = W.pop(j_local) continue else: return x, W # 情况2: 方向非零，计算最大可行步长 alpha = 1.0 for i in range(m): if i in W: continue den = A[i] @ p if den > tol: cand = (b[i] - A[i] @ x) / den alpha = min(alpha, cand) x_new = x + alpha * p # 新点处触界的约束加入工作集 touched = [ i for i in range(m) if i not in W and abs(A[i] @ x_new - b[i]) <= tol ] for i in touched: if i not in W: W.append(i) x = x_new raise RuntimeError("达到最大迭代次数，可能出现了退化或数值问题")

这里有一个实现细节值得说明：我在计算步长时，对所有不在工作集的约束都检查一遍分母 (a_i^T p)，只有为正时才可能构成限制。同时，在更新点之后统一检查哪些约束在新点处触界，一次性加入工作集，这样即使步长为1且有约束在终点恰好被碰到，也不会漏掉。

4.3 用刚才的例子验证代码

用第3节手算的例子做测试：

H = np.array([[1.0, 0.0], [0.0, 1.0]]) g = np.array([-1.0, -1.0]) # 约束顺序：x1+x2<=1, x1>=0, x2>=0 A = np.array([[1.0, 1.0], [-1.0, 0.0], [0.0, -1.0]]) b = np.array([1.0, 0.0, 0.0]) x0 = np.array([0.0, 0.5]) x_opt, W_opt = active_set_qp(H, g, A, b, x0, W0=[1]) print("最优解:", x_opt) print("最终工作集:", W_opt)

在我的环境里运行，输出是：

最优解: [0.5 0.5] 最终工作集: [0]

恰好对应手算推演的结果：(x^* = (0.5, 0.5))，最终生效的约束只有 (x_1+x_2 \le 1)。代码和手算互相验证，说明整个实现逻辑是自洽的。

5. 工程落地的实战经验

5.1 热启动是积极集法的杀手锏

真正在工程里用积极集法，最大的优势不是它原理简单，而是热启动（warm start）能力。在做MPC时，相邻两个控制周期求解的QP问题在结构上几乎一样，只是目标函数里的参考轨迹和测量值变了，而实际的约束激活模式往往变化很小。上一周期算完得到的最优解和最优点处的有效约束集合，直接作为下一周期的初始点 (x_0) 和初始工作集 (W_0)，通常只需要几次迭代就能收敛。

这个性质对嵌入式系统非常友好：每次迭代的耗时波动小，最坏情况可控，不像某些算法每轮迭代内部还有不定次数的回溯搜索。我在一个算力受限的控制器上把QP求解时间从几毫秒压到几百微秒，靠的就是热启动加提前把KKT矩阵的分解结果缓存起来。相比之下，内点法虽然迭代次数与规模关系不大，但每次迭代都要解一个更大的线性系统，而且暖启动的优势没有积极集法这么直接。

5.2 退化循环与数值容差

积极集法最经典的理论问题叫退化（degeneracy），表现出来就是算法在一个点附近反复加约束、踢约束，陷入死循环。最常见的原因是多个约束在同一个点同时触界，或者约束之间存在近似线性相关，导致乘子判定时出现平局，工作集在几条约束之间反复横跳。我在写数值实验时第一次遇到这个问题，卡在循环里出不来，打印工作集轨迹才发现某条约束被反复加进来又被踢出去。

工程上的应对手段通常有几个层次。第一，约束预处理，删掉重复或近似重复的约束，把变量和约束都做缩放，让数值在同一个量级。第二，判断规则上做文章，比如乘子为负时总选最负的那个，或者加一个很小的随机扰动打破平局。第三，设置最大迭代次数的兜底，一旦超限就触发一个降级策略。另外，(H) 只是半正定而非严格正定时，KKT矩阵可能奇异，常见的处理是给Hessian加一个很小的正则化项 (H + \epsilon I)，让问题变成严格凸，代价是解的精度会受影响，但稳定很多。

5.3 什么时候该换内点法或ADMM

积极集法不是万能的。它的迭代次数和约束切换次数强相关，如果问题规模大、约束模式变动剧烈，或者矩阵结构本身适合做稀疏分解，内点法和ADMM往往更有优势。我整理了一个简单的选型对比，直接给出结论性的经验：

比如在金融组合优化里，资产数量几千、协方差矩阵稀疏度一般，用内点法或OSQP更常见。而传统MPC领域，因为每个周期的问题小且热启动收益巨大，qpOASES这类积极集法求解器反而更流行。OSQP虽然基于ADMM，但在大规模稀疏场景下极快，也是现代MPC的热门选择。

5.4 对接现成库与两个调试技巧

如果不想自己实现，常见的现成库按接口易用度排序：quadprog是Goldfarb-Idnani型积极集法，适合中小规模稠密问题，接口极简；qpOASES专门为MPC的在线优化设计，支持热启动和动态问题；OSQP基于ADMM，适合大规模稀疏问题；scipy.optimize.minimize里的SLSQP和trust-constr也能解，但性能一般，适合原型验证。前提是确认库的底层算法和你对问题的预期匹配。

调试积极集法有两个技巧我屡试不爽。第一个是打印每一轮的工作集、(p)、(\lambda)、(\alpha)，特别是观察约束索引的序列是否在反复出现同一组索引——如果在，基本就是退化循环。第二个是调试前先把约束做归一化处理，把每条约束的 (a_i) 除以它的范数，这样步长和乘子的数值尺度一致，容差判断才不容易出偏差。我见过很多"看起来约束没问题但求解器老报不可行"的案例，最后查出来是变量单位不一致（毫米和米混用），数值差了一千倍，导致容差判断全部失效。

我一直觉得，积极集法虽然是个"老"方法，但它把约束激活模式的离散决策问题讲得最透彻。你一旦手写过一次，理解了KKT条件里互补性条件到底在说什么，后面再用内点法、ADMM或者其他任何先进求解器，遇到原理性问题时心里都有底。放在现在的嵌入式MPC项目里，我依然首选带热启动的积极集法，原因只有一个：它在实际问题上就是稳定、快速、可控。