1. 从一道经典例题切入非线性规划到底是什么如果你学过线性规划可能会觉得那套方法挺“规矩”的——目标函数和约束条件都是线性的图形是直线或平面最优解总在顶点上。但现实世界可没这么“线性”。比如你想设计一个容积最大的圆柱形易拉罐在给定表面积材料成本的限制下它的半径和高是多少这里容积是半径和高的乘积非线性表面积也是半径和高的函数非线性。再比如投资组合优化中收益可能是风险的二次函数为了平衡风险与回报约束条件里可能还有各种非线性的法规限制。这些问题线性规划就束手无策了得请出我们今天的主角非线性规划。简单说非线性规划就是研究目标函数或约束条件中至少有一个是非线性函数时的最优化问题。它的数学形式一般长这样min f(x) s.t. g_i(x) ≤ 0, i 1, ..., m h_j(x) 0, j 1, ..., p其中f(x),g_i(x),h_j(x)中至少有一个是非线性的。x是我们的决策变量向量。min f(x)表示最小化目标函数s.t.是“subject to”的缩写代表约束条件。g_i(x) ≤ 0是不等式约束h_j(x) 0是等式约束。非线性规划之所以重要且“难搞”是因为它的解空间可能非常复杂不再是单纯的多面体。最优解可能出现在可行域的内部此时约束可能不起作用称为内部解也可能出现在边界上此时一个或多个约束取等号称为边界解。更麻烦的是目标函数可能有多个“洼地”局部最优解我们的任务是找到最深的那个“洼地”全局最优解但很多算法容易陷在某个局部洼地里出不来。为了让大家有个直观感受我们从一个经典的、结构清晰的例题开始。这个例子将贯穿全文我会手把手带你用 MATLAB 把它解出来并在过程中解释每一步背后的原理和考量。例题投资组合优化简化版假设你有两种资产可供投资股票A和股票B。你希望分配你的资金比如1个单位以最大化你的期望收益但同时要控制风险。我们用一个非常经典的模型来简化决策变量x1投资于股票A的比例x2投资于股票B的比例。目标函数最大化期望收益f(x) 0.1*x1 0.15*x2。这里假设股票A的期望收益率是10%股票B是15%。注意这个目标函数是线性的但别急非线性马上就来。风险约束非线性 我们用收益的方差来衡量风险。假设两种股票收益的方差和协方差已知投资组合的方差风险为Risk 0.2*x1^2 0.1*x2^2 0.05*x1*x2。我们要求风险不能超过一个阈值比如 0.1。于是得到第一个非线性不等式约束0.2*x1^2 0.1*x2^2 0.05*x1*x2 ≤ 0.1。预算约束线性 你的投资比例之和必须为1x1 x2 1。非负约束线性 投资比例不能为负不允许做空x1 ≥ 0,x2 ≥ 0。所以完整的数学模型是max 0.1*x1 0.15*x2 s.t. 0.2*x1^2 0.1*x2^2 0.05*x1*x2 ≤ 0.1 x1 x2 1 x1 ≥ 0, x2 ≥ 0这是一个典型的目标函数线性、但包含非线性约束的非线性规划问题。我们最终要找到满足风险要求下能使收益最大化的x1和x2。2. 解题工具箱MATLAB中的fmincon函数详解面对非线性规划我们不可能每次都徒手推导解析解那会非常复杂。在数学建模和工程实践中我们高度依赖数值优化工具。MATLAB 的优化工具箱Optimization Toolbox提供了强大的求解器其中fmincon函数是解决中等规模非线性规划问题的“瑞士军刀”。fmincon这个名字可以拆解为 “Find MINimum of a CONstrained nonlinear multivariable function”即寻找有约束非线性多元函数的最小值。记住MATLAB 的优化函数默认都是最小化目标函数。如果你的问题是最大化比如我们的例题只需要将目标函数乘以 -1 即可。2.1fmincon的基本调用语法最完整的调用形式如下[x, fval, exitflag, output, lambda, grad, hessian] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, nonlcon, options)看起来参数很多别怕我们一个个拆解大部分在实际应用中并不需要同时使用。输出参数x求解得到的最优决策变量值。fval在最优解x处的目标函数值。注意对于最大化问题这是-1 * 原始目标函数在x处的值。exitflag算法终止原因的整数代码。这是最重要的诊断信息之一正数通常表示成功例如1 表示一阶最优性条件满足0 表示达到最大迭代次数或函数计算次数负数表示求解失败例如-2 表示找不到可行解。output一个结构体包含优化过程的详细信息如迭代次数、函数计算次数、算法类型、一阶最优性度量等。lambda在最优解处的拉格朗日乘子Lagrange multipliers结构体。lambda.ineqlin对应线性不等式约束lambda.eqlin对应线性等式约束lambda.ineqnonlin对应非线性不等式约束lambda.eqnonlin对应非线性等式约束。乘子的大小可以解释约束的“紧度”或“价值”。grad在最优解处目标函数的梯度可选输出。hessian在最优解处目标函数的黑塞矩阵Hessian可选输出。输入参数fun目标函数句柄。一个 MATLAB 函数文件或匿名函数输入是决策变量向量x输出是标量目标函数值f。x0初始猜测点。这是非线性规划求解的起点非常重要不同的起点可能导致算法收敛到不同的局部最优解。A, b定义线性不等式约束A*x ≤ b。如果没有用空数组[]代替。Aeq, beq定义线性等式约束Aeq*x beq。如果没有用空数组[]代替。lb, ub定义决策变量的下界和上界lb ≤ x ≤ ub。如果没有用空数组[]代替。nonlcon非线性约束函数句柄。一个 MATLAB 函数输入x输出两个向量非线性不等式约束c(x) ≤ 0和非线性等式约束ceq(x) 0。这里必须写成 ≤0 和 0 的形式options优化选项结构体用于设置算法参数如最大迭代次数、函数容差、显示输出等。通常用optimoptions(fmincon, ...)来创建。2.2 算法选择fmincon的四种内功心法fmincon内置了多种算法通过options中的Algorithm选项指定。选择合适的算法对求解效率和成功率至关重要。interior-point内点法默认算法原理通过在可行域内部构造一条路径逼近最优解。它会在目标函数中加入一个“障碍函数”来惩罚接近边界的点随着迭代障碍的影响逐渐减小路径最终收敛到边界上的最优解如果存在的话。优点对于大规模问题变量和约束多通常很有效尤其擅长处理不等式约束。它通常能很好地从可行域内部开始搜索。缺点对于某些非凸问题可能不如基于梯度的算法精确。适用场景大多数中大规模、包含不等式约束的问题的默认选择。sqp序列二次规划法原理在每一步迭代它用二次函数近似目标函数用线性函数近似约束从而构造一个二次规划子问题。求解这个子问题得到搜索方向然后沿此方向进行线搜索。优点通常能提供非常高的精度尤其对于光滑函数。它通常能很好地处理等式约束和活跃的不等式约束。缺点对于大规模问题计算每个子问题的成本可能较高。适用场景中小规模、需要高精度解的问题特别是等式约束较多的问题。active-set有效集法原理它猜测哪些不等式约束在最优解处是“活跃的”即取等号然后只把这些约束当作等式约束来处理在由这些活跃约束定义的“面”上进行搜索。如果猜测错了就调整活跃集。优点对于中小规模问题特别是解位于许多约束边界的情况可能很有效。缺点对于大规模问题维护和更新活跃集的成本很高。适用场景传统算法现在更多被interior-point和sqp取代但在某些特定问题中仍有价值。trust-region-reflective信赖域反射法原理在每一步它在当前点的一个可信赖区域内用较简单的模型如二次模型近似原问题然后在这个小区域内求解子问题。它要求目标函数能提供梯度信息。优点对于大规模无约束或边界约束问题以及某些特殊的线性等式约束问题非常高效。缺点不能处理一般的非线性约束。只能处理边界约束lb, ub和线性等式约束Aeq, beq。适用场景大规模、目标函数可求导、只有边界或线性等式约束的问题。选择建议对于初学者和大多数建模问题直接使用默认的interior-point算法即可。如果求解失败或精度不够可以尝试切换到sqp算法。如果问题只有边界约束可以尝试trust-region-reflective。3. 手把手实现例题的MATLAB代码与逐行解析现在我们回到最初的例题用 MATLAB 和fmincon来求解。我会把代码分成几个部分并详细解释每一行。3.1 第一步定义目标函数我们的原始问题是最大化0.1*x1 0.15*x2。由于fmincon默认最小化我们需要将其转换为最小化-1*(0.1*x1 0.15*x2)。我们可以使用匿名函数来简洁地定义% 定义目标函数注意fmincon求解最小值所以最大化问题要加负号 objective (x) - (0.1*x(1) 0.15*x(2));这里(x)创建了一个匿名函数输入变量是xx(1)对应x1x(2)对应x2。函数体计算-(0.1*x1 0.15*x2)。3.2 第二步定义线性约束我们的线性约束包括一个等式约束和两个不等式约束非负约束。等式约束x1 x2 1。对应Aeq * x beq其中Aeq [1, 1],beq 1。不等式约束非负x1 ≥ 0,x2 ≥ 0。这属于变量的下界约束我们用lb参数处理会更方便。所以这里A和b为空。变量下界lb [0; 0]。变量上界没有明确上界设为空[]。% 线性等式约束 Aeq * x beq Aeq [1, 1]; beq 1; % 线性不等式约束 A * x b 本例中没有用空数组 A []; b []; % 变量的下界和上界 lb [0; 0]; % x1 0, x2 0 ub []; % 无上界3.3 第三步定义非线性约束这是关键的一步。非线性约束函数需要单独写成一个函数文件或嵌套函数。它必须返回两个输出c和ceq分别表示非线性不等式约束c(x) ≤ 0和非线性等式约束ceq(x) 0。我们的风险约束是0.2*x1^2 0.1*x2^2 0.05*x1*x2 ≤ 0.1。 为了符合c(x) ≤ 0的形式我们将其改写为0.2*x1^2 0.1*x2^2 0.05*x1*x2 - 0.1 ≤ 0。 本例中没有非线性等式约束所以ceq返回空数组[]。我们创建一个独立的函数文件nonlinear_constraints.mfunction [c, ceq] nonlinear_constraints(x) % 非线性不等式约束 c(x) 0 % 风险约束: 0.2*x1^2 0.1*x2^2 0.05*x1*x2 0.1 c 0.2*x(1)^2 0.1*x(2)^2 0.05*x(1)*x(2) - 0.1; % 非线性等式约束 ceq(x) 0 本例中没有 ceq []; end或者如果你希望所有代码在一个脚本里可以使用嵌套函数或函数句柄。这里展示函数句柄的写法更简洁% 使用函数句柄定义非线性约束 nonlcon (x) deal(0.2*x(1)^2 0.1*x(2)^2 0.05*x(1)*x(2) - 0.1, []);deal函数将两个输出分别赋给c和ceq。第一个参数是c第二个是ceq。3.4 第四步设置初始点与选项并调用fmincon初始点x0的选择很重要。它应该在可行域内或附近以帮助算法更快收敛。对于本例一个简单的选择是等额投资x0 [0.5; 0.5]。我们也可以设置一些优化选项比如显示迭代过程。% 初始猜测点 (必须在可行域内或附近) x0 [0.5; 0.5]; % 假设初始投资各一半 % 设置优化选项可选但推荐 options optimoptions(fmincon, ... Display, iter, ... % 显示每次迭代的信息 Algorithm, interior-point); % 选择算法默认就是它这里显式指定 % 调用 fmincon 求解 [x_opt, fval_opt, exitflag, output] fmincon(objective, x0, A, b, Aeq, beq, lb, ub, nonlcon, options);3.5 第五步结果解读与后处理求解完成后我们需要解读输出。% 显示最优解和最优值 fprintf(最优投资比例\n); fprintf( 股票A (x1) %.4f\n, x_opt(1)); fprintf( 股票B (x2) %.4f\n, x_opt(2)); fprintf( 总投资比例和 %.4f (应等于1)\n, sum(x_opt)); % 注意fval_opt 是转换后目标函数即 -收益的最小值 % 所以最大收益 -fval_opt max_return -fval_opt; fprintf(\n最大期望收益 %.4f\n, max_return); % 验证约束是否满足 risk 0.2*x_opt(1)^2 0.1*x_opt(2)^2 0.05*x_opt(1)*x_opt(2); fprintf(实际组合风险 %.6f (约束要求 ≤ 0.1)\n, risk); fprintf(风险约束满足情况%s\n, (risk 0.11e-6) ? 满足 : 不满足); % 加一个小容差 % 检查退出标志 fprintf(\n退出标志 exitflag %d\n, exitflag); if exitflag 0 fprintf(优化成功收敛。\n); elseif exitflag 0 fprintf(达到最大迭代次数或函数计算次数。\n); else fprintf(优化未成功。\n); end % 查看输出信息 fprintf(迭代次数%d\n, output.iterations); fprintf(函数计算次数%d\n, output.funcCount);运行这段完整的代码你可能会得到类似如下的输出具体数值可能因算法和版本略有差异Iter Func-count Fval Feasibility Step Length Norm of First-order optimality 0 3 -0.125000 0.000000e00 1.000e00 0.000e00 1.250e-01 1 6 -0.136364 0.000000e00 1.000e00 2.500e-01 1.136e-01 2 9 -0.137931 0.000000e00 1.000e00 6.250e-02 1.034e-01 3 12 -0.138889 0.000000e00 1.000e00 1.562e-02 9.722e-02 ... (更多迭代) 8 27 -0.139241 0.000000e00 1.000e00 9.537e-07 9.537e-07 优化已终止: 一阶最优性度量小于 options.OptimalityTolerance 并且约束违反值小于 options.ConstraintTolerance。 最优投资比例 股票A (x1) 0.4483 股票B (x2) 0.5517 总投资比例和 1.0000 (应等于1) 最大期望收益 0.1392 实际组合风险 0.100000 (约束要求 ≤ 0.1) 风险约束满足情况满足 退出标志 exitflag 1 优化成功收敛。 迭代次数8 函数计算次数27结果分析最优解x1 ≈ 0.4483,x2 ≈ 0.5517。这意味着在风险控制下应将约44.83%的资金投入股票A55.17%投入收益更高的股票B。最大收益0.1392即13.92%介于股票A的10%和股票B的15%之间符合预期。约束验证风险恰好等于0.1说明风险约束是“活跃的”在最优解处取等号它有效地限制了我们对高收益但高风险资产B的过度配置。预算约束也严格满足。算法状态exitflag1表示成功收敛一阶最优性条件得到满足。迭代了8次计算了27次目标/约束函数效率很高。4. 避坑指南与实战经验为什么我的fmincon跑不出结果在实际使用中你可能会遇到fmincon报错、不收敛、或者结果明显不合理的情况。下面我总结了一些最常见的“坑”及其解决办法。4.1 初始点x0的选择一个好的开始是成功的一半fmincon对初始点非常敏感尤其是对于非凸问题有多个局部最优解。问题如果x0离可行域太远或者导致初始约束违反严重算法可能一开始就失败。对策物理意义法根据问题的实际背景给出一个合理的猜测。比如投资问题用均匀分配工程设计问题用经验值。随机初始化多次尝试对于复杂问题可以在可行域内或大致范围内随机生成多个初始点分别运行fmincon然后选择目标函数值最好的那个解。best_x []; best_fval inf; num_trials 20; for i 1:num_trials x0_rand rand(2,1); % 生成[0,1]随机数 x0_rand x0_rand / sum(x0_rand); % 归一化以满足等式约束 [x_temp, fval_temp] fmincon(objective, x0_rand, A, b, Aeq, beq, lb, ub, nonlcon); if fval_temp best_fval best_fval fval_temp; best_x x_temp; end end可行性优先如果找不到明显可行的点可以先用一个简单的优化问题比如最小化约束违反度来找到一个可行的初始点。4.2 非线性约束函数的定义格式错误是万恶之源这是新手最容易出错的地方。坑1不等式和等式的顺序。函数必须返回[c, ceq]c对应c(x) ≤ 0ceq对应ceq(x) 0。顺序反了约束意义就全错了。坑2约束方向。必须写成≤0和0的形式。如果你的约束是g(x) ≥ 0需要转化为-g(x) ≤ 0。坑3向量化输出。即使只有一个非线性约束c和ceq也必须是列向量。例如有两个非线性不等式c1(x)≤0和c2(x)≤0应该写成c [c1(x); c2(x)]。坑4函数未在路径中。如果使用独立的.m文件定义非线性约束函数确保该文件保存在当前工作目录或 MATLAB 搜索路径中。检查清单在调用fmincon前用初始点x0手动调用一次你的非线性约束函数检查输出是否符合预期。[c0, ceq0] nonlcon(x0); disp(初始点处的非线性不等式约束值 c(x0):); disp(c0); disp(初始点处的非线性等式约束值 ceq(x0):); disp(ceq0);确保c0的各分量 ≤ 0或接近0ceq0的各分量 0或接近0。如果c0有正的大数说明初始点严重违反约束。4.3 算法不收敛与参数调优给fmincon一点耐心和指导有时算法会因达到最大迭代次数而停止exitflag0或者根本找不到可行解exitflag-2。对策1调整options。options optimoptions(fmincon, ... Display, iter-detailed, ... % 显示更详细的迭代信息 MaxIterations, 1000, ... % 增加最大迭代次数默认通常是400 MaxFunctionEvaluations, 3000, ... % 增加最大函数计算次数 OptimalityTolerance, 1e-8, ... % 收紧一阶最优性容差要求更高精度 ConstraintTolerance, 1e-8, ... % 收紧约束违反容差 StepTolerance, 1e-10); % 收紧步长容差iter-detailed显示模式能帮你看到每次迭代的目标值、约束违反度和一阶最优性度量是诊断问题的利器。如果看到目标值很久不下降或者约束违反度降不下来可能就是问题本身或初始点有问题。对策2尝试不同算法。如果interior-point失败了切换到sqp试试。options optimoptions(fmincon, Algorithm, sqp, Display, iter);对策3提供解析导数梯度、雅可比矩阵。默认情况下fmincon用有限差分法数值估算导数这既慢又不准。如果你能为目标函数和非线性约束提供解析导数能极大提升速度、精度和稳定性。通过options设置SpecifyObjectiveGradient为true并让目标函数返回两个输出[f, gradf]。设置SpecifyConstraintGradient为true并让非线性约束函数返回四个输出[c, ceq, gradc, gradceq]。 这需要一定的数学推导但对于复杂问题收益巨大。4.4 问题尺度与数值稳定性别让计算机算“糊涂”了如果决策变量的数量级相差巨大比如x1在 1e6 级别x2在 1e-3 级别或者目标函数的数值非常大/小都会导致数值计算困难出现“病态”问题。对策尺度缩放。尽量通过变量代换让所有决策变量都在O(1)的数量级附近。例如如果x1代表以米为单位的长度范围是几千可以令x1_scaled x1 / 1000将单位变为千米。在模型中相应地调整目标函数和约束系数。求解出x_scaled后再变换回原始变量x。4.5 全局最优与局部最优跳出“洼地”fmincon找到的通常是局部最优解。对于非凸问题这可能不是全局最优。对策1多起点优化。如前所述从多个随机初始点出发求解取最好的结果。这是最实用、最常用的方法。对策2使用全局优化算法。对于确实难以找到全局最优的问题可以考虑 MATLAB 的全局优化工具箱中的函数如GlobalSearch或MultiStart。它们会在fmincon的基础上自动进行多起点搜索。problem createOptimProblem(fmincon, objective, objective, x0, x0, ... Aeq, Aeq, beq, beq, lb, lb, ub, ub, nonlcon, nonlcon); gs GlobalSearch; [x_global, fval_global] run(gs, problem);这相当于一个更系统、更智能的“多起点”方法。5. 举一反三非线性规划建模的常见类型与扩展掌握了基本解法后我们来看看非线性规划在数学建模中还有哪些常见形态和扩展思路。5.1 无约束优化fminunc与fminsearch如果问题没有约束就是无约束非线性优化。MATLAB 提供了专门的函数fminunc适用于光滑函数可以使用梯度信息。用法与fmincon类似但不需要约束参数。fminsearch使用 Nelder-Mead 单纯形法不需要导数信息适用于非光滑或求导困难的函数但效率较低变量不宜过多。5.2 最小二乘问题lsqnonlin与lsqcurvefit当目标函数是若干项的平方和时例如曲线拟合问题就构成了非线性最小二乘问题min Σ [F_i(x)]^2。lsqnonlin求解非线性最小二乘问题。lsqcurvefit专门用于数据拟合目标是使模型曲线F(x, xdata)与观测数据ydata的误差平方和最小。这类问题有特殊的结构使用这些专用函数通常比通用的fmincon更高效、更稳定。5.3 多目标优化帕累托最优解集现实中很多问题需要同时优化多个相互冲突的目标如成本最低、性能最好、风险最小。这没有单一的最优解而是一组“帕累托最优解”Pareto optimal solutions即在不使任何一个目标变差的情况下无法再改进其他目标。 MATLAB 的全局优化工具箱提供了gamultiobj函数使用遗传算法来求解多目标优化问题并绘制帕累托前沿Pareto front。5.4 整数/离散变量混合整数非线性规划如果部分决策变量必须是整数如选择设备的台数、路径规划中的节点选择问题就变成了混合整数非线性规划。这类问题难度急剧上升。MATLAB 的全局优化工具箱提供了ga遗传算法和surrogateopt代理优化等函数可以处理这类问题但它们属于启发式算法不能保证找到全局最优解。对于特定结构的问题也可以使用专门的商业求解器如 Gurobi、BARON通过 MATLAB 接口调用。5.5 从例题到实际建模思想迁移我们用的投资组合例题是高度简化的。实际中的马科维茨均值-方差模型目标函数可能是最大化夏普比率收益/风险约束可能包括行业配置上限、不允许做空lb0、最低持仓比例等形成一个更复杂的非线性规划。但核心建模步骤不变定义决策变量各资产的投资权重。建立目标函数最大化收益、最小化风险或最大化风险调整后收益。确定约束条件预算约束、风险约束、法规约束等线性和非线性。选择求解器并实现根据问题规模资产数量和约束类型选用fmincon或其他工具。结果分析与验证检查解是否满足所有约束是否合理并进行敏感性分析例如改变风险上限观察最优投资组合如何变化。这种“定义变量-建立目标与约束-数值求解-分析结果”的框架是解决绝大多数优化类建模问题的通用心法。无论是设计工程结构、调度生产资源还是训练机器学习模型其损失函数最小化本质上也是一个无约束/有约束优化问题其底层逻辑都是相通的。理解并熟练运用fmincon这样的工具就等于掌握了一把打开众多实际问题的钥匙。 MATLAB fmincon函数详解：从非线性规划原理到投资组合优化实战
1. 从一道经典例题切入非线性规划到底是什么如果你学过线性规划可能会觉得那套方法挺“规矩”的——目标函数和约束条件都是线性的图形是直线或平面最优解总在顶点上。但现实世界可没这么“线性”。比如你想设计一个容积最大的圆柱形易拉罐在给定表面积材料成本的限制下它的半径和高是多少这里容积是半径和高的乘积非线性表面积也是半径和高的函数非线性。再比如投资组合优化中收益可能是风险的二次函数为了平衡风险与回报约束条件里可能还有各种非线性的法规限制。这些问题线性规划就束手无策了得请出我们今天的主角非线性规划。简单说非线性规划就是研究目标函数或约束条件中至少有一个是非线性函数时的最优化问题。它的数学形式一般长这样min f(x) s.t. g_i(x) ≤ 0, i 1, ..., m h_j(x) 0, j 1, ..., p其中f(x),g_i(x),h_j(x)中至少有一个是非线性的。x是我们的决策变量向量。min f(x)表示最小化目标函数s.t.是“subject to”的缩写代表约束条件。g_i(x) ≤ 0是不等式约束h_j(x) 0是等式约束。非线性规划之所以重要且“难搞”是因为它的解空间可能非常复杂不再是单纯的多面体。最优解可能出现在可行域的内部此时约束可能不起作用称为内部解也可能出现在边界上此时一个或多个约束取等号称为边界解。更麻烦的是目标函数可能有多个“洼地”局部最优解我们的任务是找到最深的那个“洼地”全局最优解但很多算法容易陷在某个局部洼地里出不来。为了让大家有个直观感受我们从一个经典的、结构清晰的例题开始。这个例子将贯穿全文我会手把手带你用 MATLAB 把它解出来并在过程中解释每一步背后的原理和考量。例题投资组合优化简化版假设你有两种资产可供投资股票A和股票B。你希望分配你的资金比如1个单位以最大化你的期望收益但同时要控制风险。我们用一个非常经典的模型来简化决策变量x1投资于股票A的比例x2投资于股票B的比例。目标函数最大化期望收益f(x) 0.1*x1 0.15*x2。这里假设股票A的期望收益率是10%股票B是15%。注意这个目标函数是线性的但别急非线性马上就来。风险约束非线性 我们用收益的方差来衡量风险。假设两种股票收益的方差和协方差已知投资组合的方差风险为Risk 0.2*x1^2 0.1*x2^2 0.05*x1*x2。我们要求风险不能超过一个阈值比如 0.1。于是得到第一个非线性不等式约束0.2*x1^2 0.1*x2^2 0.05*x1*x2 ≤ 0.1。预算约束线性 你的投资比例之和必须为1x1 x2 1。非负约束线性 投资比例不能为负不允许做空x1 ≥ 0,x2 ≥ 0。所以完整的数学模型是max 0.1*x1 0.15*x2 s.t. 0.2*x1^2 0.1*x2^2 0.05*x1*x2 ≤ 0.1 x1 x2 1 x1 ≥ 0, x2 ≥ 0这是一个典型的目标函数线性、但包含非线性约束的非线性规划问题。我们最终要找到满足风险要求下能使收益最大化的x1和x2。2. 解题工具箱MATLAB中的fmincon函数详解面对非线性规划我们不可能每次都徒手推导解析解那会非常复杂。在数学建模和工程实践中我们高度依赖数值优化工具。MATLAB 的优化工具箱Optimization Toolbox提供了强大的求解器其中fmincon函数是解决中等规模非线性规划问题的“瑞士军刀”。fmincon这个名字可以拆解为 “Find MINimum of a CONstrained nonlinear multivariable function”即寻找有约束非线性多元函数的最小值。记住MATLAB 的优化函数默认都是最小化目标函数。如果你的问题是最大化比如我们的例题只需要将目标函数乘以 -1 即可。2.1fmincon的基本调用语法最完整的调用形式如下[x, fval, exitflag, output, lambda, grad, hessian] fmincon(fun, x0, A, b, Aeq, beq, lb, ub, nonlcon, options)看起来参数很多别怕我们一个个拆解大部分在实际应用中并不需要同时使用。输出参数x求解得到的最优决策变量值。fval在最优解x处的目标函数值。注意对于最大化问题这是-1 * 原始目标函数在x处的值。exitflag算法终止原因的整数代码。这是最重要的诊断信息之一正数通常表示成功例如1 表示一阶最优性条件满足0 表示达到最大迭代次数或函数计算次数负数表示求解失败例如-2 表示找不到可行解。output一个结构体包含优化过程的详细信息如迭代次数、函数计算次数、算法类型、一阶最优性度量等。lambda在最优解处的拉格朗日乘子Lagrange multipliers结构体。lambda.ineqlin对应线性不等式约束lambda.eqlin对应线性等式约束lambda.ineqnonlin对应非线性不等式约束lambda.eqnonlin对应非线性等式约束。乘子的大小可以解释约束的“紧度”或“价值”。grad在最优解处目标函数的梯度可选输出。hessian在最优解处目标函数的黑塞矩阵Hessian可选输出。输入参数fun目标函数句柄。一个 MATLAB 函数文件或匿名函数输入是决策变量向量x输出是标量目标函数值f。x0初始猜测点。这是非线性规划求解的起点非常重要不同的起点可能导致算法收敛到不同的局部最优解。A, b定义线性不等式约束A*x ≤ b。如果没有用空数组[]代替。Aeq, beq定义线性等式约束Aeq*x beq。如果没有用空数组[]代替。lb, ub定义决策变量的下界和上界lb ≤ x ≤ ub。如果没有用空数组[]代替。nonlcon非线性约束函数句柄。一个 MATLAB 函数输入x输出两个向量非线性不等式约束c(x) ≤ 0和非线性等式约束ceq(x) 0。这里必须写成 ≤0 和 0 的形式options优化选项结构体用于设置算法参数如最大迭代次数、函数容差、显示输出等。通常用optimoptions(fmincon, ...)来创建。2.2 算法选择fmincon的四种内功心法fmincon内置了多种算法通过options中的Algorithm选项指定。选择合适的算法对求解效率和成功率至关重要。interior-point内点法默认算法原理通过在可行域内部构造一条路径逼近最优解。它会在目标函数中加入一个“障碍函数”来惩罚接近边界的点随着迭代障碍的影响逐渐减小路径最终收敛到边界上的最优解如果存在的话。优点对于大规模问题变量和约束多通常很有效尤其擅长处理不等式约束。它通常能很好地从可行域内部开始搜索。缺点对于某些非凸问题可能不如基于梯度的算法精确。适用场景大多数中大规模、包含不等式约束的问题的默认选择。sqp序列二次规划法原理在每一步迭代它用二次函数近似目标函数用线性函数近似约束从而构造一个二次规划子问题。求解这个子问题得到搜索方向然后沿此方向进行线搜索。优点通常能提供非常高的精度尤其对于光滑函数。它通常能很好地处理等式约束和活跃的不等式约束。缺点对于大规模问题计算每个子问题的成本可能较高。适用场景中小规模、需要高精度解的问题特别是等式约束较多的问题。active-set有效集法原理它猜测哪些不等式约束在最优解处是“活跃的”即取等号然后只把这些约束当作等式约束来处理在由这些活跃约束定义的“面”上进行搜索。如果猜测错了就调整活跃集。优点对于中小规模问题特别是解位于许多约束边界的情况可能很有效。缺点对于大规模问题维护和更新活跃集的成本很高。适用场景传统算法现在更多被interior-point和sqp取代但在某些特定问题中仍有价值。trust-region-reflective信赖域反射法原理在每一步它在当前点的一个可信赖区域内用较简单的模型如二次模型近似原问题然后在这个小区域内求解子问题。它要求目标函数能提供梯度信息。优点对于大规模无约束或边界约束问题以及某些特殊的线性等式约束问题非常高效。缺点不能处理一般的非线性约束。只能处理边界约束lb, ub和线性等式约束Aeq, beq。适用场景大规模、目标函数可求导、只有边界或线性等式约束的问题。选择建议对于初学者和大多数建模问题直接使用默认的interior-point算法即可。如果求解失败或精度不够可以尝试切换到sqp算法。如果问题只有边界约束可以尝试trust-region-reflective。3. 手把手实现例题的MATLAB代码与逐行解析现在我们回到最初的例题用 MATLAB 和fmincon来求解。我会把代码分成几个部分并详细解释每一行。3.1 第一步定义目标函数我们的原始问题是最大化0.1*x1 0.15*x2。由于fmincon默认最小化我们需要将其转换为最小化-1*(0.1*x1 0.15*x2)。我们可以使用匿名函数来简洁地定义% 定义目标函数注意fmincon求解最小值所以最大化问题要加负号 objective (x) - (0.1*x(1) 0.15*x(2));这里(x)创建了一个匿名函数输入变量是xx(1)对应x1x(2)对应x2。函数体计算-(0.1*x1 0.15*x2)。3.2 第二步定义线性约束我们的线性约束包括一个等式约束和两个不等式约束非负约束。等式约束x1 x2 1。对应Aeq * x beq其中Aeq [1, 1],beq 1。不等式约束非负x1 ≥ 0,x2 ≥ 0。这属于变量的下界约束我们用lb参数处理会更方便。所以这里A和b为空。变量下界lb [0; 0]。变量上界没有明确上界设为空[]。% 线性等式约束 Aeq * x beq Aeq [1, 1]; beq 1; % 线性不等式约束 A * x b 本例中没有用空数组 A []; b []; % 变量的下界和上界 lb [0; 0]; % x1 0, x2 0 ub []; % 无上界3.3 第三步定义非线性约束这是关键的一步。非线性约束函数需要单独写成一个函数文件或嵌套函数。它必须返回两个输出c和ceq分别表示非线性不等式约束c(x) ≤ 0和非线性等式约束ceq(x) 0。我们的风险约束是0.2*x1^2 0.1*x2^2 0.05*x1*x2 ≤ 0.1。 为了符合c(x) ≤ 0的形式我们将其改写为0.2*x1^2 0.1*x2^2 0.05*x1*x2 - 0.1 ≤ 0。 本例中没有非线性等式约束所以ceq返回空数组[]。我们创建一个独立的函数文件nonlinear_constraints.mfunction [c, ceq] nonlinear_constraints(x) % 非线性不等式约束 c(x) 0 % 风险约束: 0.2*x1^2 0.1*x2^2 0.05*x1*x2 0.1 c 0.2*x(1)^2 0.1*x(2)^2 0.05*x(1)*x(2) - 0.1; % 非线性等式约束 ceq(x) 0 本例中没有 ceq []; end或者如果你希望所有代码在一个脚本里可以使用嵌套函数或函数句柄。这里展示函数句柄的写法更简洁% 使用函数句柄定义非线性约束 nonlcon (x) deal(0.2*x(1)^2 0.1*x(2)^2 0.05*x(1)*x(2) - 0.1, []);deal函数将两个输出分别赋给c和ceq。第一个参数是c第二个是ceq。3.4 第四步设置初始点与选项并调用fmincon初始点x0的选择很重要。它应该在可行域内或附近以帮助算法更快收敛。对于本例一个简单的选择是等额投资x0 [0.5; 0.5]。我们也可以设置一些优化选项比如显示迭代过程。% 初始猜测点 (必须在可行域内或附近) x0 [0.5; 0.5]; % 假设初始投资各一半 % 设置优化选项可选但推荐 options optimoptions(fmincon, ... Display, iter, ... % 显示每次迭代的信息 Algorithm, interior-point); % 选择算法默认就是它这里显式指定 % 调用 fmincon 求解 [x_opt, fval_opt, exitflag, output] fmincon(objective, x0, A, b, Aeq, beq, lb, ub, nonlcon, options);3.5 第五步结果解读与后处理求解完成后我们需要解读输出。% 显示最优解和最优值 fprintf(最优投资比例\n); fprintf( 股票A (x1) %.4f\n, x_opt(1)); fprintf( 股票B (x2) %.4f\n, x_opt(2)); fprintf( 总投资比例和 %.4f (应等于1)\n, sum(x_opt)); % 注意fval_opt 是转换后目标函数即 -收益的最小值 % 所以最大收益 -fval_opt max_return -fval_opt; fprintf(\n最大期望收益 %.4f\n, max_return); % 验证约束是否满足 risk 0.2*x_opt(1)^2 0.1*x_opt(2)^2 0.05*x_opt(1)*x_opt(2); fprintf(实际组合风险 %.6f (约束要求 ≤ 0.1)\n, risk); fprintf(风险约束满足情况%s\n, (risk 0.11e-6) ? 满足 : 不满足); % 加一个小容差 % 检查退出标志 fprintf(\n退出标志 exitflag %d\n, exitflag); if exitflag 0 fprintf(优化成功收敛。\n); elseif exitflag 0 fprintf(达到最大迭代次数或函数计算次数。\n); else fprintf(优化未成功。\n); end % 查看输出信息 fprintf(迭代次数%d\n, output.iterations); fprintf(函数计算次数%d\n, output.funcCount);运行这段完整的代码你可能会得到类似如下的输出具体数值可能因算法和版本略有差异Iter Func-count Fval Feasibility Step Length Norm of First-order optimality 0 3 -0.125000 0.000000e00 1.000e00 0.000e00 1.250e-01 1 6 -0.136364 0.000000e00 1.000e00 2.500e-01 1.136e-01 2 9 -0.137931 0.000000e00 1.000e00 6.250e-02 1.034e-01 3 12 -0.138889 0.000000e00 1.000e00 1.562e-02 9.722e-02 ... (更多迭代) 8 27 -0.139241 0.000000e00 1.000e00 9.537e-07 9.537e-07 优化已终止: 一阶最优性度量小于 options.OptimalityTolerance 并且约束违反值小于 options.ConstraintTolerance。 最优投资比例 股票A (x1) 0.4483 股票B (x2) 0.5517 总投资比例和 1.0000 (应等于1) 最大期望收益 0.1392 实际组合风险 0.100000 (约束要求 ≤ 0.1) 风险约束满足情况满足 退出标志 exitflag 1 优化成功收敛。 迭代次数8 函数计算次数27结果分析最优解x1 ≈ 0.4483,x2 ≈ 0.5517。这意味着在风险控制下应将约44.83%的资金投入股票A55.17%投入收益更高的股票B。最大收益0.1392即13.92%介于股票A的10%和股票B的15%之间符合预期。约束验证风险恰好等于0.1说明风险约束是“活跃的”在最优解处取等号它有效地限制了我们对高收益但高风险资产B的过度配置。预算约束也严格满足。算法状态exitflag1表示成功收敛一阶最优性条件得到满足。迭代了8次计算了27次目标/约束函数效率很高。4. 避坑指南与实战经验为什么我的fmincon跑不出结果在实际使用中你可能会遇到fmincon报错、不收敛、或者结果明显不合理的情况。下面我总结了一些最常见的“坑”及其解决办法。4.1 初始点x0的选择一个好的开始是成功的一半fmincon对初始点非常敏感尤其是对于非凸问题有多个局部最优解。问题如果x0离可行域太远或者导致初始约束违反严重算法可能一开始就失败。对策物理意义法根据问题的实际背景给出一个合理的猜测。比如投资问题用均匀分配工程设计问题用经验值。随机初始化多次尝试对于复杂问题可以在可行域内或大致范围内随机生成多个初始点分别运行fmincon然后选择目标函数值最好的那个解。best_x []; best_fval inf; num_trials 20; for i 1:num_trials x0_rand rand(2,1); % 生成[0,1]随机数 x0_rand x0_rand / sum(x0_rand); % 归一化以满足等式约束 [x_temp, fval_temp] fmincon(objective, x0_rand, A, b, Aeq, beq, lb, ub, nonlcon); if fval_temp best_fval best_fval fval_temp; best_x x_temp; end end可行性优先如果找不到明显可行的点可以先用一个简单的优化问题比如最小化约束违反度来找到一个可行的初始点。4.2 非线性约束函数的定义格式错误是万恶之源这是新手最容易出错的地方。坑1不等式和等式的顺序。函数必须返回[c, ceq]c对应c(x) ≤ 0ceq对应ceq(x) 0。顺序反了约束意义就全错了。坑2约束方向。必须写成≤0和0的形式。如果你的约束是g(x) ≥ 0需要转化为-g(x) ≤ 0。坑3向量化输出。即使只有一个非线性约束c和ceq也必须是列向量。例如有两个非线性不等式c1(x)≤0和c2(x)≤0应该写成c [c1(x); c2(x)]。坑4函数未在路径中。如果使用独立的.m文件定义非线性约束函数确保该文件保存在当前工作目录或 MATLAB 搜索路径中。检查清单在调用fmincon前用初始点x0手动调用一次你的非线性约束函数检查输出是否符合预期。[c0, ceq0] nonlcon(x0); disp(初始点处的非线性不等式约束值 c(x0):); disp(c0); disp(初始点处的非线性等式约束值 ceq(x0):); disp(ceq0);确保c0的各分量 ≤ 0或接近0ceq0的各分量 0或接近0。如果c0有正的大数说明初始点严重违反约束。4.3 算法不收敛与参数调优给fmincon一点耐心和指导有时算法会因达到最大迭代次数而停止exitflag0或者根本找不到可行解exitflag-2。对策1调整options。options optimoptions(fmincon, ... Display, iter-detailed, ... % 显示更详细的迭代信息 MaxIterations, 1000, ... % 增加最大迭代次数默认通常是400 MaxFunctionEvaluations, 3000, ... % 增加最大函数计算次数 OptimalityTolerance, 1e-8, ... % 收紧一阶最优性容差要求更高精度 ConstraintTolerance, 1e-8, ... % 收紧约束违反容差 StepTolerance, 1e-10); % 收紧步长容差iter-detailed显示模式能帮你看到每次迭代的目标值、约束违反度和一阶最优性度量是诊断问题的利器。如果看到目标值很久不下降或者约束违反度降不下来可能就是问题本身或初始点有问题。对策2尝试不同算法。如果interior-point失败了切换到sqp试试。options optimoptions(fmincon, Algorithm, sqp, Display, iter);对策3提供解析导数梯度、雅可比矩阵。默认情况下fmincon用有限差分法数值估算导数这既慢又不准。如果你能为目标函数和非线性约束提供解析导数能极大提升速度、精度和稳定性。通过options设置SpecifyObjectiveGradient为true并让目标函数返回两个输出[f, gradf]。设置SpecifyConstraintGradient为true并让非线性约束函数返回四个输出[c, ceq, gradc, gradceq]。 这需要一定的数学推导但对于复杂问题收益巨大。4.4 问题尺度与数值稳定性别让计算机算“糊涂”了如果决策变量的数量级相差巨大比如x1在 1e6 级别x2在 1e-3 级别或者目标函数的数值非常大/小都会导致数值计算困难出现“病态”问题。对策尺度缩放。尽量通过变量代换让所有决策变量都在O(1)的数量级附近。例如如果x1代表以米为单位的长度范围是几千可以令x1_scaled x1 / 1000将单位变为千米。在模型中相应地调整目标函数和约束系数。求解出x_scaled后再变换回原始变量x。4.5 全局最优与局部最优跳出“洼地”fmincon找到的通常是局部最优解。对于非凸问题这可能不是全局最优。对策1多起点优化。如前所述从多个随机初始点出发求解取最好的结果。这是最实用、最常用的方法。对策2使用全局优化算法。对于确实难以找到全局最优的问题可以考虑 MATLAB 的全局优化工具箱中的函数如GlobalSearch或MultiStart。它们会在fmincon的基础上自动进行多起点搜索。problem createOptimProblem(fmincon, objective, objective, x0, x0, ... Aeq, Aeq, beq, beq, lb, lb, ub, ub, nonlcon, nonlcon); gs GlobalSearch; [x_global, fval_global] run(gs, problem);这相当于一个更系统、更智能的“多起点”方法。5. 举一反三非线性规划建模的常见类型与扩展掌握了基本解法后我们来看看非线性规划在数学建模中还有哪些常见形态和扩展思路。5.1 无约束优化fminunc与fminsearch如果问题没有约束就是无约束非线性优化。MATLAB 提供了专门的函数fminunc适用于光滑函数可以使用梯度信息。用法与fmincon类似但不需要约束参数。fminsearch使用 Nelder-Mead 单纯形法不需要导数信息适用于非光滑或求导困难的函数但效率较低变量不宜过多。5.2 最小二乘问题lsqnonlin与lsqcurvefit当目标函数是若干项的平方和时例如曲线拟合问题就构成了非线性最小二乘问题min Σ [F_i(x)]^2。lsqnonlin求解非线性最小二乘问题。lsqcurvefit专门用于数据拟合目标是使模型曲线F(x, xdata)与观测数据ydata的误差平方和最小。这类问题有特殊的结构使用这些专用函数通常比通用的fmincon更高效、更稳定。5.3 多目标优化帕累托最优解集现实中很多问题需要同时优化多个相互冲突的目标如成本最低、性能最好、风险最小。这没有单一的最优解而是一组“帕累托最优解”Pareto optimal solutions即在不使任何一个目标变差的情况下无法再改进其他目标。 MATLAB 的全局优化工具箱提供了gamultiobj函数使用遗传算法来求解多目标优化问题并绘制帕累托前沿Pareto front。5.4 整数/离散变量混合整数非线性规划如果部分决策变量必须是整数如选择设备的台数、路径规划中的节点选择问题就变成了混合整数非线性规划。这类问题难度急剧上升。MATLAB 的全局优化工具箱提供了ga遗传算法和surrogateopt代理优化等函数可以处理这类问题但它们属于启发式算法不能保证找到全局最优解。对于特定结构的问题也可以使用专门的商业求解器如 Gurobi、BARON通过 MATLAB 接口调用。5.5 从例题到实际建模思想迁移我们用的投资组合例题是高度简化的。实际中的马科维茨均值-方差模型目标函数可能是最大化夏普比率收益/风险约束可能包括行业配置上限、不允许做空lb0、最低持仓比例等形成一个更复杂的非线性规划。但核心建模步骤不变定义决策变量各资产的投资权重。建立目标函数最大化收益、最小化风险或最大化风险调整后收益。确定约束条件预算约束、风险约束、法规约束等线性和非线性。选择求解器并实现根据问题规模资产数量和约束类型选用fmincon或其他工具。结果分析与验证检查解是否满足所有约束是否合理并进行敏感性分析例如改变风险上限观察最优投资组合如何变化。这种“定义变量-建立目标与约束-数值求解-分析结果”的框架是解决绝大多数优化类建模问题的通用心法。无论是设计工程结构、调度生产资源还是训练机器学习模型其损失函数最小化本质上也是一个无约束/有约束优化问题其底层逻辑都是相通的。理解并熟练运用fmincon这样的工具就等于掌握了一把打开众多实际问题的钥匙。