简介：面向MATLAB优化学习者和工程应用人员的二次规划（QP）求解源码包，适合希望从代码层面理解约束优化算法实现的读者。资源围绕二次规划的标准形式展开，涉及Hessian矩阵、梯度向量及不等式/等式约束的设置，并展示quadprog函数与有效集法、路径跟踪法等核心求解思路，可直接用于工程优化、经济模型、信号处理等典型场景。压缩包共4个文件，全部为.m源文件，整体大小2KB，文件虽少但模块划分清晰，便于逐行研读和二次开发。已有293人浏览学习。通过研读这份源码，读者可掌握二次规划建模与算法迭代细节，了解自定义优化选项（如迭代限制、精度控制）对结果的影响，为后续扩展和迁移到实际问题提供扎实基础。 年初做投资组合权重优化的时候，我又一次被“二次规划”卡住了。需求说简单也简单：在总仓位等于1、不允许做空的条件下，找一组风险最小的资产权重。说难也难，因为打开MATLAB文档准备调用quadprog时，发现自己手里的数学式子跟函数要求的输入形式根本对不上。后来我把二次规划的数学模型、quadprog的输入输出，以及源代码转换这一步从头到尾捋了一遍，才把问题彻底搞明白。这篇笔记适合刚接触MATLAB优化的人，也适合已经调过quadprog但总在结果上栽跟头的工程师。核心就三件事：二次规划的标准形式是什么、怎么把实际问题写成MATLAB代码、以及求解失败时该从哪里下手。
1. 二次规划到底在解什么：先把数学形式对上号
1.1 那个 0.5 的约定，最容易栽跟头
二次规划的标准形式长这样：
min 0.5 * x' * H * x + f' * x s.t. A * x <= b Aeq * x = beq lb <= x <= ub其中H是二次项系数矩阵，f是一次项系数向量，A和b是不等式约束，Aeq和beq是等式约束，lb和ub是决策变量的上下界。这些都好理解，真正坑人的是目标函数开头的那个0.5。
quadprog约定的是最小化0.5*x'*H*x，但很多教材、论文里写的二次规划是min x'*H*x + f'*x，没有0.5。如果你从那些资料里直接抄一个H矩阵塞进quadprog，等于把目标函数放大了两倍。对于无约束问题，最优解位置可能不变，但一旦加上不等式或等式约束，目标函数缩放就会改变约束与目标之间的权衡，结果就偏了。
判断方法很简单：看H的来源。如果原始问题本身是min x'Σx这种形式，那转成quadprog时H要写成2*Σ；如果原始问题写成了min 0.5*x'Σx + c'x，H就是Σ。我建议所有代码里都在注释里写明“这里的H对应的是带0.5的标准形式还是不带0.5的形式”，不然过两周回来看代码铁定犯迷糊。
1.2 约束矩阵怎么对上号
二次规划的应用场景很广，投资组合里常见的是权重约束，模型预测控制里常见的是状态量、控制量的上下限约束，工程优化里常见的是资源总量限制。所有这些约束最终都要落到四种类型里：不等式、等式、上下界、以及没有约束。
写代码前先把下面这张对应表想清楚：
| 问题里的描述 | 需要转换成的形式 | 转换方式 |
|---|---|---|
| 所有变量之和等于1 | Aeq * x = beq | Aeq = ones(1, n), beq = 1 |
| 某个线性组合不小于阈值 | A * x <= b | 两边取负号，变成 -组合系数*x <= -阈值 |
| 变量在0到1之间 | lb <= x <= ub | lb = zeros(n, 1), ub = ones(n, 1) |
| 没有约束 | 传空数组 | 对应位置写 [] |
新手最容易出问题的就是“大于等于”约束。quadprog只接受小于等于，所以遇到mu' * w >= targetRet这种约束，必须两边乘负号，变成-mu' * w <= -targetRet。这一步丢失了符号，后面怎么调都白搭。
维度也要反复确认。A必须是m×n，b必须是m×1；Aeq是p×n，beq是p×1。MATLAB对维度不匹配的报错还算友好，但有些版本会把维度错误静默地转成另一种解释，导致结果不对还不报错，这个后面调试章节会细说。
2. 源代码实战：从最小二乘到带约束的组合优化
2.1 最小二乘问题怎么改写成二次规划
很多人不知道，最普通的最小二乘拟合其实就是二次规划的特例。考虑min ||C*x - d||^2，展开之后：
||C*x - d||^2 = x' * (C'*C) * x - 2 * d' * C * x + d' * d最后一项d'*d是常数，不影响最优解，可以丢掉。对照quadprog标准形式：
- H = 2 * (C' * C)
- f = -2 * (C' * d)
如果原始问题里写的是min 0.5 * ||C*x - d||^2，那H和f都要跟着缩放，但缩放后最优解x不变。这里可以写一小段代码验证：
% 最小二乘转二次规划验证 C = [1 2; 3 4]; d = [1; 1]; H = 2 * (C' * C); f = -2 * (C' * d); % 解析解：x = (C'*C)^{-1} * C' * d x_analytic = -H \ f; % quadprog 求解 options = optimoptions('quadprog', 'Display', 'off'); x_qp = quadprog(H, f, [], [], [], [], [], [], [], options); disp([x_analytic, x_qp]);如果这两个结果不一致，说明H和f构造有问题，而不是求解器有问题。这个例子我每次教别人调quadprog都会先用一遍，几秒钟就能定位出90%的模型转换错误。
2.2 一个完整可运行的组合权重优化示例
回到投资组合的例子。假设有5个资产，已知期望收益向量mu和协方差矩阵Sigma，想求解风险最小的权重组合，约束条件有三个：全仓、不允许做空、组合收益不低于目标收益率。
% 模拟数据：5个资产的期望收益 mu = [0.10; 0.12; 0.08; 0.15; 0.11]; % 模拟协方差矩阵（确保半正定） rng(0); A_rand = randn(5) * 0.2; Sigma = A_rand * A_rand' + eye(5) * 0.01; % 目标：最小化 w' * Sigma * w % 转成 quadprog 标准形式：H = 2*Sigma, f = 0 n = length(mu); H = 2 * Sigma; f = zeros(n, 1); % 全仓约束：sum(w) = 1 Aeq = ones(1, n); beq = 1; % 不允许做空：w >= 0，同时每只股票权重不超过1 lb = zeros(n, 1); ub = ones(n, 1); % 收益约束：mu' * w >= targetRet % 转成 A*w <= b： -mu' * w <= -targetRet targetRet = 0.10; A = -mu'; b = -targetRet; % 求解 options = optimoptions('quadprog', ... 'Algorithm', 'interior-point-convex', ... 'Display', 'iter'); [x, fval, exitflag, output] = quadprog(H, f, A, b, Aeq, beq, lb, ub, [], options);注意两点。第一，H为什么是2*Sigma，因为原始目标函数是w'*Sigma*w，quadprog标准形式里自带0.5，所以必须乘2。第二，收益约束为什么要对A和b取负号，因为quadprog只能处理小于等于不等式，取负号只是移项，不是改变约束本质。
跑完以后，x就是最优权重，fval是0.5*x'*H*x + f'*x，由于f是零向量，fval就等于x'*Sigma*x，也就是组合方差。想要组合波动率就取sqrt(fval)。这里是最容易误读的地方，quadprog返回的fval不是原始目标函数w'*Sigma*w本身吗？严格说fval是标准形式下的目标函数值，但因为H=2Σ，0.5*x'Hx = x'Σx，恰好等于原始目标，所以直接读fval没毛病。如果你把H写成Σ，那fval就不对应原始风险值了。
2.3 exitflag与output：结果到底靠不靠谱
常看exitflag的取值，大概分几种情况：
| exitflag | 含义 | 常见处理方式 |
|---|---|---|
| 1 | 找到最优解，满足收敛条件 | 正常 |
| 0 | 达到最大迭代次数，可能还没收敛 | 调大MaxIterations，或换算法 |
| -2 | 找不到可行点，约束之间矛盾 | 检查约束是否过紧、是否自相矛盾 |
| -3 | 问题无界 | 检查lb/ub，以及H是否半正定 |
| -4 | 求解过程遇到NaN或Inf | 检查矩阵中是否有异常值，数据是否脏 |
output结构体里值得多看两个字段：iterations和firstorderopt。前者告诉你了迭代多少步，后者是一阶最优性条件的度量，理论上接近0才是好的。如果exitflag=0但firstorderopt已经很小，说明其实离最优解很近了，只是没有达到默认的OptimalityTolerance，这时候把容差稍微放松一点就能拿到结果。
3. 结果离谱时，完整的调试排查链路
3.1 第一件事：用解析解验证模型，而不是怀疑求解器
quadprog本身是经过大量测试的，绝大多数“结果不对”都是模型到代码的转换问题。所以我的调试顺序很固定：先跑一个没有约束的小问题，拿解析解和quadprog结果对比。比如min 0.5*x'Hx + f'x，无约束最优解就是x = -H\f。
这一步能快速判断H和f是否构造正确。如果解析解对上了，再去检查约束部分；如果连解析解都对不上，问题一定出在H/f转换，和后端求解器没有关系。别一上来就调各种tolerance，那是在错误方向上浪费时间。
3.2 迭代不收敛：先看算法，再看迭代上限
exitflag=0对应的是“迭代没跑完但触顶了”。我见过很多人在这种情况下直接把MaxIterations从默认值调大，但结果还是老样子。更有效的做法是先把Display设为'iter'，看每一轮的目标函数值和一阶最优性指标有没有在下降。如果每轮都在下降，只是下降得慢，那调大MaxIterations有意义；如果几轮之后目标函数基本不动了，说明算法已经到极限，这时候应该调整的是约束的数值尺度，或者换一个算法。
quadprog在凸二次规划上常用的算法是interior-point-convex，它默认自带一些预处理。老版本或者特殊场景下也可以选active-set。active-set在中小规模问题上迭代次数通常更稳定，但它要求H必须是正定的，如果有零特征值或者负特征值，会遇到问题。遇到迭代不收敛，我可以先跑一个关闭约束的版本，再把约束逐步加回来，定位到底是哪个约束引起的。
3.3 数值尺度不一致：一个小权重，一把辛酸泪
实际项目里的数据不会像官方示例那么整洁。比如协方差矩阵的元素数量级可能只有10^-4，而收益率约束是0.1，目标函数和约束条件之间尺度差了好几个数量级。interior-point算法对尺度很敏感，矩阵条件数太差时，即使问题本身有解，也可能算出非常离谱的权重。
处理办法有两个。一是对数据做归一化，比如把收益和协方差都乘同一个倍数，让它们的数量级落在1附近。二是设置合理的ConstraintTolerance和OptimalityTolerance。默认的tolerance不一定适配你的问题尺度，但也不要随便放大，放大太多会得到一个“看着可行其实违反约束”的伪最优解。
另外，如果协方差矩阵是用历史数据算出来的，经常会有很小的负特征值，这是数值噪声导致的，不是真的非凸。先eig(Sigma)看一眼，如果负特征值接近机器精度，可以直接用Sigma = (Sigma + Sigma') / 2 + 1e-8 * eye(n)做对称化加微扰，保证半正定。
3.4 无可行解：约束自相矛盾的排查顺序
exitflag=-2是最让人头疼的，因为它不是算法不好，而是可行域本身就是空的。也就是说，你给出的约束条件里，根本没有一个点能满足所有限制。
排查原则：先从最简单的约束组合开始，逐步加约束。比如先只保留全仓约束和上下界，看有没有解；再加上收益约束，再看有没有解。哪一步开始报不可行，问题就出在哪一步。还可以用一个最小可行性的线性规划来验证：
% 用 linprog 检查可行域是否为空 f_feas = zeros(n, 1); [x_feas, ~, exitflag_feas] = linprog(f_feas, A, b, Aeq, beq, lb, ub);如果exitflag_feas不是1，说明约束本身就矛盾了。最常见的矛盾是哪几种？lb和ub重叠，例如要求权重在0到0.5之间，又要求权重之和等于1，显然不可能。还有就是“大于等于”约束取负号时没转干净，比如本该是-mu'*w <= -targetRet，写成了mu'*w <= -targetRet，方向反了，可行域自然可能变成空集。
4. 进阶边界：什么时候该换其他求解思路
4.1 大规模稀疏问题：quadprog的极限在哪里
quadprog对于几千个变量、几千个约束的稠密问题基本能扛住，但规模再往上走，或者约束矩阵具有明显稀疏结构，就要注意内存和计算效率了。interior-point-convex算法内部需要解一个大型线性方程组，H和A如果不是稀疏存储，内存会先撑不住。
一个有效做法是显式把矩阵转成稀疏类型再传给quadprog：
H_sparse = sparse(H); A_sparse = sparse(A);有了稀疏格式，quadprog内部会走稀疏线性代数路径，能处理的规模会大很多。再往上，如果变量上万，约束结构又很复杂，MATLAB生态里还可以接Gurobi、Mosek这类专业求解器，它们在预处理、并行计算方面更强。纯MATLAB环境下，也有osqp这样的ADMM求解器，但需要自己处理QP标准形式的转换，学习成本另算。
4.2 非凸二次规划与整数变量：quadprog解决不了的部分
quadprog默认要求问题至少是凸的，更准确地说，H必须是对称半正定矩阵。判断方法是：
eig(H)只要出现一个明显的负特征值，就不适合直接扔给quadprog。此时quadprog可能报错，也可能给出某个局部解，但无法保证是全局最优。对小型非凸问题，可以考虑枚举或分支定界；对中型问题，可以试试fmincon换个非线性求解器，但对于二次目标它也只是找局部解，需要根据实际场景判断这个局部解接不接受。
还有一类更特殊的问题，变量要求取整数，比如“哪些资产纳入组合，用0/1变量表示”。这种带整数约束的二次规划叫MIQP，MATLAB自带的intlinprog只能处理线性目标函数，quadprog又不接受整数约束，两边都搭不上。到了这一步就需要上商业求解器或者专门的启发式算法了。所以遇到MIQP，第一反应不应该是找“matlab怎么解”，而是先评估问题规模，再决定是用Gurobi这类工具还是自己写启发式。
最后分享一个我自己的习惯。无论问题看起来多简单，我都会先用一个可以手算的小样例跑通，再换真实数据。quadprog把矩阵运算和迭代细节封装得很干净，也正因如此，一旦结果不对，问题往往出在“标准形式”的转换上，而不是求解器本身。拿这篇文章里的源代码框架，换成你自己的C、d、A、b参数，基本半小时就能跑出第一版结果。真正的难点从来不是调用函数，而是建模时把每一个约束都翻译成quadprog认识的矩阵和向量。
本文还有配套的精品资源，点击获取
