一种基于阻力剖面的高升阻比高超平稳滑翔再入制导方法与流程

本发明涉及了一种基于阻力剖面的高升阻比高超平稳滑翔再入制导方法，属于航天技术、武器技术领域。

背景技术：

高超声速滑翔飞行器一般由助推火箭送入预定高度，然后无动力从轨道或亚轨道再入大气层。为了获得远程飞行能力和强横向机动能力，飞行器通常采用大升阻比的升力体结构。

再入飞行是一个多约束问题，具有任务多变、飞行时间长、速度和位置变化大、飞行环境变化剧烈等特点。这给滑翔式飞行器再入制导设计带来了困难和挑战。标准轨道制导律通过跟踪标准轨迹进行制导，其特点是计算量小，易于实现，通常以阻力加速度剖面为参考剖面，通过滑翔距离解析解求解名义弹道预测剩余航程，调整pid增益控制倾侧角跟踪参考剖面。但多数滑翔距离解析解由于忽略了地球自转引起的惯性力，在旋转地球背景下精度不高，存在跟踪误差。

技术实现要素：

本发明的目的在于提供一种基于阻力剖面的高升阻比高超平稳滑翔再入制导方法，针对高超声速飞行器平稳滑翔上述问题，在旋转地球背景下，证明了线化地球自转引起的惯性力方法的精确性。在此基础上，推导出精度更高的滑翔距离解析解，并根据基于能量的解析解规划参考阻力剖面进行跟踪。针对再入弹道存在的长周期弱阻尼问题，采用考虑地球自转的弹道阻尼技术进行抑制。

本发明一种基于阻力剖面的高升阻比高超平稳滑翔再入制导方法，其特征在于：它包括以下几个步骤(如图12所示)：

步骤一：高超声速飞行器建模

在旋转地球背景下，采用标准大气模型和平方反比引力场模型，建立高超声速飞行器的六自由度动力学方程如下：

其中，d/dt表示对时间求导，t是时间，h是海拔高度，r是射程，λ是经度，φ是纬度，v是飞行器相对于地球的速率，γ是弹道倾角，ψ是飞行器航向角，以当地北向为基准，σ是倾侧角，re是地球半径，取值6378.137km，ωe是地球自转角速度，l＝ρv²scl/2m为升力加速度，d＝ρv²scd/2m为阻力加速度。

步骤二：滑翔距离解析解推导；

1.积分形式的射程预测公式

在旋转地球背景下，考虑由地球自转引起的惯性力，基于步骤一中的动力学模型推导积分形式的射程预测公式。

(1)基于速度的射程预测公式

由射程导数公式(2)与式(5)，可得到射程对速度的导数。假设cosγ≈1，sinγ≈0，得到考虑地球自转的基于速度的射程预测公式

(2)基于能量的射程预测公式

由于飞行器在再入过程中是无动力的，能量逐渐衰减，且终端条件也可以由终端速度或终端能量定义，这里使用能量作为独立变量，其好处是不用考虑再入过程开始、结束时间的差异，且可以统一合理地描述不同轨迹段。定义能量为

e＝gh+1/2v²(9)

将上式对高度求导可得

根据式(1)，(2)，(5)及式(10)，经过一系列代数化简，假设cosγ≈1，sinγ≈0，以及量级很小，可得基于能量的考虑地球自转的射程积分公式

2.惯性力线化及证明

由于地球自转引起的牵连加速度和科氏加速度均是时变的，且具有较强的非线性，无法直接积分，故对其进行线化处理。由于平稳滑翔段在再入飞行过程中占很大比例，主要对平稳滑翔段的惯性力进行线化处理。

在平稳滑翔段，弹道倾角很小，令式(6)为0，cosγ≈1，得到平稳滑翔段阻力满足：

上式结合式(8)，式(11)，考虑可近似为常数，故仅需要线化cos²φ，sinφcosφcosψ和cosφsinψ项。cos²φ与sinφcosφcosψ前面的系数(即)很小，约为10^-4级别，所以在此粗略的用能量的线性函数来拟合，拟合公式如下：

sinφcosφcosψ＝k11e+k10(13)

cos²φ＝k21e+k20(14)

不用拟合cosφsinψ是因为其值近似为常数，严谨的证明过程如下：

因为cosφsinψ与速度无关，所以可忽略速度影响，视飞行器仅受重力作用做圆周运动。由于重力对ze轴力矩为0，故可得

在ned坐标系下，由

转到地心赤道旋转坐标系中可得

所以

又因为飞行器仅受重力做匀速圆周运动，且re+h≈re为定值，cosγ≈1，可得

cosφsinψ≈const(20)

故cosφsinψ可视为常值，代入初始值φ0和ψ0求解此常值。

3.平稳滑翔段滑翔距离解析解推导

在推导基于能量的滑翔距离解析解时，式(12)中的速度一次项难以处理，故在此将改速度项线化为能量的一次函数如下：

在此基础上，经过一系列代数变换，得到常值纵向升阻比的平稳滑翔段滑翔距离解析解如下：

(1)基于速度的滑翔距离解析解

另

若δ≥0

若δ＜0

其中，

a2＝(k11+k21)ωe²(re+h)²/2+1(24)

(2)基于能量的滑翔距离解析解

其中，

b1＝2+2kv1ωe(re+h)cosφsinψ+k21ωe²(re+h)²+k11ωe²(re+h)(29)

对比分析可知，基于能量的滑翔距离解析解精度更高，形式更为简洁，故以能量为自变量进行再入轨迹规划。

步骤三：再入轨迹规划与跟踪制导；

1.基准攻角剖面规划

对于再入机动飞行器而言，射程和热防护是攻角剖面设计的主要考虑因素。针对cav-h(highperformancecommonaerovehicle)模型，这里设计攻角曲线为能量分段函数：

其中，es1和es2表示攻角分段函数的分界点。在再入初始下降阶段，为尽量拉高弹道，减小热流密度，采用较大攻角α1飞行在平稳滑翔段，为了发挥最大机动能力，飞行器以最大升阻比攻角滑翔，即α2＝(l/d)max。出于控制系统跟踪性和气动稳定性的考虑，设计es1和es2之间的线性段来保证攻角剖面平稳过渡。

2.阻力-能量剖面规划

由于初始段阻力加速度走廊狭窄，且走廊上边界可以近似为二次曲线，故将初始段的参考d(e)曲线设计成两个关于能量e的二次多项式，形式如下：

其中，a11，a12，a13，a21，a22，a23为常值。

再入走廊中部范围较大，考虑平稳滑翔条件，设计滑翔段参考d(e)曲线为一条直线。根据式(6)，得：

即为关于能量e的一次曲线，比较剩余射程和预测射程大小，确定其与二次多项式阶段交点。其中，l1＝lcosσ，为升力在铅垂平面的分量，l1/d为纵向升阻比，可由步骤二中平稳滑翔段滑翔距离解析解得到。

由于再入走廊末端较为宽阔，在线性段末端规划一次曲线段，以满足过程约束和终端约束。

d＝k3e+b3(e＜e3)(34)

其中，k3，b3为常值。

3.参考弹道跟踪指令

通过对公式(1)和阻力公式求导，进行一系列代数化简，忽略阻力系数cd的一阶、二阶导数，得到如下一个公式，其可以将参考阻力转换为纵向升阻比l1/d。

进而通过调节基准倾侧角σbsl来跟踪纵向升阻比曲线，实现倾侧反转。

飞行器横向运动通过航向误差门限进行控制。

步骤四：基于地球自转模型的平稳滑翔弹道倾角推导；

弹道阻尼技术以平稳滑翔弹道倾角作为反馈项进行控制，可以成功抑制再入弹道存在的长周期弱阻尼振荡。在平稳滑翔段，相比于传统的pid反馈控制，弹道阻尼控制技术跟踪精确性更高，形式简单，便于应用。在旋转地球背景下，根据式(6)推导考虑地球自转的平稳滑翔弹道倾角γ'sg，满足

d2γ'sg²+d1γ'sg+d0＝0(37)

各项系数为：

代入弹道阻尼技术中作为反馈项添加到基准攻角和倾侧角上进行控制：

αcmd＝αbsl+cosσbslkγ(γ'sg-γ)(41)

用上述弹道阻尼控制技术可精确跟踪平稳滑翔弹道，有效抑制弹道振荡。

通过上述四个步骤，可使弹道稳定在平稳滑翔弹道附近，最终得到了一种基于阻力剖面的高升阻比高超平稳滑翔再入制导方法。

本发明一种基于阻力剖面的高升阻比高超平稳滑翔再入制导方法，优点在于：

(1)线化由地球自转引起的惯性力，并给出了严格证明。获得了精度更高的基于能量的滑翔距离解析解，形式简单，便于应用。

(2)以能量为自变量规划阻力剖面，在旋转地球背景下将阻力剖面转换为纵向升阻比，然后调整倾侧角在最大升阻比下进行跟踪，横向运动由航向误差门限控制。

(3)针对再入弹道的振荡问题，推导基于地球自转模型的平稳滑翔弹道倾角，并采用以此为反馈项的弹道阻尼技术抑制振荡。

附图说明

图1是地心赤道旋转坐标系(ger)和当地北东下坐标系(ned)示意图。

图2是滑翔距离解析解对比曲线。

图3(a)是再入高度走廊。

图3(b)是再入阻力加速度走廊。

图4是规划的阻力加速度剖面。

图5是航向误差门限。

图6(a)是反馈控制对比高度曲线。

图6(b)是反馈控制对比阻力剖面跟踪曲线。

图7是七种目标状态的经纬度曲线。

图8是七种目标状态高度射程曲线。

图9是七种目标状态速度时间曲线。

图10是七种目标状态弹道倾角曲线。

图11(a)是目标东南方向高度曲线。

图11(b)是目标东南方向速度时间曲线。

图11(c)是目标东南方向弹道倾角时间曲线。

图11(d)是目标东南方向倾侧角时间曲线。

图11(e)是目标东南方向航向角曲线。

图11(f)是目标东南方向航向角误差曲线。

图12为本发明方法流程图。

上述图中，涉及到的符号、代号说明如下：

图1中，表示地心赤道旋转坐标系(ger)，oxyz表示当地北东下坐标系(ned)。h是海拔高度，r是射程，λ是经度，φ是纬度，v是飞行器相对于地球的速率，γ是弹道倾角，ψ是飞行器航向角，以当地北向为基准，σ是倾侧角。图1中，trajectorysimulation表示弹道仿真结果，agrfree表示基于能量考虑地球自转的滑翔距离解析解，agrfres表示基于速度考虑地球自转的滑翔距离解析解，agrfnree表示基于能量无地球自转的滑翔距离解析解，agrfnres表示基于速度无地球自转的滑翔距离解析解。图4中，e1，e2，e3表示规划的阻力加速度剖面分段点。图5中，withoutcontrol表示不加反馈控制，pidcontrol表示添加pid控制，tdct表示弹道阻尼控制。图9中trajectorysimulation表示弹道仿真结果，referencetrajectory表示参考剖面，entrycorridor表示高度走廊，los表示视线航向角。所有图例中，height、altitude表示高度，range表示射程，longitude表示经度，latitude表示纬度，velocity表示速度，energy表示能量，time表示时间，headingangle表示航向角，headingerror表示航向误差，flight-pathangle表示弹道倾角，downrange表示纵程，dragacceleration表示阻力加速度。

具体实施方式

下面将结合附图和实施案例对本发明作进一步的详细说明。

针对高超声速飞行器平稳滑翔上述问题，在旋转地球背景下，证明了线化地球自转引起的惯性力方法的精确性。在此基础上，推导出精度更高的滑翔距离解析解，并根据基于能量的解析解规划参考阻力剖面进行跟踪。针对再入弹道存在的长周期弱阻尼问题，采用考虑地球自转的弹道阻尼技术进行抑制。

本发明一种基于阻力剖面的高升阻比高超平稳滑翔再入制导方法，其特征在于：它包括以下几个步骤：

步骤一：再入模型建立；

1.再入运动方程

如图1所示，这里在地心赤道旋转坐标系(ger)和当地北东下坐标系(ned)下观察飞行器运动。

其中，d/dt表示对时间求导，t是时间，h是海拔高度，r是射程，λ是经度，φ是纬度，v是飞行器相对于地球的速率，γ是弹道倾角，ψ是飞行器航向角，以当地北向为基准，σ是倾侧角，re是地球半径，取值6378.137km，ωe是地球自转角速度，l＝ρv²scl/2m为升力加速度，d＝ρv²scd/2m为阻力加速度。

2.过程约束

飞行器滑翔飞行中所要满足的热流、动压、过载与平衡滑翔条件等约束共同限制了滑翔弹道的可行范围，即滑翔飞行再入走廊。高超声速飞行器滑翔段的过程约束如下：

其中，kq是热流密度相关常数，是最大热流密度，nmax是飞行器允许得最大过载，qmax为最大来流动压。

3.终端约束

终端约束主要为终端速度约束，终端高度约束，终端能量约束，弹目距离约束，，航向误差约束，终端倾侧角约束以及约束弹道倾角为一0°附近的小量。

步骤二：滑翔距离解析解；

1.积分形式的射程预测公式

(1)基于速度的射程预测公式

由射程导数公式

结合式(4)，可得到射程对速度的导数

由于弹道倾角接近于0，假设cosγ≈1，sinγ≈0，则上述公式可以简化为

(2)基于能量的射程预测公式

由于飞行器在再入过程中是无动力的，能量逐渐衰减，且终端条件也可以由终端速度或终端能量定义，这里使用能量作为独立变量，其好处是不用考虑再入过程开始、结束时间的差异，且可以统一合理地描述不同轨迹段。定义能量为

e＝gh+1/2v²(13)

将上式对高度求导可得

结合式(1)和式(4)可得

将式(15)代入式(14)中，g项被抵消，可得

结合式(1)和式(10)可得

将式(17)代入式(16)中可得

得

考虑弹道倾角很小，即cosγ≈1，sinγ≈0，且量级很小，射程可由下式估算：

2.惯性力线化

由于地球自转引起的牵连加速度和科氏加速度均是时变的，且具有较强的非线性，无法直接积分。故对其进行线化处理。由于平稳滑翔段在再入飞行过程中占很大比例，主要对平稳滑翔段的惯性力进行线化处理。

在平稳滑翔段，弹道倾角很小，令式(5)为0，cosγ≈1，得到平稳滑翔段阻力满足：

上式结合式(12)和式(20)，考虑可近似为常数，故仅需要线化cos²φ，sinφcosφcosψ和cosφsinψ项。cos²φ与sinφcosφcosψ前面的系数(即)很小，约为10^-4级别，所以在此粗略的用能量的线性函数来拟合，拟合公式如下：

sinφcosφcosψ＝k11e+k10(22)

cos²φ＝k21e+k20(23)

不用拟合cosφsinψ是因为其值近似为常数，严谨的证明过程如下：

因为cosφsinψ与速度无关，所以可忽略速度影响，视飞行器仅受重力作用做圆周运动。由于重力对ze轴力矩为0，故可得

在ned坐标系下，由

转到地心赤道旋转坐标系中可得

所以

又因为飞行器仅受重力做匀速圆周运动，且re+h≈re为定值，cosγ≈1，可得

cosφsinψ≈const(29)

故cosφsinψ可视为常值，代入初始值φ0和ψ0求解此常值。

3.平稳滑翔段滑翔距离解析解推导

下面试推导平稳滑翔段滑翔距离解析解。在推导基于能量的滑翔距离解析解时，式(21)中的速度一次项难以处理，故在此将改速度项线化为能量的一次函数如下：

在此基础上，经过一系列代数变换，得到常值纵向升阻比的平稳滑翔段滑翔距离解析解如下表1。

表1

在旋转地球环境下通过仿真对比以上四种滑翔距离解析解的精度。

以下表1为解析解与弹道仿真结果比较

表2

由图2和以下表1可以看出，忽略地球自转影响会大幅降低滑翔距离解析解的精度，因此，推导滑翔距离解析解时考虑地球自转是十分必要的。尽管粗略地线化地球自转引起的惯性力造成一定的误差，但仿真结果表明，线化后得到的基于能量的解析解精度仍然很高。此外，仿真结果表明，无论是否考虑地球自转，基于能量的滑翔距离解析解精度都高于基于速度的滑翔距离解析解。对比式(11)与式(12)，发现二者的区别在于基于速度的滑翔距离解析解简化时忽略了项，而基于能量的滑翔距离解析解忽略了项。可见前者忽略的势能的影响，并且易知前者的数量级远大于后者。在弹道倾角量级较大，即sinγ不严格等于0的情况下，以能量为自变量来预测射程相较以速度为自变量更为准确。同时，对比式(11)与式(12)，易看出基于能量的滑翔距离解析解形式更为简单，故选取基于能量的滑翔距离解析解来进行射程预测。

步骤三：再入轨迹规划与跟踪制导；

1.基准攻角剖面规划

对于再入机动飞行器而言，射程和热防护是攻角剖面设计的主要考虑因素。针对cav-h模型，这里设计攻角曲线为能量分段函数：

在再入初始下降阶段，为尽量拉高弹道，减小热流密度，采用较大攻角α1飞行在平稳滑翔段，为了发挥最大机动能力，飞行器以最大升阻比攻角滑翔，即α2＝(l/d)max。出于控制系统跟踪性和气动稳定性的考虑，设计es1和es2之间的线性段来保证攻角剖面平稳过渡。

2.再入走廊设计

再入走廊上边界由式(5)在和σ＝0的情况下求解得到。再入走廊下边界通过求解过程约束(7)-(9)得到，如图3所示。

3.阻力-能量剖面规划

由于初始段阻力加速度走廊狭窄，且走廊上边界可以近似为二次曲线，故将初始段的参考d(e)曲线设计成两个关于能量e的二次多项式，形式如下：

其中，a11，a12，a13，a21，a22，a23为常值。

再入走廊中部范围较大，考虑平稳滑翔条件，设计滑翔段参考d(e)曲线为一条直线。根据式(5)，得：

即为关于能量e的一次曲线，比较剩余射程和预测射程大小，确定其与二次多项式阶段交点。其中，l1＝lcosσ，为升力在铅垂平面的分量，l1/d为纵向升阻比，可由平稳滑翔段滑翔距离解析解得到。

由于再入走廊末端较为宽阔，在线性段末端规划一次曲线段，以满足过程约束和终端约束。

d＝k3e+b3(e＜e3)(0.34)

其中，k3，b3为常值。规划的阻力加速度剖面如图4所示。

4.参考弹道跟踪指令

通过对公式(1)和阻力公式求导，进行一系列代数化简，忽略阻力系数cd的一阶、二阶导数，得到如下一个公式，其可以将参考阻力转换为纵向升阻比l1/d。

进而通过调节基准倾侧角σbsl来跟踪纵向升阻比曲线，实现倾侧反转。

5.倾侧角反转

侧向制导方法利用预先设计的航向误差门限控制倾侧角反转，将飞行器的侧向运动限定再预先设定的走廊内，当航向误差到达走廊边界时，令倾侧角反向，使横向运动向航向角误差减小的方向进行。

定义航向误差为

δψ＝ψ-ψlos(0.37)

其中，ψlos为视线方位角，即为在飞行器当前铅垂面内，飞行器指向目标的向量与当地北向的夹角。

为了保证δψ能够进入航向误差走廊，再入初始阶段飞行器应朝终端目标点飞行。进入再入阶段后，倾侧角符号由下式确定：

其中，sign[φ(tk-1)]表示上一时刻由侧向制导算法给出的倾侧角符号，δψh(e)为航向误差门限。通常根据横程大小、倾侧角反转次数和终端误差等多方面约束将δψh(e)设计为能量的分段线性函数。按照终端航向角误差要求，将航向误差门限设计为如图5所示的漏斗形。

步骤四：基于地球自转模型的平稳滑翔弹道倾角推导；

采用基准攻角和基准倾侧角跟踪参考阻力-能量剖面，对比传统的pid反馈控制和yu提出的平稳滑翔弹道倾角作为反馈项进行控制的弹道阻尼控制方法。所得到的弹道如图6所示，可以看出，弹道阻尼控制技术和pid控制都可以成功抑制弹道长周期弱阻尼振荡。但图6(b)显示，在平稳滑翔段，pid控制的阻力加速度跟踪曲线存在较大的跟踪误差，而弹道阻尼控制技术则可以较为准确的跟踪参考剖面，优势明显，故选用弹道阻尼技术抑制弹道振荡。

在旋转地球背景下，根据式(5)推导考虑地球自转的平稳滑翔弹道倾角γ'sg，满足

d2γ'sg²+d1γ'sg+d0＝0(39)

各项系数为：

代入弹道阻尼技术中作为反馈项添加到基准攻角和倾侧角上进行控制：

αcmd＝αbsl+cosσbslkγ(γsg-γ)(43)

用上述弹道阻尼控制技术可精确跟踪平稳滑翔弹道，有效抑制弹道振荡。

实施案例：

为了验证算法的精度，以cav-h作为再入飞行器模型在七种不同目标状态下进行仿真。再入初始状态均为h0＝80km，v0＝7000m/s，λ0＝0°，θ0＝45°，γ0＝0°。ψ0依据目标方向取值。弹道阻尼反馈项系数kγ＝5。终端约束为vf＝2000m/s，hf＝25km，弹目距离rtm≤5km。终端航向角误差|δψf|≤5°。仿真结果如图7-图10。

表3

图7中，按照逆时针方向从左至右，分别是七种不同的目标位置，即从西北方位到东北方位。图8-图10和表3(终端状态量)显示了弹道阻尼控制技术能控制飞行器在不同目标状态下成功抑制弹道振荡，精准命中目标。以东南方向的目标位置为例具体分析参考剖面跟踪情况。

由图11(a)可明显看出，所得弹道处于狭窄的高度走廊内，满足过程约束。同时，相较于原先再入弹道，弹道振荡的情况得到了极大的缓解，弹道的平稳性大大增强。但在初始下降阶段，弹道仍存在一个波峰的振荡，可通过增大弹道阻尼系数来消除。图11(b)中，实际速度基本跟踪上了参考速度，达到了末端速度要求。图11(c)中弹道倾角精确跟踪参考弹道倾角，基本稳定在0附近，满足要求。图11(d)为倾侧角时间曲线，结合图11(e)和图11(f)，可看出在仿真过程中，倾侧角翻转了6次，末端横向误差|δψf|≤5°，基本命中目标点。

综上所述，通过上述步骤，推导出了本发明方法，即一种基于阻力剖面的高升阻比高超平稳滑翔再入制导方法，案例仿真结果表明本发明方法能够精准跟踪平稳滑翔弹道，成功抑制弹道振荡，提高导弹的命中精度，综合性能优异。