第四章 最优制导律的数值求解框架:两点边值问题、打靶法、多重打靶法
好,咱们进入正题。前面几章我们把最优控制的理论基础打了一遍,从变分法到庞特里亚金极值原理,再到各种性能指标下的解析解。说实话,解析解看着漂亮,但实际工程中能用上的情况真不多。为什么?因为大多数真实场景下的制导问题,边界条件复杂、约束多、模型非线性,解析解根本写不出来。
那怎么办?数值求解呗。这一章我就带你看看,怎么用计算机把最优制导律“算”出来。核心就是解决一个叫两点边值问题的东西,简称TPBVP。
4.1 两点边值问题(TPBVP)——最优制导的数学骨架
先说说什么是两点边值问题。你想想看,最优制导律的求解,本质上是在找一条轨迹,从初始状态飞到终端状态,同时让某个性能指标最小。这个过程由状态方程和协态方程共同描述,状态有初始条件,协态有终端条件。
这就麻烦了。状态方程从初始时刻往前推,协态方程从终端时刻往后推。一个向前,一个向后,在中间某个点碰头。这就是典型的两点边值问题。
TPBVP的标准形式:
状态方程: dx/dt = f(x, u, t), x(t0) = x0
协态方程: dλ/dt = -∂H/∂x, λ(tf) = ∂φ/∂x(tf)
控制方程: ∂H/∂u = 0 (最优性条件)
我在做某个拦截弹的制导律设计时,遇到过这样一个情况:状态变量有6个,协态变量又有6个,总共12个微分方程要同时求解。边界条件一半在起点,一半在终点。当时我盯着这堆方程看了半天,心想这玩意儿怎么解?
嗯,后来我明白了,解决TPBVP的核心思路就一个——把边值问题转化成初值问题。怎么转化?打靶法就是干这个的。
4.2 打靶法——最朴素的思路
打靶法这个名字挺形象的。你想想看,打靶的时候,你瞄准一个方向开一枪,子弹飞出去,看落点偏了没。偏了就调整角度,再打一枪,直到命中靶心。
打靶法求解TPBVP也是这个道理:
- 猜测一组缺失的初始条件(比如协态变量的初值)
- 正向积分整个微分方程组到终端时刻
- 检查终端条件是否满足
- 调整猜测值,重复直到收敛
说白了,就是不断试错。但试错也得有方法,不能瞎猜。常用的方法是牛顿法或者拟牛顿法来更新猜测值。
我个人习惯:打靶法的初始猜测很关键。我一般先用解析解或者工程经验给一个大概的值,别从零开始猜。否则收敛慢不说,还可能发散。
下面给一个简单的打靶法代码框架,用Python写的:
import numpy as np
from scipy.integrate import solve_ivp
from scipy.optimize import fsolve
def shooting_method(x0, tf, guess_lambda0, system_dynamics, terminal_condition):
"""
打靶法求解TPBVP
x0: 初始状态
tf: 终端时间
guess_lambda0: 协态初值猜测
"""
def f_to_solve(lambda0):
# 将状态和协态拼在一起积分
def augmented_dynamics(t, y):
x = y[:n]
lam = y[n:]
# 计算最优控制
u = compute_optimal_control(x, lam)
# 返回状态和协态的导数
dx = state_dynamics(x, u)
dlam = costate_dynamics(x, lam, u)
return np.concatenate([dx, dlam])
# 初始条件
y0 = np.concatenate([x0, lambda0])
# 正向积分
sol = solve_ivp(augmented_dynamics, [0, tf], y0, method='RK45')
# 提取终端状态
xf = sol.y[:n, -1]
lamf = sol.y[n:, -1]
# 返回终端条件误差
return terminal_condition(xf, lamf)
# 用fsolve求解
lambda0_opt = fsolve(f_to_solve, guess_lambda0)
return lambda0_opt
这段代码看着简单,但实际用起来坑不少。我曾经在某个项目里用打靶法求解一个高超声速飞行器的再入制导问题,结果死活不收敛。后来发现是终端条件太敏感,稍微改一点猜测值,终端状态就飞了。
注意:打靶法对初始猜测很敏感。如果系统是非线性的,或者积分时间很长,一个小小的误差就会被放大,导致终端条件严重偏离。这就是所谓的“敏感性”问题。
4.3 多重打靶法——把长距离拆成短距离
既然打靶法怕长积分区间,那我把一段长路径拆成多段短路径,不就行了?这就是多重打靶法的核心思想。
多重打靶法把整个时间区间分成N段,每段都当成一个独立的打靶问题。段与段之间通过连续性条件连接起来——上一段的终点状态等于下一段的起点状态。
这样做的好处很明显:
- 每段积分时间短,数值稳定性好
- 对初始猜测的敏感度降低
- 可以并行计算各段
当然,代价也大——未知变量数量暴增。原来只需要猜一组协态初值,现在每段都要猜状态和协态的初值。求解规模从n维变成了N×2n维。
下面我用SVG画一张图,帮你理解多重打靶法的结构:
你看这张图,时间轴被分成了5段,每段独立积分。段与段之间用连续性条件约束。这样即使某一段的积分误差较大,也不会扩散到整个区间。
4.4 两种方法的对比与选择
到底用打靶法还是多重打靶法?我个人的经验是:
| 对比项 | 打靶法 | 多重打靶法 |
|---|---|---|
| 变量数量 | 少(n个) | 多(N×2n个) |
| 数值稳定性 | 差(长积分区间) | 好(短积分区间) |
| 对初值敏感度 | 高 | 低 |
| 实现难度 | 简单 | 中等 |
| 适用场景 | 短时间、弱非线性 | 长时间、强非线性 |
如果你只是做个简单的制导律验证,积分时间短,打靶法就够了。但如果你要处理的是像再入飞行器那样的长航时、强非线性问题,我建议你直接上多重打靶法。别问我怎么知道的——我曾经在打靶法上浪费了两周时间,换成多重打靶法两天就搞定了。
避坑指南:多重打靶法的段数怎么选?我一般按这个原则:每段积分时间不超过系统最小时间常数的5倍。这样既能保证稳定性,又不会让变量数量爆炸。
4.5 数值求解的工程要点
最后,我总结几个实际编程中容易踩的坑:
- 雅可比矩阵一定要给:用fsolve或者牛顿法时,如果能提供解析的雅可比矩阵,收敛速度会快很多。数值差分算雅可比虽然省事,但精度差,容易导致不收敛。
- 缩放很重要:状态变量和协态变量的量级可能差好几个数量级。我习惯把所有变量归一化到[0,1]或者[-1,1]区间,这样求解器工作起来舒服得多。
- 积分器选型:对于刚性问题,用隐式积分器(比如Radau方法);对于非刚性问题,RK45就够了。别问我怎么知道的——有一次我用RK45解一个刚性方程,结果步长小到几乎算不动。
- 终端条件要松弛:有时候终端条件太严格,数值上根本达不到。我一般会加一个小的容忍度,比如终端位置误差在1米以内就算满足。
好了,这一章的内容就到这儿。数值求解TPBVP是制导律设计的核心技能,打靶法和多重打靶法是最常用的两种手段。你想想看,掌握了这些,以后遇到再复杂的制导问题,至少知道从哪儿下手了。
公众号:蓝海资料掘金营,微信deep3321