学习最优控制时哈密尔顿函数法是我接触最早、也最容易被“表面上的严格”误导的工具。教材里的标准推导通常长这样给出性能指标和状态方程定义哈密尔顿函数写出协态方程和 (\partial H/\partial u0)再配合横截条件然后结束。这套流程确实优雅但放到真实项目中十有八九会遇到“标准假设不满足”的情况控制量有上下限、终端状态不是完全固定而是落在某个流形上、终端时间也是未知的、甚至在最优弧段上哈密尔顿函数对控制变量是线性关系导致无法直接求解。这篇补充想把那些课本里往往一句带过、但在数值求解和工程实践中绕不过去的内容系统讲一遍适合已经学过一遍最优控制、但动手算例时卡壳的读者。1. 标准哈密尔顿函数法的推导边界驻值条件不等于最优1.1 经典五件套的推导路径先快速回顾一下标准问题。性能指标写为[ J \Phi(x(t_f)) \int_{t_0}^{t_f} L(x,u,t),dt ]状态方程[ \dot{x} f(x,u,t) ]定义哈密尔顿函数[ H L \lambda^T f ]一阶必要条件通常包括状态方程(\dot{x} \partial H/\partial \lambda)协态方程(\dot{\lambda} -\partial H/\partial x)最优性条件(\partial H/\partial u 0)横截条件根据终端状态是否固定来决定 (\lambda(t_f))这套推导的数学基础是变分法。你在控制输入上给一个任意小量 (\delta u)它会引起状态变化 (\delta x)再引起性能指标变化 (\delta J)令一阶变分为零就得到上述条件。这里必须强调一个容易忽略的点这套必要条件本质上是“驻值条件”。它就像有限维优化里的 (\nabla f0)只有函数具备凸性等条件时驻点才等价于极小值点。如果哈密尔顿函数对 (u) 不是凸的(\partial H/\partial u0) 可能给出的是鞍点甚至局部极大值。1.2 教材没说的四个实战变体教材例题几乎都集中在“终端固定、控制无约束、终端时间固定、状态无约束”的理想情形。但实际工程项目里我遇到的大多数问题至少破坏其中一条假设。按我的经验需要补充的四类场景可以这么划分变体工程诱因对标准条件的影响终端状态不是完全给定只要求落在目标附近、末端速度不要求精确横截条件出现 (\lambda(t_f)) 的额外自由度控制量有界电机饱和、阀门开度限制、推力上限(\partial H/\partial u0) 失效要用庞特里亚金极小值原理状态路径约束位置不能进入禁区、温度不能超限需要引入额外乘子激活区间和内点边界条件终端时间自由最快到达、最小时间轨迹出现终端时间横截条件哈密尔顿函数守恒变得关键下面四个章节就按照这四个变体展开。每一部分都会先讲理论修正点再给一个可以直接照着算的小例子。2. 横截条件协态变量终点的真正来源2.1 从末端代价和终端约束推导统一形式一旦终端状态不是“完全固定”横截条件就成了整个问题最容易出错的地方。考虑更一般的性能指标带终端代价 (\Phi(x(t_f))) 和终端等式约束 (\psi(x(t_f))0)。为了推导把终端约束乘上拉格朗日乘子 (\mu)增广指标为[ J_a \Phi(x(t_f)) \mu^T\psi(x(t_f)) \int_{t_0}^{t_f} \left[L \lambda^T(f-\dot{x})\right]dt ]对 (x) 做变分经过分部积分后终端处的项是[ -\lambda(t_f)^T\delta x(t_f) \Phi_x\delta x(t_f) \mu^T\psi_x\delta x(t_f) ]其中 (\Phi_x) 和 (\psi_x) 是在 (t_f) 处求值的雅可比矩阵。由于 (x(t_f)) 的变分不是完全自由的它必须满足终端约束的一阶近似 (\psi_x\delta x(t_f)0)。最终得到统一形式的横截条件[ \lambda(t_f) \Phi_x^T \psi_x^T\mu ]这里 (\mu) 是终端约束对应的乘子需要和 TPBVP 一起解出来。如果问题没有终端代价项也没有终端约束且终端状态完全自由那么 (\lambda(t_f)0)。如果 (x(t_f)x_f) 固定则 (\lambda(t_f)) 是自由的由状态约束给足边界条件。不同教材的符号习惯可能把 (\lambda(t_f)) 写成 (-\Phi_x)这个差异来自积分方向或哈密尔顿函数符号定义本质物理含义不变。2.2 协态变量的灵敏度含义把哈密尔顿函数法看成有约束最优化后(\lambda(t)) 其实是“状态方程约束”的影子价格。严谨一点说(\lambda(t_0)) 等于性能指标对初始状态 (x(t_0)) 的灵敏度[ \lambda(t_0) \frac{\partial J^*}{\partial x(t_0)} ]这在数值求解里非常有用。很多时候你猜不准 (\lambda(t_0))是因为没有意识到这个参数对解的敏感程度。终端固定时(\lambda(t_0)) 的一个微小扰动会导致终端状态偏出十万八千里所以直接用单打靶法很容易发散。理解了这层灵敏度含义就会明白为什么工程上推荐使用配点法而不是简单打靶。2.3 算例带终端代价的最小能量控制拿一个最简单的标量系统展示横截条件的用法。系统[ \dot{x}u,\quad x(0)x_0 ]性能指标[ J\frac{1}{2}\int_0^T u^2dt\frac{1}{2}\alpha\left(x(T)-x_d\right)^2 ]目标希望终端尽可能靠近 (x_d)但允许有一点偏差(\alpha) 是终端权重。哈密尔顿函数[ H\frac{1}{2}u^2\lambda u ]协态方程[ \dot{\lambda}-\frac{\partial H}{\partial x}0 ]所以 (\lambda) 是常数。最优性条件[ \frac{\partial H}{\partial u}u\lambda0 \Rightarrow u-\lambda ]状态响应[ x(T)x_0-\lambda T ]横截条件[ \lambda(T)\alpha\left(x(T)-x_d\right)\alpha\left(x_0-\lambda T-x_d\right) ]整理后[ \lambda\frac{\alpha(x_0-x_d)}{1\alpha T} ]于是最优控制为[ u-\frac{\alpha(x_0-x_d)}{1\alpha T} ]注意看如果让 (\alpha\to\infty)也就是终端代价无穷大就退化成了固定终端问题控制量约等于 ((x_d-x_0)/T)这个结果和用固定终端条件直接推完全一致。这个小例子能清楚说明横截条件在“终端状态非固定”时如何参与求解也能看出协态变量和终端代价之间的耦合关系。3. 有控制约束时从驻值条件到极小值条件3.1 为什么教材里的求导条件会被打破实际物理系统的控制量几乎没有不受限的电机有最大电压、火箭有最大推力、阀门有最大开度。当可行控制域是紧集时(\partial H/\partial u0) 不是必要条件因为最优控制可能落在可行域边界上边界处导数不为零。庞特里亚金极小值原理的处理方式是把“驻值”改成“极小值”[ u^(t)\arg\min_{u\in\mathcal{U}} H(x^,\lambda^*,u,t) ]这个条件和哈密尔顿函数法是兼容的只是把无约束最优化替换成了约束最优化。如果控制约束是不等式 (u_{min}\le u\le u_{max})且哈密尔顿函数对 (u) 光滑那么可以先用 KKT 条件判断是内部解还是边界解再分段拼接。3.2 开关函数与 bang-bang 控制很多工程问题中哈密尔顿函数关于控制是线性的。举一个典型结构[ H 1 \lambda^T(AxBu) ]此时求导得到[ \frac{\partial H}{\partial u}B^T\lambda ]这个表达式里根本没有 (u)无法通过导数为零解出控制量。最优控制完全由开关函数 (s(t)B^T\lambda(t)) 的符号决定[ u^*-U_{max},\text{sign}(s(t)) ]这就是 bang-bang 控制控制量永远在全开和全关之间切换。系统状态和协态共同决定切换时刻。工程上更关心的是开关次数和切换条件。对于线性时不变系统协态方程是线性的开关函数通常是若干指数项的组合在有限时间内最多只有有限次切换。实现时要注意开关函数接近零的时刻不能靠简单比较符号判断数值积分步长不够小时容易漏掉一次切换导致后续轨迹完全错误。3.3 状态路径约束的引入如果问题还带状态路径约束比如[ g(x,t)\le 0 ]那么需要在哈密尔顿函数中引入额外的非负乘子 (v(t))[ \tilde{H}L\lambda^T fv g ]这个乘子只在 (g0) 的激活区间非零在未激活区间为零并且满足互补松弛条件[ v(t)g(x,t)0 ]棘手的地方在于路径约束激活区间的切换点位置一开始是未知的。在每个未激活区间系统按正常 TPBVP 积分在激活区间系统被约束在 (g0) 的流形上运动控制量由约束的二阶时间导数决定。两个区间之间的状态和协态需要满足内点边界条件这些边界条件本质上是哈密尔顿函数连续性和协态连续性的一部分。我自己做轨迹优化时处理路径约束最稳妥的方式是先不直接上精确约束而是用一个大惩罚项逼近路径约束等求解器收敛到合理区域再逐步增大惩罚系数把解逼到真正的约束边界。直接让算法从含路径约束的 TPBVP 开始迭代极少能一次收敛。4. 终端时间自由哈密尔顿函数多出来一个守恒量4.1 终端时间变分带来的额外必要条件终端时间 (t_f) 也可以作为优化变量。例如“尽快到达目标”“最小时间轨迹”这类问题。此时性能指标对终端时间的导数也必须为零。推导后得到的终端时间横截条件为[ H(t_f)\frac{\partial \Phi}{\partial t_f}0 ]如果性能指标里没有显式终端时间项即 (\partial \Phi/\partial t_f0)那么终端时刻哈密尔顿函数值必须为零。如果系统是自治的也就是 (f) 和 (L) 都不显含时间那么哈密尔顿函数沿最优轨迹满足[ \frac{dH}{dt}\frac{\partial H}{\partial t}0 ]所以 (H) 是常数。结合终端时间横截条件这个常数通常被固定为某个确定值。这个性质极其实用它是检验数值解是否正确的最强工具之一。4.2 算例最短时间双积分器考虑一个最简单的双积分器[ \dot{x}_1x_2,\quad \dot{x}_2u,\quad |u|\le1 ]目标是从初始状态转移到原点并且终端时间 (t_f) 最小。哈密尔顿函数[ H1\lambda_1 x_2\lambda_2 u ]因为是对 (u) 线性最优控制为 bang-bang 形式[ u-\text{sign}(\lambda_2) ]协态方程[ \dot{\lambda}_10,\quad \dot{\lambda}_2-\lambda_1 ]所以 (\lambda_2) 随时间线性变化每个控制弧段上最多只有一次符号变化。整个最优轨迹最多经历一次切换。在相平面里切换曲线是[ x_1-\frac{1}{2}x_2|x_2| ]这条曲线就是著名的 bang-bang 开关线。如果初始状态在这条曲线上直接用一个控制量到达原点不在曲线上先用一个控制量切换到曲线再沿曲线到达原点。这个例子是解析解和数值解互相验证的绝佳对象因为你能清楚地看到哈密尔顿函数常数、终端时间横截条件、开关次数这些抽象概念分别对应什么。4.3 用哈密尔顿函数常数校验数值解在处理终端时间自由问题时我养成的习惯是常把时间变量做无量纲化处理。令 (\taut/t_f)把 (t_f) 作为未知参数加入边界条件中。这样做的好处是积分区间固定为 ([0,1])适合标准求解器。但是引入了额外边界条件初值猜测难度也增加了。有一个非常可靠的校验办法每隔一段时间检查一次哈密尔顿函数值看它是否在整条轨迹上保持常数。如果状态和协态都收敛得很漂亮但哈密尔顿函数明显不是常数基本可以断定解有问题。我在解带多个切换点的轨迹优化问题时遇到过好几次这种“看起来很好、实际上错误”的情况最终都是靠哈密尔顿函数守恒找出来的。这个技巧不花额外成本强烈建议做数值求解前先在代码里写好这个检查项。5. 奇异弧段哈密尔顿函数法最容易被忽略的场景5.1 奇异最优控制的出现原因当哈密尔顿函数关于控制变量是线性时最优性条件 (\partial H/\partial u0) 退化成了一个纯函数等式[ s(t)\frac{\partial H}{\partial u}0 ]这个方程里没有 (u)所以无法直接解出控制量。如果 (s(t)) 在某个区间上恒等于零而不是只在孤立时刻过零就进入了奇异弧段。这种情况在教材例题中很少出现因为教学例题都喜欢用 (\frac{1}{2}u^2) 这样的二次型控制代价让最优性条件对 (u) 可解。但实际工程中最小燃料、某些能量管理、交通流模型里控制变量往往是以线性方式进入状态方程的奇异弧段并不罕见甚至是最优解中非常关键的一段。5.2 奇异弧段的求解套路奇异弧段虽然不能直接解出 (u)但可以借助“奇异弧上控制变量的系数函数必须保持为零”这一事实对时间逐阶求导。假设[ s\frac{\partial H}{\partial u} ]在奇异弧上有[ s0,\quad \frac{ds}{dt}0,\quad \frac{d^2s}{dt^2}0,\dots ]一直求导到 (u) 显式出现在某个表达式中。设最低阶数为 (q)那么奇异弧通常可以用[ s,\dot{s},\dots,s^{(2q-1)} ]这些量共同描述。用 (s^{(2q)}0) 可以解出奇异控制量 (u_s)但它只是必要条件还需要检查广义勒让德-克莱伯条件来判断这个弧段是不是局部最优[ (-1)^q\frac{\partial}{\partial u}\left(\frac{d^{2q}}{dt^{2q}}s\right)\ge0 ]这个条件说起来复杂实际应用时可以通过引入一个很小的正则项来绕开。例如把性能指标改成[ J_\epsilon J \frac{\epsilon}{2}\int u^2dt ]当 (\epsilon) 很小时问题退化为非奇异控制求解器可以正常工作。求出解以后把 (\epsilon) 逐渐减小再求解外推得到奇异最优解的近似值。我在很多机械臂和移动小车问题中使用过这个办法比直接推导奇异弧条件省事得多尤其适合工程节拍比较快的场景。5.3 奇异弧段在数值求解中为什么特别难奇异弧段的本质问题是沿着弧段哈密尔顿函数的二阶变化率中等可能退化导致 TPBVP 的雅可比矩阵接近奇异。直接打靶法基本会失败因为初值 (s) 的微小扰动会让切换函数在奇异弧附近反复振荡积分过程极不稳定。配点法通常能稍微好一点但也需要给定一个足够好的初始猜测。我的经验是除非问题非常简单否则不推荐在奇异弧上直接手推解析公式而是用正则化手段把奇异弧变成非奇异问题求解最后再看正则项趋于零时的趋势。这种“先软化再逼近”的策略在很多非线性最优控制问题上都适用不只是奇异问题。6. 从必要条件到实际算例打靶法与避坑经验6.1 把哈密尔顿函数法变成可求解的两点边值问题无论标准情形还是各种补充变体哈密尔顿函数法最后都会把最优控制问题化成如下形式[ \dot{x}H_\lambda,\quad \dot{\lambda}-H_x ]再加上边界条件和横截条件构成一个两点边值问题。最简单直观的解法是单打靶法。思路是把未知的协态初值 (\lambda(t_0)s) 当作待求参数从 (t_0) 正向积分到 (t_f)然后看终端残差是否满足[ r(s)\lambda(t_f)-\lambda_{\text{target}}(x(t_f))0 ]如果 (r(s)) 不为零就调整 (s) 重新积分。调整方法可以用二分法或牛顿法。出于教学目的我经常给学生看一个非常小的 Python 框架import numpy as np from scipy.integrate import solve_ivp from scipy.optimize import fsolve def dynamics(t, y, params): x1, x2, lam1, lam2 y u -np.sign(lam2) # bang-bang控制 dx1 x2 dx2 u dlam1 0.0 dlam2 -lam1 return [dx1, dx2, dlam1, dlam2] def terminal_residual(s): y0 [x1_0, x2_0, s[0], s[1]] sol solve_ivp(dynamics, [0, t_f], y0, t_eval[t_f]) x1_f, x2_f, lam1_f, lam2_f sol.y[:, -1] return [x1_f, x2_f] s_guess [1.0, 1.0] s_sol fsolve(terminal_residual, s_guess)这里面的fsolve天然是局部算法初始猜测差太远就会失败。所以我不建议在没有把握的情况下直接上 fsolve更推荐先把问题无量纲化再用解析解或一个简化模型给出初值。6.2 直接调库还是手写积分手写单打靶法虽然直观但非线性较强的多输入系统几乎必炸。实际工程里更推荐直接使用成熟的边界值问题求解器比如 MATLAB 的bvp4c或 Python 的scipy.integrate.solve_bvp。它们内部用的配点法会把整个时间域离散成网格同时求解状态、协态和未知参数比单纯正向积分加牛顿迭代要稳定得多。不过我在这里提醒一句配点法对网格密度和初始猜测也有要求尤其是 bang-bang 控制这类非光滑问题切换点附近必须把网格放密否则控制切换会被抹平解的精度差一个量级。能用解析式表达控制律的问题尽量把控制表达式写进方程里而不是交给求解器去离散。6.3 数值调试中最常踩的五个坑结合我自己的项目经验以下五个问题几乎覆盖了最优控制数值求解的绝大多数“翻车”现场横截条件写错符号。这是出现频率最高的问题。建议拿到一个新问题先找一个有解析解的极限情形核对符号再开始调参数。控制约束的饱和逻辑写错。特别是多个输入同时饱和时判断顺序不对会让积分结果在某个时刻突跳。开关函数接近零时没做精细处理。工程上经常遇到切换函数在很长一段时间内都接近零看似奇异弧其实不是盲目按奇异处理会让弧段结构出错。量纲不一致导致数值病态。状态量、控制量、时间尺度如果差几个数量级求解器收敛性会急剧下降。初始猜测全是零。如果从零猜很多最优控制问题的协态积分会直接退化到平凡解最后得到一条不满足约束的假轨迹。调试时我的固定动作是先检查哈密尔顿函数是否满足应有的守恒特性再检查终端横截条件残差最后检查控制序列在每个时间点上是否真的让哈密尔顿函数达到最小。这三步走完绝大多数问题都能定位到具体原因。我个人养成的一个习惯是拿到一个最优控制问题后不急着上数值工具先把哈密尔顿函数写出来然后逐条对照标准假设——终端状态是固定还是自由、控制量是否有界、是否存在路径约束、终端时间是否固定、哈密尔顿函数对控制是否严格凸。把这些关键点确认完再决定使用哪套补充条件。很多问题之所以解不出来并不是缺少软件或计算量而是一开始最优性条件就没有写完整横截条件漏了一项、奇异弧段没有识别、哈密尔顿函数常数没有利用这些小问题几乎是我见过最多的翻车原因。希望这篇围绕哈密尔顿函数法补充展开的内容能帮你把这些坑提前避开。