1. 清华大学 航天航空学院, 北京 100084;
2. 北京航空航天大学 自动化科学与电气工程学院, 北京 100191
收稿日期：2019-03-06
基金项目：国家自然科学基金资助项目（11672146）
作者简介：程林(1989-), 男, 博士后
通信作者：蒋方华, 副教授, E-mail:jiangfh@tsinghua.edu.cn
摘要：针对飞行器再入制导问题，该文引入控制变量参数化、积分问题转换和在线模型辨识等技术，提出一种跨周期迭代的可行轨迹预测校正算法，并结合标称轨迹跟踪算法形成一套多约束复合制导方案。利用一种复合高度-速度（height velocity，HV）飞行走廊，将再入轨迹规划问题简化为单调函数寻根问题。为提高射程预测计算效率，引入Gauss-Legendre求积公式，将积分问题转化为函数计算问题。采用递推最小二乘估计方法，收集历史预测信息，实现模型在线辨识功能，并采用跨周期Newton-Raphson方法完成高度权重系数的在线修正。在标称轨迹跟踪器设计的基础上，开展飞行器数值仿真试验，结果表明：基于在线模型辨识的复合制导方法具有显著的速度优势，且具有优异的自主性和自适应能力。
关键词：再入制导复合飞行走廊Gauss求积法递推最小二乘估计模型辨识跨周期参数校正
Multi-constraint compound reentry guidance based on onboard model identification
CHENG Lin1, ZHANG Qingzhen2, JIANG Fanghua1
1.School of Aerospace Engineering, Tsinghua University, Beijing 100084, China;
2.School of Automation Science and Electrical Engineering, Beihang University, Beijing 100191, China
Abstract: A period-crossing feasible trajectory planning algorithm for reentry guidance was developed based on control variable parameterization, integral transformations, and onboard model identification. A compound height velocity (HV) corridor simplifies the reentry guidance problem into a root-searching problem. A Gauss integral is introduced to improve the time efficiency of the range prediction with the original integral problem converted into a function calculation problem. The recursive least squares estimation method was used to develop functions for on-board information mining and model identification. The reliable, explicit solution model can easily correct the weight coefficients using the period-crossing Newton-Raphson method. Numerical simulations show that the reentry guidance method based on on-board model identification is much faster, more autonomous and more adaptable than the reference trajectory tracking design method.
Key words: reentry guidancecompound altitude-velocity corridorGaussian integral methodrecursive least squares estimation methodmodel identificationperiod-crossing parameter correction
“再入制导”是指高超声速飞行器按照一定的导引规律将飞行器从初始投放点引导到目标状态点，期间满足多种过程约束、控制约束和终端约束。从本质上讲，轨迹优化设计问题是一个复杂的非线性最优控制问题，它的实时性规划一直是难题。传统的标称轨迹制导方法采用离线轨迹设计与在线标称轨迹跟踪相结合的制导方案，此举规避了轨迹在线设计的难题，使得再入制导工程化成为可能[1]。然而，标称轨迹制导方法的离线轨迹设计的方案决定了它只能应用于特定场景下的特定任务，从根本上限制了制导的鲁棒性以及任务的灵活性，随着近年来任务需求的不断增强，它的弊端不断突显，包括耗费人力物力、鲁棒性有限、任务灵活性低等。
20世纪90年代末，在新一代可重复使用运载飞行器发展背景下，美国NASA牵头启动了先进制导与控制计划(AG & C)[2], 旨在为下一代飞行器提供更加自主、可靠与省钱的再入制导，掀起了新型再入制导方法研究的热潮。新型再入制导技术的发展主要围绕2个方向，包括标称轨迹在线重构[3-7]以及预测校正制导方法[8-12]。2种方法中，可行轨迹在线快速设计都是技术关键。1997年，Lu等[13-14]引用伪谱法概念以阻力为自变量，以总吸热量最小为优化目标，采用成熟的序列二次规划方法实现快速阻力加速度生成。2003年，Shen等[10]提出了拟平衡滑翔条件，基于约束生成的倾侧角走廊实现特定射程的倾侧角剖面快速生成。2014年，Lu[15]在总结近年来轨迹设计结果的基础上，进一步提出了一套通用的轨迹快速规划方案，并率先提出采用时标分离理论对轨迹进行修正，可以得到比较平缓的再入轨迹。在文[16]的工作基础上，雍恩米等[17]进一步研究考虑地球旋转影响和大横程的三维参考轨迹快速生成。以上方法取得不同技术突破，但是由于设计参数与设计目标之间没有显式的函数关系，算法大多采取周期内多次弹道积分，并依据轨迹预测信息提取出来的梯度信息迭代到目标参数值。这种周期内多次轨迹预测和迭代的方法给实际制导工程化带来严峻的挑战。与直接打靶法相比，伪谱法同时离散状态和控制变量，不需要对状态积分，具有求解速度更快的显著优势，在轨迹快速规划方面具有广泛应用前景，受到研究者的持续关注[18-22]。
结合以上研究现状分析和问题描述，本文提出一种基于在线模型辨识的可行轨迹实时规划算法被提出。该算法将不断产生于周期间的轨迹预测信息收集起来，并采用递推最小二乘估计算法对求解模型在线辨识。辨识得到的模型可以在周期间传递，进而应用于可行轨迹的参数在线修正。不断修正的可行轨迹与经典的标称轨迹跟踪制导方法结合起来，可以构成一种新型复合制导方案。制导多任务仿真验证了整个方案的有效性和可行性。
1 再入制导问题1.1 再入运动模型本论文以高超声速滑翔式再入飞行器为研究对象，其考虑地球自转的3自由度再入运动模型表示如下[23]：
| $\left\{ \begin{array}{l}\dot r = v\sin \theta ,\\\dot \lambda = \frac{v}{{r\cos \phi }}\cos \theta \sin \psi ,\\\dot \phi = \frac{v}{r}\cos \theta \cos \psi ,\\{{\dot S}_e} = \frac{v}{r}\cos \theta ,\\v = - \frac{D}{{{m_0}}} - g\sin \theta + \omega _e^2r\cos \phi \left( {\sin \theta \cos \phi - \cos \theta \sin \phi \cos \psi } \right),\\\dot \theta = \frac{1}{v}\left[ {\frac{{L\cos \sigma }}{{{m_0}}} + \left( {\frac{{{v^2}}}{r} - g} \right)\cos \theta + 2{\omega _e}v\cos \phi \sin \psi + } \right.\\\;\;\;\;\;\left. {\omega _e^2r\cos \phi \left( {\cos \theta \cos \phi + \sin \theta \sin \phi \cos \psi } \right)} \right],\\\dot \psi = \frac{1}{v}\left[ {\frac{{L\sin \sigma }}{{{m_0}\cos \theta }} + \frac{{{v^2}}}{r}\cos \theta \sin \psi \tan \phi + \frac{{\omega _e^2r}}{{\cos \theta }}\cos \psi \sin \phi \cos \phi - } \right.\\\;\;\;\;\;\;\left. {2{\omega _e}v\left( {\cos \phi \cos \psi \tan \theta - \sin \phi } \right)} \right].\end{array} \right.$ | (1) |
| ${\alpha _{{\rm{ref}}}} = \left\{ \begin{array}{l}{\alpha _{\max }},\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;v > {v_1};\\\frac{{{\alpha _{\max }} - {\alpha _{\min }}}}{{{v_1} - {v_2}}}\left( {v - {v_2}} \right),\;\;\;\;\;{v_2} < v < {v_1};\\{\alpha _{\min }},\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;v < {v_2}.\end{array} \right.$ | (2) |
1.2 再入轨迹约束条件再入过程约束主要包括热峰值约束、过载约束和动压约束[16]。
| $\dot Q\left( \tau \right) = \frac{{{C_1}}}{{\sqrt {{R_{\rm{d}}}} }}{\left( {\frac{\rho }{{{\rho _0}}}} \right)^{0.5}}{\left( {\frac{v}{{{V_C}}}} \right)^{3.15}} \le {\dot Q_{\max }},$ | (3) |
| $n\left( \tau \right) = \sqrt {{{\bar L}^2} + {{\bar D}^2}} = q\sqrt {C_L^2 + C_D^2} \frac{{{S_{{\rm{ref}}}}}}{{{m_0}{g_0}}} \le {n_{\max }},$ | (4) |
| $q\left( \tau \right) = \frac{1}{2}\rho {v^2} \le {q_{\max }}.$ | (5) |
考虑到与终端能量管理端的交接，再入段轨迹设计的终端需满足给定高度、速度条件以及位置条件，即
| $\begin{array}{*{20}{c}}{h\left( {{\tau _f}} \right) = {h_f},v\left( {{\tau _f}} \right) = {v_f},}\\{\lambda \left( {{\tau _f}} \right) = {\lambda _f},\phi \left( {{\tau _f}} \right) = {\phi _f}.}\end{array}$ | (6) |
1.3 仿真对象为验证本论文提出算法的有效性，引用美国波音公司1998年设计的带控制翼的锥形体再入机动飞行器CAV-L作为仿真对象，总体参数如表 1所示，其他总体气动参数详见文[17]和[23]。
表 1 初始再入状态与约束设置
| 参数 | 取值 |
| 初始高度/km | 100 |
| 初始航迹角/(°) | -2 |
| 终端速度/(m·s-1) | 1 800 |
| 目标经度/(°) | 225 |
| 1 200 | |
| qmax/kPa | 200 |
| 初始速度/(m·s-1) | 7 200 |
| 初始航向角/(°) | 55 |
| 终端高度/km | 20 |
| 目标纬度/(°) | 5 |
| nmax | 4 |
表选项
2 飞行器控制参数化在初始投放阶段，大气密度稀薄，飞行器机动能力有限，因此飞行器不作轨迹调整，沿固定倾侧角飞行。整个再入制导过程分为初始下滑段和拟平衡滑翔段，初始下滑段固定倾侧角设计方法和截止条件详见文[9]，在此不再赘述。
飞行器进入拟平衡滑翔段，将有足够的气动能力进行轨迹调整，所以本文的可行轨迹预测校正算法针对的是平衡滑翔段。
2.1 复合HV走廊飞行器高度-速度(height velocity, HV)走廊是轨迹设计中最常用的一种飞行走廊，它具有求解简单，物理含义明确的优点，且在实际轨迹设计中应用方便。飞行器复合走廊由常规HV走廊HCUB和HCLB、考虑末端状态的末段飞行包络HTUB和HTLB复合而成，具体设计过程详见文[23]，如图 1所示。复合走廊上界记为HUb，下界记为HLb。复合走廊为飞行器提供了可行的飞行空间，原则上，在走廊内部，飞行器满足滑翔条件、过程约束和终端高度、速度约束。复合走廊的构建为飞行器可行轨迹实时解算奠定了基础。
 |
| 图 1 (网络版彩图)复合HV走廊 |
| 图选项 |
2.2 射程与高度的关系根据式(1)，并忽略地球自转项。可以得到射程与速度的简化导数关系为
| $\frac{{{\rm{d}}{S_{\rm{e}}}}}{{{\rm{d}}v}} = \frac{{v\cos \theta /r}}{{ - D/{m_0} - g\sin \theta }}.$ | (7) |
| $\frac{{{\rm{d}}{S_{\rm{e}}}}}{{{\rm{d}}v}} = \frac{{{m_0}}}{{ - \frac{1}{2}\rho v{C_L}{S_{{\rm{ref}}}}{R_0}}}.$ | (8) |
2.3 控制剖面参数化可行轨迹快速规划的本质是如何在有限时间内规划一条满足多种约束条件的可行剖面。本文在复合走廊构建的基础上，提出将复合走廊上下界加权得到高度-速度剖面的可行轨迹生成方法，表示为
| ${h_{{\rm{ref}}}} = \omega \cdot {H_{{\rm{UB}}}} + \left( {1 - \omega } \right) \cdot {H_{{\rm{LB}}}}.$ | (9) |
3 轨迹在线迭代预测校正算法3.1 基于Gauss积分法的射程快速预测鉴于权重系数ω值大小可实现设计剖面href的高低调整，从而影响射程，所以ω与射程之间可以建立一定的量化函数关系。
已知高度剖面，射程的推导公式为
| $\begin{array}{*{20}{c}}{{S_{{\rm{pre}}}} = \int_{{v_{{\rm{now}}}}}^{{v_f}} {\frac{{v\cos \theta /\left( {{R_0} + {h_{{\rm{ref}}}}} \right)}}{{ - D - g\sin \theta }}{\rm{d}}v} = }\\{\int_{{v_{{\rm{now}}}}}^{{v_f}} {f\left( {{h_{{\rm{ref}}}},v} \right){\rm{d}}v} .}\end{array}$ | (10) |
| $\sin \theta = \frac{{ - D \cdot \left( {{\rm{d}}{h_{{\rm{ref}}}}/{\rm{d}}v} \right)}}{{{m_0}\left( {g \cdot \left( {{\rm{d}}{h_{{\rm{ref}}}}/{\rm{d}}v} \right) + v} \right)}}.$ | (11) |
由式(10)可知，射程的估算需要射程积分，而飞行器远程滑翔的特性决定了射程积分的代价是很高的，因此本文引进Lagrange多项式对被积函数进行拟合，配点选作LGL点，它们是多项式
| ${P_K}\left( \tau \right) = \frac{1}{{{2^K}K!}}\frac{{{\rm{d}}k}}{{{\rm{d}}{\tau ^K}}}\left[ {{{\left( {{\tau ^2} - 1} \right)}^K}} \right],K = 0,1, \cdots ,n.$ | (12) |
| $f\left( {{h_{{\rm{ref}}}},v} \right) \approx \sum\limits_{i = 1}^K {{{\tilde L}_i}\left( \tau \right)f\left( {{H_{{\rm{ref}}}}\left( {{\tau _i}} \right),{v_i}} \right)} .$ | (13) |
最终，射程预测公式可以改进为
| ${S_{{\rm{pre}}}} = \int_{{v_{{\rm{Now}}}}}^{{v_f}} {\sum\limits_{i = 1}^K {{{\tilde L}_i}\left( \tau \right)f\left( {{H_{{\rm{ref}}}}\left( {{\tau _i}} \right),{v_i}} \right){\rm{d}}v} } .$ | (14) |
| ${S_{{\rm{pre}}}} = \frac{{{v_f} - {v_{{\rm{Now}}}}}}{2}\sum\limits_{i = 1}^K {{w_k}f\left( {{H_{{\rm{ref}}}}\left( {{\tau _i}} \right),{v_i}} \right)} .$ | (15) |
对比式(10)和(15)可以看出，通过对积分函数f(href, v)采用Lagrange多项式进行拟合，最终将积分问题转换为两个向量内积问题，从而大大提高算法运行速度。
3.2 基于递推最小二乘估计的模型在线辨识传统的轨迹校正方法是单个周期内采用优化手段确定满足约束的控制参数，经典的如Lu等[15]采用Gauss-Newton求解满足射程的初值倾侧角。传统的轨迹校正方法存在两点不足。首先，周期与周期之间虽然参数可以传递，但是模型信息是分割的，也就是说之前求解的信息没有充分利用起来。其次，也正是每个周期迭代开始时，对模型了解不充分，参数求解需要多次迭代，才能结合梯度或者模型信息，求解出高精度的解，这是由问题的属性决定的。
结合式(9)和(15)，可以建立走廊权重系数w与预测射程的函数关系为
| ${S_{{\rm{pre}}}} = f\left( {\omega ,{x_{{\rm{Now}}}}} \right).$ | (16) |
基于以上理由，需要采用一套可靠、在线运行的模型辨识算法对模型进行持续辨识与传递，并加以利用。考虑Spre与ω为非线性关系，在当前状态下按照小扰动理论展开为
| ${S_{{\rm{pre}}}} = f\left( {\omega ,{x_{{\rm{Now}}}}} \right) \approx a \cdot \omega + b.$ | (17) |
考虑线性离散观测模型为
| $\mathit{\boldsymbol{h}}_i^{\rm{T}}\mathit{\boldsymbol{x}} + {v_i} = {z_i}.$ | (18) |
| $\mathit{\boldsymbol{H}}{ = _k} + {\mathit{\boldsymbol{V}}_k} = {\mathit{\boldsymbol{Z}}_k}.$ | (19) |
| $\begin{array}{*{20}{c}}{{\mathit{\boldsymbol{H}}_k} = \left[ {\begin{array}{*{20}{c}}{\mathit{\boldsymbol{h}}_1^{\rm{T}}}\\ \vdots \\{\mathit{\boldsymbol{h}}_k^{\rm{T}}}\end{array}} \right] \in {R^{k \times 2}},{\mathit{\boldsymbol{V}}_k} = \left[ {\begin{array}{*{20}{c}}{{v_1}}\\ \vdots \\{{v_k}}\end{array}} \right] \in {R^k},}\\{{\mathit{\boldsymbol{Z}}_k} = \left[ {\begin{array}{*{20}{c}}{{z_1}}\\ \vdots \\{{z_k}}\end{array}} \right] \in {R^k}.}\end{array}$ | (20) |
| ${\mathit{\boldsymbol{\hat x}}_k} = {\left( {\mathit{\boldsymbol{H}}_k^{\rm{T}}\mathit{\boldsymbol{Q}}_k^{ - 1}{\mathit{\boldsymbol{H}}_k}} \right)^{ - 1}}\mathit{\boldsymbol{H}}_k^{\rm{T}}\mathit{\boldsymbol{Q}}_k^{ - 1}{\mathit{\boldsymbol{Z}}_k}.$ | (21) |
在线参数估计中，新的信息不断得到，本文引进递推的最小二乘估计来实现参数的在线修正。假设在时刻得到新的观测值zk+1，即
| $\mathit{\boldsymbol{h}}_{k + 1}^{\rm{T}}\mathit{\boldsymbol{x}} + {v_{k + 1}} = {z_{k + 1}},$ | (22) |
| ${\mathit{\boldsymbol{P}}_k} = {\left( {\mathit{\boldsymbol{H}}_k^{\rm{T}}\mathit{\boldsymbol{Q}}_k^{ - 1}{\mathit{\boldsymbol{H}}_k}} \right)^{ - 1}}.$ | (23) |
| ${\mathit{\boldsymbol{\hat x}}_k} = {\mathit{\boldsymbol{P}}_k}\mathit{\boldsymbol{H}}_k^{\rm{T}}\mathit{\boldsymbol{Q}}_k^{ - 1}{\mathit{\boldsymbol{Z}}_k}.$ | (24) |
| $\mathit{\boldsymbol{Q}}_{k + 1}^{ - 1} = \left[ {\begin{array}{*{20}{c}}{\mathit{\boldsymbol{Q}}_k^{ - 1}}&{}\\{}&{R_{k + 1}^{ - 1}}\end{array}} \right].$ | (25) |
| ${\mathit{\boldsymbol{K}}_{k + 1}} = {\mathit{\boldsymbol{P}}_k}{\mathit{\boldsymbol{h}}_{k + 1}}{\left( {{R_{k + 1}} + \mathit{\boldsymbol{h}}_{k + 1}^{\rm{T}}{\mathit{\boldsymbol{P}}_k}{\mathit{\boldsymbol{h}}_{k + 1}}} \right)^{ - 1}},$ | (26) |
| ${\mathit{\boldsymbol{P}}_{k + 1}} = \left( {\mathit{\boldsymbol{I}} - {\mathit{\boldsymbol{K}}_{k + 1}}\mathit{\boldsymbol{h}}_{k + 1}^{\rm{T}}} \right){\mathit{\boldsymbol{P}}_k}.$ | (27) |
| ${\mathit{\boldsymbol{\hat x}}_{k + 1}} = {\mathit{\boldsymbol{\hat x}}_k} + {\mathit{\boldsymbol{K}}_{k + 1}}\left( {{z_{k + 1}} - \mathit{\boldsymbol{h}}_{k + 1}^{\rm{T}}{{\mathit{\boldsymbol{\hat x}}}_k}} \right).$ | (28) |
根据参数的最优估计值为
| $\left[ {\begin{array}{*{20}{c}}{\hat a}\\{\hat b}\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}{{{\mathit{\boldsymbol{\hat x}}}_{k + 1}}\left( 1 \right)}\\{{{\mathit{\boldsymbol{\hat x}}}_{k + 1}}\left( 2 \right)}\end{array}} \right],$ | (29) |
| ${S_{{\rm{pre}}}} \approx \hat a \cdot \omega + \hat b.$ | (30) |
| ${\omega _{k + 1}} = {\omega _k} + \left( {{S_{{\rm{pre}}}}\left( {{\omega _k}} \right) - {S_{{\rm{go}}}}} \right)/\hat a.$ | (31) |
| ${h_{{\rm{ref}}}}\left( {k + 1} \right) = {\omega _{k + 1}} \cdot {\mathit{\boldsymbol{H}}_{{\rm{UB}}}} + \left( {1 - {\omega _{k + 1}}} \right) \cdot {\mathit{\boldsymbol{H}}_{{\rm{LB}}}}.$ | (32) |
 |
| 图 2 轨迹在线预测-校正示意图 |
| 图选项 |
4 跟踪控制器设计本文设计的可行剖面是高度剖面，考虑到高度剖面变化太快会引起倾侧角跟踪指令的大幅跳跃，因此在轨迹跟踪器之前对高度进行滤波，滤波算法为
| $\begin{array}{*{20}{c}}{h_{{\rm{cmd}}}^k = \left( {1 - {k_{{\rm{Pf}}}}} \right) \cdot h_{{\rm{cmd}}}^{k - 1} + {k_{{\rm{Pf}}}} \cdot h_{{\rm{ref}}}^k + }\\{{k_{{\rm{If}}}}\int {\left( {h_{{\rm{cmd}}}^{k - 1} - h_{{\rm{ref}}}^{k - 1}} \right){\rm{d}}t} .}\end{array}$ | (33) |
高度的跟踪是通过调整倾侧角实现的。首先对高度进行连续二次求导，得到
| $\ddot h = \left[ {\frac{{L\cos \sigma }}{m} + \left( {\frac{{{v^2}}}{r} - g} \right)} \right]\cos \theta - \left( {\frac{D}{m} + g\sin \theta } \right)\sin \theta .$ | (34) |
| $\Delta \ddot h + 2 \cdot {\lambda _h}\Delta \dot h + \lambda _h^2\Delta h = 0.$ | (35) |
5 数字仿真测试为验证本文提出算法的有效性，开展两组试验，分别测试算法的实时性、可靠性和标准环境下的任务自主能力。
5.1 算法实时性测试此部分试验主要对比3种参数迭代方法的求解时间，3种方法分别为基于积分射程预测+Gauss-Newton方法求解(I & GN)，基于Gauss积分法+Gauss-Newton方法求解(G & GN)，以及基于Gauss积分公式+跨周期New-Raphson方法求解(G & PCNR)。
从求解时间对比来看，I & GN与G & GN求解时间的对比主要是显示Lagrange插值多项式对积分对象拟合的效果，对比显示求解速度提高约46倍。G & GN与G & PCNR的时间对比主要显示Gauss-Newton方法与本文提出的基于模型在线辨识的跨周期New-Raphson方法的对比，前后信息的综合利用可以提高算法约6.2倍。3种方法的实时性在表 2中更加直观。
表 2 3种方法实时性对比
| 射程/km | 求解时间/s | ||
| I & Gn | G & Gn | I & PCNR | |
| 2 000 | 8.50 | 0.187 | 0.027 1 |
| 3 000 | 8.17 | 0.182 | 0.027 1 |
| 4 000 | 7.48 | 0.171 | 0.027 1 |
| 5 000 | 6.62 | 0.125 | 0.027 1 |
| 6 000 | 6.73 | 0.156 | 0.027 1 |
| 7 000 | 7.87 | 0.171 | 0.027 1 |
| 7 862 | 6.93 | 0.187 | 0.027 1 |
表选项
5.2 标准环境下多任务仿真试验本文随机选取4个典型再入过程，再入点分别为case1(160, 5)、case2(180, 5)、case3(175, 10)、case4(175, 0)。
图 3—6给出了4个单独任务的仿真图，从图中可以看出：首先，从制导效果来看，飞行器能实现精确地再入飞行，横侧向控制良好，且图 4给出的高度-速度剖面可以看出，飞行器平稳再入滑翔，虽然存在轨迹超过上界，短时间内不满足平衡滑翔的状况，但是软约束不会影响飞行的安全，且所有飞行轨迹高于走廊下界，说明任务再入过程都严格满足过程约束限制，此外终端高度和速度约束也较好地得到满足。其次，从算法迭代稳定性来看，算法能够迅速收敛到目标值，且随着射程需求的改变，ω能够稳定调整，从高度-速度剖面上看，不同射程需求对应不同的飞行高度剖面。再次，为消除横向机动带来的误差损失，图 5中射程系数ω随着飞行慢慢变大，这也验证了制导算法在线自主调整能力。图 6给出了算法的射程调整精度和调节时间统计。可以看出，整个再入过程，飞行器通过对过去信息的挖掘及建模，可实现非常准确的跨周期参数迭代，迭代精度和实时性都非常高。
 |
| 图 3 (网络版彩图)经纬度 |
| 图选项 |
 |
| 图 4 (网络版彩图)高度-速度剖面 |
| 图选项 |
 |
| 图 5 (网络版彩图)权重系数ω的收敛曲线 |
| 图选项 |
 |
| 图 6 (网络版彩图)射程偏差和求解时间 |
| 图选项 |
6 结论高超声速再入制导轨迹设计是一个复杂的多约束优化问题，本文提出一套跨周期运行的轨迹预测校正算法，并结合轨迹跟踪控制器设计，实现飞行器自主和高精度再入复合制导。
本文基于复合HV走廊将复杂轨迹设计问题转化为一元函数求根问题，在此基础上，引入Lagrange多项式拟合被积函数，实现轨迹射程的快速预测，并采用递推最小二乘估计方法实现在线信息挖掘与模型辨识，最终实现可靠求解模型的传递。基于已知求解模型，提出跨周期迭代算法，实现求解参数的跨周期、高精度迭代求解。综合以上技术手段，最终将再入制导轨迹预测校正算法控制在严格制导周期以内(MATLAB仿真下小于0.1 s)。结合轨迹跟踪和多任务仿真验证了算法的实时性和自主性。
参考文献
| [1] | HARPOLD J C, GRAVES JR C A. Shuttle entry guidance[C]// Proceedings of the 25th American Astronautical Society, Anniversary Conference. Houston, USA: NASA, 1978: 35-40. |
| [2] | HANSON J. Advanced guidance and control project for reusable launch vehicles[C]// Proceedings of the AIAA Guidance, Navigation, and Control Conference and Exhibit. Dever, USA: AIAA, 2000: 1-10. |
| [3] | MEASE K D, CHEN D T, TEUFEL P, et al. Reduced-order entry trajectory planning for acceleration guidance[J]. Journal of Guidance, Control, and Dynamics, 2002, 25(2): 257-266. |
| [4] | KANG Y L, CHENG L, ZHANG Q Z, et al. Data-driven RLV multi-objective reentry trajectory optimization based on new QABC algorithm[J]. The International Journal of Advanced Manufacturing Technology, 2016, 84(1-4): 453-471. DOI:10.1007/s00170-015-8124-9 |
| [5] | FAHROO F, ROSS I M. Direct trajectory optimization by a Chebyshev pseudospectral method[J]. Journal of Guidance, Control, and Dynamics, 2002, 25(1): 160-166. |
| [6] | BOLLINO K P, LEWIS L R, SEKHAVAT P, et al. Pseudospectral optimal control: a clear road for autonomous intelligent path planning[C]// Proceedings of the AIAA InfoTech at Aerospace Conference and Exhibit. Rohnert Park, USA: AIAA, 2007: 1228-1241. |
| [7] | REA J R. A Legendre pseudospectral method for rapid optimization of launch vehicle trajectories[D]. Cambridge: Massachusetts Institute of Technology, 2001. |
| [8] | LU P. Entry guidance: a unified method[J]. Journal of Guidance, Control, and Dynamics, 2014, 37(3): 713-728. DOI:10.2514/1.62605 |
| [9] | SHEN Z J, LU P. Onboard generation of three-dimensional constrained entry trajectories[J]. Journal of Guidance, control, and Dynamics, 2003, 26(1): 111-121. |
| [10] | 赵頔, 沈作军. 基于在线轨迹迭代的自适应再入制导[J]. 北京航空航天大学学报, 2016, 42(7): 1526-1535. ZHAO D, SHEN Z J. Adaptive reentry guidance based on on-board trajectory iterations[J]. Journal of Beijing University of Aeronautics and Astronautics, 2016, 42(7): 1526-1535. (in Chinese) |
| [11] | ROSS I M, SEKHAVAT P, FLEMING A, et al. Pseudospectral feedback control: Foundations, examples and experimental results[C]// Proceedings of the 2006 AIAA Guidance, Navigation, and Control Conference. Keystone, USA: AIAA, 2006: 78-99. |
| [12] | LU P. Entry guidance and trajectory control for reusable launch vehicle[J]. Journal of Guidance, Control, and Dynamics, 1997, 20(1): 143-149. |
| [13] | LU P, HANSON J M. Entry guidance for the X-33 vehicle[J]. Journal of Spacecraft and Rockets, 1998, 35(3): 342-349. DOI:10.2514/2.3332 |
| [14] | LU P, BRUNNER C W, STACHOWIAK S J, et al. Verification of a fully numerical entry guidance algorithm[J]. Journal of Guidance, Control, and Dynamics, 2017, 40(2): 230-247. |
| [15] | 雍恩米.高超声速滑翔式再入飞行器轨迹优化与制导方法研究[D].长沙: 国防科学技术大学, 2008: 3-6. YONG E M. Study on trajectory optimization and guidance approach for hypersonic glide-reentry vehicle[D]. Changsha: National University of Defense Technology, 2008: 3-6. (in Chinese) |
| [16] | FAHROO F, DOMAN D B, NGO A D. Footprint generation for reusable launch vehicles using a direct pseudospectral method[C]// Proceedings of the 2003 American Control Conference. Denver, USA: IEEE, 2003: 2163-2168. |
| [17] | LU P. Closed-form control laws for linear time-varying systems[J]. IEEE Transactions on Automatic Control, 2000, 45(3): 537-542. DOI:10.1109/9.847739 |
| [18] | JORRIS T R, MCCRACKEN B. Aircraft system identification using pseudospectral parameter optimization with adaptive nodes[C]// Proceedings of AIAA Atmospheric Flight Mechanics Conference. Portland, USA: AIAA, 2011: 178-201. |
| [19] | BITTNER M, FISCH F, HOLZAPFEL F, et al. A multi-model Gauss pseudospectral optimization method for aircraft trajectories[C]// Proceedings of AIAA Atmospheric Flight Mechanics Conference. Minneapolis, USA: AIAA, 2012: 329-340. |
| [20] | TAYLOR C P. Optimization study of a trans-Atlantic abort for the U.S. space shuttle using a pseudospectral Legendre method[D]. Cambridge: Massachusetts Institute of Technology, 2003. |
| [21] | STANTON S A. Optimal orbital transfer using a Legendre pseudospectral method[D]. Cambridge: Massachusetts Institute of Technology, 2003. |
| [22] | CHENG L, WANG Z B, CHENG Y, et al. Multi-constrained predictor-corrector reentry guidance for hypersonic vehicles[J]. Proceedings of the Institution of Mechanical Engineers, Part G: Journal of Aerospace Engineering, 2018, 232(16): 3049-3067. DOI:10.1177/0954410017724185 |
| [23] | PHILLIPS T H. A Common Aero Vehicle (CAV) model, description, and employment guide[R]. Schafer Corporation for AFRL and AF-SPC, 2003. |
