企业官网建设流程全解析

航天器末端追逃博弈仿真，本质上就是在问一个问题：当攻防双方都在用脑子计算怎么对付对方时，最后会演化成什么状态？追方想脱靶量尽可能小，逃方想脱靶量尽可能大，两边都有自己的控制能力和燃料约束。这不是一个“单方面最优”问题，而是双方都在动态博弈。复现这类仿真的难点，并不在于把相对运动方程写出来，而在于如何把一个“两边都不吃亏”的策略组找出来。严格意义上的纳什均衡在连续非线性系统中几乎不可能解析求解，所以工程实践里退而求其次，找ε-纳什均衡。下面我完整复盘一下这个项目的建模、求解和Matlab复现过程。

这种仿真的核心价值在于把博弈论中的均衡概念工程化。你不需要追求“绝对最优”，只需要验证“任何一方单方面改变策略，额外收益不会超过ε”。这个思路非常贴近真实工程：信息不全、模型不精确、对手不按常理出牌，追求精确解没有实际意义。这篇文章适合正在做航天器相对运动、微分博弈、或者想入门“博弈+最优控制”仿真的人参考，代码结构可以直接套到类似问题上。

1. 项目整体设计与为什么选择ε-纳什均衡

1.1 末端追逃博弈的工程背景

末端追逃不是科幻片里的专属场景。交会对接的最后逼近段、在轨服务中的非合作目标接近、空间碎片清除的末制导环节，都可以抽象成追逃博弈。以在轨服务为例：服务航天器要接近目标星，但目标星如果具备机动能力，它的最优策略恰恰是“别被抓住”。于是服务星不能只按照预先规划好的轨迹飞，它必须实时根据目标星的机动调整自己的逼近策略，同时目标星也在实时观测跟踪星的机动意图。

这种双向耦合的动态决策问题，用单边最优控制解不了。单边最优控制的逻辑是“目标按已知或可预测的方式运动，我方最优”，但博弈里对手的行为是策略性的，它会利用你的预测误差。所以必须把双方的决策放在同一个框架里同时考虑，这就是微分博弈的基本思想。

1.2 精确均衡在工程里的困境

微分博弈的精确解通常对应Hamilton-Jacobi-Isaacs(HJI)方程或两点边值问题。状态维度低、模型线性时还能解，可一旦引入非线性相对运动、控制饱和、状态约束，HJI方程基本没法解析求解，数值求解也要面对“维数灾难”和收敛性难题。

ε-纳什均衡做了一件务实的事：把“最优”的标准放宽。数学定义也不复杂，假设追方策略为u_p，逃方策略为u_e，如果对任意可行策略u_p'和u_e'都有：

• J_p(u_p, u_e) ≤ J_p(u_p', u_e) + ε
• J_e(u_p, u_e) ≤ J_e(u_p, u_e') + ε

那就称策略组(u_p, u_e)是ε-纳什均衡。第一个式子说追方单方面换策略，最多能让自己的代价再降低ε；第二个式子说逃方单方面换策略，最多能让自己的代价再降低ε。当双方都找不到能显著改善自己处境的反制策略时，这个局面在工程上就是稳定的。

打个比方，就像两个人下棋，谁都不能保证每一步都是绝对最优，但只要任何一方换一套走法，胜率提升不超过一个很小的百分比，双方就会认可当前局面“够均衡了”。ε就是这个“可接受的提升空间”。

1.3 整体技术路线

整个复现分成三层：

1. 模型层：建立追方与逃方的相对运动方程，定义双方性能指标和控制约束。
2. 求解层：用SDRE（状态依赖Riccati方程）在每个时间步近似求解局部最优反馈策略。
3. 验证层：对求得的策略组做策略偏离扫描，计算最大单方面收益改善量，判断是否满足预设ε阈值。

选择Matlab的原因很直接：矩阵运算、控制相关函数（care、ode45）齐全，画图方便，特别适合这种“状态维度不高但逻辑嵌套较深”的仿真项目。

2. 系统建模与博弈问题数学表达

2.1 相对运动状态方程

末端追逃一般建立在相对运动坐标系下。为了兼顾可读性和实用性，我以二维平面内的相对运动为例，状态向量取为：

x = [Δr_x, Δr_y, Δv_x, Δv_y]^T

前两个分量是追方相对逃方的位置偏差，后两个是相对速度偏差。如果参考轨道是圆轨道，经典模型是Hill方程，也就是CW方程。追方和逃方的控制输入分别记为u_p和u_e，作用在加速度层。但纯粹CW方程是线性时不变系统，直接上LQR就够了，体现不出SDRE的必要性。所以我在状态矩阵里加了一项与相对距离相关的非线性修正，模拟真实末端接近时重力差项和视线角变化的耦合影响。

写出来大概是：

x' = A(x)x + B_p u_p + B_e u_e

其中A(x)中包含线性部分A_lin和与当前相对位置相关的非线性项ΔA(x)。具体到代码里，就是根据当前状态实时组装矩阵，这也是“状态依赖矩阵”这个名字的来源。

2.2 双方性能指标与控制约束

追方的目标是让终端相对距离尽可能小，同时控制能量不要太大。性能指标为：

J_p = ∫(x^T Q_p x + u_p^T R_p u_p) dt + x(T)^T Q_fp x(T)

逃方正好反过来，要让终端相对距离尽量大。但为了数值上可解，我把逃方的终端目标设置为一个距离追方较远的参考状态，而不是直接对终端距离取负。逃方指标形式类似：

J_e = ∫(u_e^T R_e u_e) dt + (x(T) - x_ref)^T Q_fe (x(T) - x_ref)

注意这里的x_ref不一定需要是固定点，在末端博弈里可以理解成“逃方希望到达的相对位置区域”。

控制约束方面，双方都有加速度上限。工程里直接用饱和函数处理，代码里用saturate函数把指令加速度限制在各自上限内。这个细节非常关键，因为无约束最优控制的解往往远超物理可实现范围，加饱和后曲线才真实。

2.3 反馈策略的表示

无论追方还是逃方，我都把策略限定为线性状态反馈：

u_p = -K_p x u_e = -K_e x

K矩阵由SDRE在每个时间步实时更新。这样做的工程理由很充分：星载计算机需要简单、确定、可预计算的策略，非线性最优反馈难以在线实现，而线性反馈配合状态依赖增益可以在每个控制周期内快速计算出来。

3. SDRE求解原理与ε均衡验证流程

3.1 为什么用SDRE而不是直接求解HJI

微分博弈的精确解需要求解HJI偏微分方程，状态维度一上去就非常困难。SDRE的思路是绕开全局偏微分方程，把非线性系统在每一个采样时刻“伪线性化”，然后求解一个标准的代数Riccati方程。

具体来说，对非线性系统x' = f(x) + g(x)u，如果能写成x' = A(x)x + B(x)u的形式，那么在当前状态x_k处，就把A(x_k)当作线性系统的系统矩阵，求解对应的Riccati方程得到P，然后令K = R^{-1}B^T P。下一时刻状态变了，再重新组装A(x)再解一次。这就是“状态依赖”的含义。

SDRE不是精确最优解，但它有几个非常实用的优点：实现简单、每步计算量小、能够跟随非线性特性的变化，只要A(x)的组装方式合理，闭环稳定性在大多数工程场景下都有保障。

3.2 单步Riccati方程求解

对追方，在每个时间步求解：

A^T P_p + P_p A - P_p B R_p^{-1} B^T P_p + Q_p = 0

对逃方同理，只是Q和R换成自己的权重。Matlab里的care函数可以直接返回代数Riccati方程的解，这也是用Matlab复现这类问题最顺手的地方。

这里有个容易踩坑的点：care函数求解时需要保证能控性和权重矩阵的正定性。如果设置的R矩阵有问题，或者状态矩阵A在某个状态附近出现退化，care会直接报错。我个人建议在每次调用care前加一个判断，检测B的能控性矩阵的秩是否满足要求，不满足时跳过该步并沿用上一步的增益。这个保护机制在仿真里能省掉很多头疼事。

3.3 策略迭代与ε均衡校验流程

SDRE解出的策略组只是“双方各自认为自己局部最优”的策略，并不天然满足ε均衡的条件。所以必须加一步策略验证，流程如下：

1. 用SDRE求出初始策略组(u_p^0, u_e^0)。
2. 固定u_e^0，对u_p做扰动扫描，生成若干候选策略u_p^i，逐一计算J_p(u_p^i, u_e^0)，找出相对J_p(u_p^0, u_e^0)的最大改善量ΔJ_p。
3. 固定u_p^0，对u_e做同样的扰动扫描，找出最大改善量ΔJ_e。
4. 取ε_est = max(ΔJ_p, ΔJ_e)，与预设阈值ε_th比较。
5. 如果ε_est > ε_th，调整权重矩阵重新求解SDRE，重复整个流程。
6. 如果ε_est ≤ ε_th，认为当前策略组达到ε-纳什均衡，输出并绘制结果。

策略扰动扫描不用太复杂，通常把增益矩阵乘以0.8、1.2、1.5等系数，再叠加一个小幅偏置项，就能覆盖主要偏离方向。真正做项目时，可以在扫描集合里加入几个有代表性的机动模式，比如逃方的“横向脉冲规避”“连续正弦机动”“末端急转弯”，这比纯粹随机扰动更有工程意义。

4. Matlab代码实现与参数配置

4.1 项目文件结构

我习惯把仿真拆成几个独立文件，方便调试：

• main.m：主程序，负责参数初始化、主循环、结果汇总
• computeA.m：根据当前状态组装状态依赖矩阵
• solve_riccati.m：求解追方和逃方的Riccati方程
• dynamics.m：相对运动动力学函数
• rk4_step.m：四阶Runge-Kutta积分步
• saturate.m：控制量饱和处理
• plot_results.m：绘图

拆开写的好处是，想单独测试某一环时不用把整个主程序跑一遍。尤其是computeA这种核心函数，单独写一个测试脚本验证矩阵组装是否正确，比在完整仿真里肉眼找bug要快得多。

4.2 主循环与积分器实现

主循环的核心结构如下：

% 参数初始化 dt = 0.05; T = 30; t = 0:dt:T; N = length(t); % 初始相对状态 [dx, dy, dvx, dvy] x0 = [50; 60; -5; 2]; % 权重矩阵 Q_p = diag([10, 10, 0.5, 0.5]); R_p = 2 * eye(2); Q_e = diag([5, 5, 0.2, 0.2]); R_e = 3 * eye(2); Q_fp = diag([30, 30, 1, 1]); Q_fe = diag([20, 20, 1, 1]); % 控制饱和上限 u_p_max = 5; u_e_max = 3; x = x0; x_history = zeros(4, N); x_history(:, 1) = x; for k = 1:N-1 A = computeA(x); [K_p, K_e] = solve_riccati(A, Q_p, R_p, Q_e, R_e, B); u_p = -K_p * x; u_e = -K_e * x; u_p = saturate(u_p, u_p_max); u_e = saturate(u_e, u_e_max); x = rk4_step(@dynamics, x, u_p, u_e, dt); x_history(:, k+1) = x; end

四阶Runge-Kutta积分器是自己写的，只有十几行：

function x_next = rk4_step(f, x, u_p, u_e, dt) k1 = f(x, u_p, u_e); k2 = f(x + 0.5*dt*k1, u_p, u_e); k3 = f(x + 0.5*dt*k2, u_p, u_e); k4 = f(x + dt*k3, u_p, u_e); x_next = x + (dt/6) * (k1 + 2*k2 + 2*k3 + k4); end

为什么不用ode45？因为每个时间步的增益矩阵K都要根据当前状态重新计算，ode45内部的每个函数求值点都需要重新计算K，误差控制逻辑反而会拖慢计算并增加复杂度。固定步长RK4在这个场景下足够精确，dt取0.05秒时，30秒仿真只要600步，计算量完全不是问题。

4.3 状态依赖矩阵的计算

computeA函数是最容易写错的地方。我踩过的坑是：A(x)线性项和非线性项加错位置，导致仿真在某个边界状态突然跳变。推荐做法是先算线性CW部分，再叠加修正项，并加一个随距离渐变的加权因子：

function A = computeA(x) dr_x = x(1); dr_y = x(2); r = sqrt(dr_x^2 + dr_y^2); % 参考轨道角速度，近地圆轨道约0.001 rad/s omega = 0.001; % CW方程线性部分 A_lin = [0, 0, 1, 0; 0, 0, 0, 1; 3*omega^2, 0, 0, 2*omega; 0, 0, -2*omega, 0]; % 非线性修正项，r越小影响越明显 mu = 3.986e14; R0 = 6771000; delta = mu / R0^3; r_ref = max(r, 10); A_nonlin = [0, 0, 0, 0; 0, 0, 0, 0; -delta * (R0/r_ref)^3, 0, 0, 0; 0, -delta * (R0/r_ref)^3, 0, 0]; A = A_lin + A_nonlin; end

这里r_ref下限取10米，是为了避免末端距离趋于零时产生数值奇异。这个下限值对结果有影响，建议根据仿真尺度和控制精度要求调整。

4.4 参数配置的工程心得

Q和R矩阵的选取在这个项目里既是科学也是手感。我调试时的经验是：

• Q_p相对Q_e越大，追方就越“激进”，会在早期消耗更多燃料压着逃方打。
• R_p越小，追方的控制代价越便宜，末端速度会变得很高。
• R_e越大，逃方越畏手畏脚，容易被追方逼入死角。
• 两组权重不能悬殊过大，否则SDRE解出的策略会触发饱和限幅，仿真曲线出现高频抖动。

我给出一组能跑通的参数供参考，但强烈建议自己调一遍。调参的直观方法是：先固定初始状态，只改一个权重参数，对比相对距离曲线和控制量曲线的变化，理解每个权重的作用后再组合调整。这个项目的“手感”非常重要，权重选得合理，后期ε校验会容易通过；选得不合理，策略偏离扫描的结果会很难看。

5. 仿真结果解析与均衡定量验证

5.1 典型场景与初始状态

仿真场景设置如下：追方初始相对位置在(50, 60)米，相对速度为(-5, 2)米/秒，即追方正在向逃方方向靠近。仿真时长30秒，双方控制加速度上限分别为5米/秒²和3米/秒²，这个不对称模拟的是追方机动能力强于逃方的情形。

5.2 相对轨迹与脱靶量分析

运行代码后，最直观的结果是相对距离随时间的变化曲线。如果模型和权重设置合理，曲线不会单调下降，而是会出现若干段“僵持”区域：追方逼近，逃方横向机动拉开距离，追方再修正方向继续逼近。这个过程就体现了博弈的动态耦合，也是这类仿真最有趣的地方。

我复现的典型结果是终端相对距离约在8到15米之间。追方尽管机动能力占优，但逃方通过横向机动有效避免了脱靶量被压到零。如果把追方权重调到特别大，终端距离能压到3米以内，但代价是追方控制量全程贴近饱和，且ε校验指标明显变差——因为逃方只要稍微改变策略就能大幅提升自身收益，策略组离均衡很远。

5.3 控制量曲线特征

追方的控制量曲线通常呈“阶梯状”，初期接近时有一个明显的持续推力段，中途遇到逃方机动产生方向修正，末端再次加大制动。逃方控制量则更多是“脉冲式”，因为逃方不需要持续大推力对抗，只需要在关键时机做横向机动，打破追方的预测，让自己脱离追杀窗口。

如果控制量曲线出现高频来回振荡，大概率是权重设置过小或者SDRE增益在某个状态附近发生跳变。遇到这种情况我会先检查A(x)的组装在对应状态附近是否连续，其次再考虑增大R矩阵的惩罚力度。

5.4 ε均衡的定量验证结果

我用增益扰动扫描方法对SDRE策略组做验证，结果大致如下：

当ε阈值设为5%时，大部分策略偏离都在可接受范围内，但逃方如果改用正弦横向机动，其收益改善达到3.7%，虽然没超过5%，但已经说明当前策略组存在被针对性突破的可能。如果把阈值收紧到3%，就需要重新调整权重并迭代求解。

这恰好体现了ε均衡和精确均衡的差异：精确均衡要求任何偏离都无法改善收益，工程上基本做不到；而ε均衡给出一个可以接受的范围，只要所有合理偏离的收益改善都不超过阈值，策略组就可以投入使用。

6. 复现避坑与调参经验实录

6.1 仿真发散问题

最常遇到的发散原因是A(x)矩阵在r接近零时出现大数值项。我在computeA里用r_ref=max(r,10)做了下限保护，这个处理非常关键。另一个常见原因是dt取得太大，RK4的稳定性不够，SDRE增益更新跟不上系统状态快速变化。建议把dt从0.01逐步增加到0.1试一遍，观察曲线是否稳定。

如果仿真中后期突然发散，优先检查控制量是否触发饱和。饱和函数的处理会引入非线性，如果单纯在代码里用min/max截断，可能在切换点产生增益突变。推荐用平滑饱和函数过渡，或者至少在主循环里增加一个控制量变化率限制。

6.2 care函数报错排查

care报错一般只有两类原因：一是R矩阵不是正定矩阵，二是能控性条件不满足。检查方法很简单，用eig(R)看特征值是否全大于零，用ctrb(A,B)算能控性矩阵的秩。在SDRE里还有一个不容易发现的问题：A(x)在某个时间步变成奇异矩阵，导致care内部算法无法收敛。

我的做法是在solve_riccati里包一层try-catch，cat检测到异常时沿用上一步的增益矩阵。虽然这会让策略暂时失去最优性，但至少仿真不会中断。对于博弈验证来说，偶尔一步的保守增益完全可以接受。

6.3 ε均衡校验始终不达标怎么办

如果策略偏离扫描出来的最大收益改善总是超过阈值，不要急着增加扫描范围，先检查双方性能指标的权重是否合理匹配。最常见的问题是追方权重优势过大，导致逃方始终能找到大幅改善自身收益的偏离策略。这时候把Q_p适当降低、Q_e适当提高，让双方“势均力敌”，通常两三轮迭代后ε就降下来了。

另外还要注意性能指标中状态权重和控制权重的量纲匹配。追逃问题里状态量是米和米/秒，控制量是米/秒²，量纲差了两个数量级，直接设置的话Q和R数值上很难协调。我习惯先把状态做归一化处理，或者通过量纲折算系数把R设成与Q同一数量级，这样得到的K矩阵物理意义更清晰，ε校验也更容易通过。

6.4 几个实用的调试小技巧

第一，把每个时间步的P矩阵、K矩阵、控制量都存下来，不要只存状态历程。后续做策略验证和分析时，这些中间量比状态曲线更有价值。第二，调试时先不加饱和，观察无约束解的形态，再逐步加入约束，可以很快定位问题是出在策略求解还是出在约束处理。第三，在绘制相对距离曲线的同时，把双方的控制量也画在同一张图上，用双坐标轴显示，往往能一眼看出博弈的“攻防转折点”在哪个时刻。

这套方法不只适用于这个项目。我后来在另一个接近导引策略验证项目里也用到了相同的框架：SDRE求解策略组，然后通过策略偏离扫描来判断策略的鲁棒性。差异只是换了一套状态方程和性能指标，代码框架几乎原封不动就能复用。

我个人在实际调试中最大的体会是：这类博弈仿真最难的从来不是数学推导，也不是代码实现，而是你得先想清楚“双方到底在争什么”。性能指标的定义方式直接决定仿真的结果走向。你把逃方的终端目标设为远离追方的某个参考点，和直接把终端项设为负的脱靶量，得到的两套策略会有本质差别。所以动手写代码之前，先把双方的“输赢标准”想透，后面的工作都是顺水推舟。

最后分享一个小技巧：在做ε均衡验证时，不要只扫描增益系数。加入几种有物理意义的对抗机动模式，比如逃方的正弦横向机动、脉冲规避、末端急转弯，测试结果的参考价值高得多。增益扫描说明的是“策略参数没调好”，而物理机动模式测试说明的是“策略结构能不能扛住针对性反制”，后者才是工程上真正关心的均衡质量。代码跑通了只是开始，能把均衡验证做到让人信服，才算把这个项目吃透了。