CN104217041A - 一种多约束在线高斯伪谱末制导方法 - Google Patents

Abstract

本发明涉及一种多约束在线高斯伪谱末制导方法，该方法包括如下步骤：（1）采用高斯伪谱法对非线性状态方程进行末端轨迹优化；（2）在一个制导周期内，采用时钟保持机制将最优控制指令带入实际飞行过程；（3）在同一制导周期上，利用上一周期的全弹道信息及当前状态量进行制导指令的在线重构；（4）在下一制导周期上，重复步骤（2）和（3），直到接近目标为止；（5）接近目标后，应用比例导引方法击中目标。能够对速度、攻角、攻击方向等非线性较强的状态量的终端值进行控制，能够很好的消除飞行过程中各种干扰引起的偏差，并且能够很快的重构最优控制指令。

Description

一种多约束在线高斯伪谱末制导方法

技术领域

本发明涉及一种多约束在线高斯伪谱末制导方法，应用于高超声速飞行器的末制导过程中，属于飞行器的制导控制技术领域。

背景技术

经过近半个世纪的发展，常规飞行器的末制导技术已经趋于成熟，但随着高超声速飞行器的发展，尤其是全球打击任务的迫切需要，要求飞行器以高超声速的速度精确打击目标的同时还需保证终端的角度、速度等指标在一定的范围之内，这类指标的控制能够保证飞行器获得最佳的毁伤效果。

尽管现在有许多末段制导律能够保证飞行器精确的打击目标，但却很少有制导律能够导引飞行器以特定的角度攻击目标，能够控制末段攻角的制导律则更少。1973年，M.Kim和Kelly Grider利用线化后的Riccati方程求解线性二次问题获得了保证终端飞行弹道倾角的次优制导律。之后，Kenney在八十年代末期做了大量关于求解线性二次问题获得末段制导律的研究。1988年，Rusank提出了一种新的求解Riccati方程问题的方法，该方法成为一种处理最优末制导问题的有效方法。1995年，Idan获得了具有初始和终端弹道倾角约束的末制导问题的解析解。之后，Kim提出了一种保证终端弹道倾角的偏置比例导引。2008年，LuPing基于控制比例导引方程，推导了一种保证终端角度的自适应末段制导律，该方法通过逻辑关系将不同制导平面的弹道导引到能够保证飞行器终端角度收敛的状态上去。上面的工作都是基于末端飞行角度的制导律，并没有包括末段的攻角控制。2010年，邢强通过分段的策略推导了既能保证终端飞行角度又能保证末段攻角的分段比例导引。

上面提到的方法都是基于线性模型或部分非线性方程解析解推导得到的制导律。由于速度具有较强的非线性和多耦合的特点，对于飞行器的末端速度控制一直都是制导领域中的一个难点。

近年来，随着数值解算非线性规划问题方法的发展，能够很快的求解一个光滑的非线性最优控制问题。1995年，Elnager将求解流体偏微分方程问题的勒让德伪谱法引入到求解最优控制问题之中。之后，I.Ross提出了关于非线性最优控制问题的一阶条件与非线性规划问题的KKT条件之间等价的协态映射理论，揭示了伪谱法求解非线性最优控制问题的高精度性质。2006年，DavidBeson在勒让德伪谱法的基础上，引入高斯节点对非线性最优控制问题进行离散，提出了高斯伪谱法，并发现该方法的非线性最优控制问题的一阶条件与非线性规划问题的KKT条件之间具有更好的映射完整性，因此具有更高的效率。与此同时，采用高斯节点也能够保证高斯型积分公式具有更好的精度。高效率的算法为在线应用铺平了道路，2006年，I.Ross研究了勒让德伪谱法在线应用的收敛性，只要解算足够快，就能保证在线应用的有效收敛。

本发明受高斯伪谱法以及勒让德伪谱在线应用研究的启发，提出了一种多约束在线高斯伪谱末制导方法。该方法结合了高斯伪谱法、基于时钟信号的指令保持技术以及初值猜测策略，在一个制导周期内采用基于时钟信号的指令保持将最优控制指令带入实际飞行过程，并采用上一步的求解结果作为初始猜测，通过高斯伪谱法求解当前制导周期下的最优控制问题。本发明的优势在于对飞行器末制导具有很高的计算效率，只需较少的离散点就能获得较高的精度，能够消除飞行过程中存在的各种干扰因素（模型误差、风干扰、大气环境偏差、气动偏差等）引起的偏差，能够对终端状态（速度、飞行角度、攻角等）进行控制。本发明通过蒙塔卡罗打靶试验进行测试，仿真结果表明该方法能够很好的消除过程干扰引起的偏差，并对飞行器的终端状态进行有效的控制。

发明内容

本发明是一种多约束在线高斯伪谱末制导方法，该方法利用飞行器末段弹道的光滑性质，采用高斯伪谱法对全状态的非线性动力学方程形成的最优控制问题进行求解，获得快速的弹道解算。该方法的实施流程如图1所示，包括：

（1）采用高斯伪谱法对非线性状态方程进行末端轨迹优化；

以初始状态量为初始条件，考虑下面末段飞行器的动力学方程和性能指标所构成的最优控制问题。

$\stackrel{\·}{x}=f(x,u,t)---\left(1\right)$

$J[x(\·),u(\·)]={(x\left(t_f\right)-x_f)}^{T}K(x\left(t_f\right)-x_f)+\underset{t_0}{\overset{t_f}{\∫}}{u}_1^2+u_2^2\mathrm{dt}---\left(2\right)$

采用高斯伪谱法使用下面的公式可以将上述最优控制问题转化为时间区域在τ∈[-1,1]的问题。

$t=\frac{t_f-t_0}{2}\mathrm{\τ}+\frac{t_f-t_0}{2}---\left(3\right)$

高斯伪谱法和勒让德伪谱法、切比雪夫伪谱法一样，都是采用插值多项式来逼近状态变量和控制变量。在LG节点的每个Lagrange多项式的微分方程可以用一个微分矩阵来逼近，这个矩阵D∈R^N×(N+1)，其元素可以表示为

$D_{\mathrm{kl}}={\stackrel{\·}{L}}_i\left({\mathrm{\τ}}_k\right)=\underset{l=0}{\overset{N}{\mathrm{\Σ}}}\frac{\underset{j=0,j\≠i,l}{\overset{N}{\mathrm{\Π}}}({\mathrm{\τ}}_k-{\mathrm{\τ}}_j)}{\underset{j=0,j\≠i}{\overset{N}{\mathrm{\Π}}}({\mathrm{\τ}}_k-{\mathrm{\τ}}_j)}---\left(4\right)$

其中k＝1,...,N且i＝0,...,N。这样状态微分方程就可以通过微分矩阵转换为代数约束，其约束表示为

$\underset{l=0}{\overset{N}{\mathrm{\Σ}}}{D}_{\mathrm{kl}}{X}_l-\frac{(t_f-t_0)}{2}f(X_k,U_k,\mathrm{\τ};t_0,t_f)=0---\left(5\right)$

连续的性能指标函数也可以通过积分公式表示为

$J=\mathrm{\Φ}(x(-1),x\left(1\right),t_0,t_f)+\frac{(t_f-t_0)}{2}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\ω}}_{k}g(X_k,U_k,\mathrm{\τ};t_0,t_f)---\left(6\right)$

其中的ω_k是高斯权函数。

这样就把连续的Bolza问题转换为非线性规划（NLP）问题，通过成熟的非线性规划求解器就可以找到该问题的最优解。

（2）在一个制导周期内，采用时钟保持机制将最优控制指令带入实际飞行过程；

在线高斯伪谱法是通过数值解算技术对基本的最优控制问题进行求解，获得离散的时变反馈控制，反馈的目的是消除模型的不确定性，外部干扰以及估计误差造成的影响。指令的保持是根据上一步求得的以时间为自变量的控制指令，采用图2的形式进行保持。

根据文献证明：如果下式

${\mathrm{\τ}}_i\≤\frac{W\left(r\right)}{\mathrm{Lip}{f}_x}---\left(7\right)$

成立，则有||x_R(t_i+1)-x_M(t_i+1)||≤ε，即能将误差控制在很小的范围之内，其中，W(r)是Lambert-W函数，r＝ε/δ。因此，成功运用反馈控制的关键是能够快速的得到一条开环控制，高斯伪谱法就是一种对光滑问题就有指数收敛性质的直接法。

（3）在同一制导周期上，利用上一周期的全弹道信息及当前状态量进行制导指令的在线重构；

当飞行器在当前制导周期上时，根据所受各种干扰的当前状态值，重新构造最优控制问题。假如没有干扰存在，上一周期的结果就为本次最优控制问题的最优解，假定干扰属于小量，则干扰存在情况下的最优解也在上一周期最优解的附近。因此，对于非线性规划问题，为了减少寻优算法的计算时间，利用上一周期制导计算产生的结果作为初值，估值策略如下

$t_i=\frac{t_f-t_0}{2}{\mathrm{\τ}}_i+\frac{t_f-t_0}{2}---\left(8\right)$

$L_i\left(t\right)=\underset{j=0,j\≠i}{\overset{N}{\mathrm{\Π}}}\frac{t-t_j}{t_i-t_j}---\left(9\right)$

$\left(t\right)\underset{i=0}{\overset{N}{\mathrm{\Σ}}}\left(t_i\right)L_i\left(t\right)---\left(10\right)$

其中τ_i∈[-1,1]，是LG节点，L_i(t)为拉格朗日插值多项式，X(t_i)为上次优化的状态变量和控制变量。

插值时间为

$t_j=\frac{t_f-t_0-t_c}{2}{\mathrm{\τ}}_j+\frac{t_f-t_0-t_c}{2}---\left(11\right)$

插值结果为

$X\left(t_j\right)=\underset{j=0}{\overset{N}{\mathrm{\Σ}}}X\left(t_i\right)L_i\left(t_j\right)---\left(12\right)$

其中，t_c为制导周期，X(t_j)作为初值带入到非线性规划问题求解器进行计算，由仿真结果可以看出计算效率很高。

（4）在下一制导周期上，重复步骤（2）和（3），直到接近目标为止；

通过步骤（2）和（3）可以实现飞行器制导指令的在线重构，每次修正都能消除该飞行器之前的所有干扰影响，实现实时情况下的最优制导，该方法对于非线性较强的终端状态量（如速度、攻角、飞行角度）具有很好的控制效果。

（5）接近目标后，应用比例导引方法击中目标。

在飞行的最后时刻，为了提高打击精度，可以将制导律切换为比例导引方法。速度、飞行角度等状态量的控制都属于长周期控制过程，在飞行的前段时间可以通过（1）～（4）流程实现对该类状态量的控制，在飞行的后段时间则可以通过比例导引方法提高打击精度。

本发明的优点在于：

（1）采用全状态的非线性动力学方程，不进行线性假设，对于飞行动力学模型具有很好的数学描述；

（2）能够对终端状态如速度、攻角、飞行角度等非线性较强的状态量进行控制；

（3）能够很好的消除飞行过程中各种干扰引起的偏差，并能够很快的重构最优控制。

附图说明

使用以下说明和附图会更好地理解本申请，其中：

图1为实施流程图；

图2是指令保持示意图；

图3是空间弹道；

图4是速度曲线；

图5是弹道倾角曲线；

图6是攻角曲线；

图7是飞行时间；

图8是X方向误差概率分布情况；

图9是Y方向误差分布情况；

图10是速度误差分布情况；

图11是弹道倾角误差分布情况；

图12是攻角误差分布情况。

具体实施方式

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

将本发明方法应用于高超声速飞行器的末段制导，保证飞行器以确定的弹道倾角攻击目标，并且在打击时刻的速度和攻角都在一定的范围之内。建立飞行器在飞行过程中的扰动模型，并考虑末段关机影响进行蒙特卡洛打靶试验。所有工作都在3.3Ghz处理器的台式电脑上使用Matlab2008b版本实现，更为高效的编译方式还能提高计算效率。

1动力学方程模型

1.1动力学方程

考虑大气相对于地球相对静止的无动力三自由度质点弹道运动学方程如下

其中，x,y,z为位置坐标，V为速度，γ为弹道倾角，为航向角，α为攻角，σ为侧倾角，g＝μ/r²是重力加速度，μ为地球重力加速度常数，升力L和阻力D为

$L=\frac{q{S}_{\mathrm{ref}}\mathrm{Cl}}{m};D=\frac{{\mathrm{qS}}_{\mathrm{ref}}\mathrm{Cd}}{m}---\left(14\right)$

其中，是动压，Cl是升力系数，Cd是阻力系数，S_ref是参考面积，m为质量，ρ为大气密度，大气密度的通用公式为

ρ＝ρ₀exp(-βH)      (15)

其中，ρ₀是海平面大气密度，H为高度，β为密度常数。

1.2约束条件

从结构和热防护的角度出发，要求高超声速飞行器在飞行过程中满足动压约束以及过载约束，其约束的表达式为

$q=\frac{1}{2}\mathrm{\ρ}{V}^2\≤q_{\max };\mathrm{nz}=\frac{L}{\mathrm{mg}}\≤{\mathrm{nz}}_{\max }---\left(16\right)$

其中q_max和nz_max可以根据具体飞行器选取。

考虑到高超声速飞行器在末段飞行过程中的恶劣环境，要求控制攻角α，侧倾角σ保持在一定的范围之内，并且其变化率也必须保持在一定的变化范围之内

$\left\{\begin{array}{c}\left|\mathrm{\α}\right|\≤{\mathrm{\α}}_{\max };\left|u_1\right|\≤u_{1\max }\\ \left|\mathrm{\σ}\right|\≤{\mathrm{\σ}}_{\max };\left|u_2\right|\≤u_{2\max }\end{array}\right.---\left(17\right)$

1.3优化性能指标

为使优化得到的控制变量足够光滑，选取各个控制变量的平方和最小这一性能指标作为优化性能指标。在线反馈的过程中，为了增加寻优算法的收敛区域，增强制导算法的鲁棒性，可以把末端状态变量的约束加入到性能指标中，具体形式为

$J[x(\·),u(\·)]={(x\left(t_f\right)-x_f)}^{T}K(x\left(t_f\right)-x_f)+\underset{t_0}{\overset{t_f}{\∫}}(u_1^2+u_2^2)\mathrm{dt}---\left(18\right)$

其中，K为加权系数，根据变量数量级选取。

对于某些特殊的擦边界最优控制问题，由于原问题在擦边界的过程中，不再是光滑的最优控制问题，伪谱法对其进行解算不会具有指数收敛性质，将失去解算的快速性，这时通过改变性能指标的方式可以将该问题改变为光滑问题，进而具有指数收敛的性质以进行在线优化。

2仿真结果及分析

2.1飞行器模型

本文采用的飞行器模型为美国在“猎鹰”计划支持下，发展的新一代以高超声速为载体的远程快速精确打击武器系统——通用航空飞行器（CAV），该飞行器具有乘波构型，依靠气动力在大气层内滑行，具有升组比大，飞行马赫数高，空间环境恶劣，对控制系统要求较高等特点，根据NASA公布的公开资料，CAV的总体参数如表1所示。

表1CAV总体参数

低层稠密大气气动拟合公式如下

$\left\{\begin{array}{c}\mathrm{Cl}={\mathrm{Cl}}_1\mathrm{\α}+{\mathrm{Cl}}_2{e}^{{\mathrm{Cl}}_3\mathrm{Ma}}+{\mathrm{Cl}}_0\\ \mathrm{Cd}={\mathrm{Cd}}_1{\mathrm{\α}}^2+{\mathrm{Cd}}_2{e}^{{\mathrm{Cd}}_3\mathrm{Ma}}+{\mathrm{Cd}}_0\end{array}\right.---\left(19\right)$

其中气动数据参照文献。

2,2任务要求

飞行器的末端任务是在初始拉偏和过程拉偏情况下满足以下任务要求

攻角、侧倾角变化范围分别为

-15°≤α≤25°;-270°≤σ≤-90°      (21)

控制变量约束为

任务要求的主要难点在于速度控制，飞行器的速度对升力阻力非常敏感，也是长周期的控制过程，调节能力非常有限，要求飞行器在末端能量管理段进行精确的控制。

2.3标称情况下的制导结果

考查飞行器对于不同的初始偏差所形成的标准弹道，从而验证制导律对于不同工况作用下的适应性。表2为考虑5种工况情况的初始状态，初速均为1718m/s，表3为终端状态要求。

表2初始状态

表3终端状态

图3为空间弹道曲线，可以看出飞行器在弹道形式上能够满足多种工况的需要，并且以相似的弹道形式击中目标。

图4和图5分别为不同工况下的速度曲线和弹道倾角曲线，不同的工况由于末端状态相同，在飞行的后期差异不大，对于飞行弹道的调整主要集中在飞行的前期，此时的速度和弹道倾角差异较大，都能满足末速3马赫，末弹道倾角-80度的要求。

图6为不同工况下的攻角变化曲线，攻角在飞行后期差异不大，在飞行的前期差异很大，飞行器对于弹道的调节也主要是集中在飞行的前能量管理段，攻角在末端都能满足0度的要求。

图7为在线制导的计算时间，五种工况均在主频为3.3GHz的计算机下使用Matlab7.0进行仿真，时间均集中在0.1秒左右，由于不存在过程干扰，优化时间主要集中在克服LG阶段之间的龙格现象。同样表明，0.8s的制导周期完全满足指令解算的需要。

表4为不同工况下的飞行标准轨迹制导结果，5种工况的在脱靶量（Miss）和速度都能保证在0.01的数量级上，命中角度能保证在1E-4数量级上。由此可知，角度更容易得到保证。

表格4末端状态

2.4扰动模型

建立飞行器的温度，密度，风干扰随高度变化的拉偏模型，其具体表达形式为

其中，W(H),P(H),T(H)为风速、密度、温度随高度变化的函数。W(H)最大能达到50m/s，P(H)最大能够达到9%，T(H)最大能够达到30度。

飞行器的气动拉偏为标准气动升阻力的基础上加以偏差值，其具体表达形式为

由于飞行器需要进行速度控制，到达窗口很小，这需要飞行器在巡航段准确的进行关机控制。假设飞行器在巡航段匀速的均匀飞行，考虑从初始状态开始利用高斯伪谱法寻找最优关机时刻，并以此作为初值带入末端反馈制导中。

2.5过程扰动存在下的蒙特卡洛打靶结果

为了测试制导律的鲁棒性，下面对考虑关机情况下的飞行器进行1000次蒙特卡罗打靶实验，拉偏分布情况如表5所示。

表5蒙特卡洛打靶拉偏分布

图8和图9为位置方向的概率分布，各个方向的误差都在0.1m以内，落点概率偏向负方向。图10和图11为速度和弹道倾角的概率分布，速度和弹道倾角大部分满足所需要的范围，并且速度分布偏低，弹道倾角分布较为均匀。

图12为末端攻角的概率分布，方向由伪谱反馈对准下，攻角会很快的收敛到机动过载为0的攻角附近，从分布图上也能看出，主要集中在-0.7度附近，满足所需要求。

表6为末端状态的数理统计结果，飞行器在末端的脱靶量，速度，弹道倾角以及攻角的均值都很小，完全满足要求，从命中概率上看，脱靶量的落入概率为100%，速度的落入概率相对较低，为81.2%，弹道倾角和攻角分别为93.3%和99.7%，也能反映出速度项最难保证。蒙特卡罗打靶结果也表明伪谱反馈制导能够在恶劣的飞行环境下具有很好的制导能力和鲁棒性。

表6末端状态统计结果

Claims (1)

1.一种多约束在线高斯伪谱末制导方法，该方法包括如下步骤：

（1）采用高斯伪谱法对非线性状态方程进行末端轨迹优化；

（2）在一个制导周期内，采用时钟保持机制将最优控制指令带入实际飞行过程；

（3）在同一制导周期上，利用上一周期的全弹道信息及当前状态量进行制导指令的在线重构；

（4）在下一制导周期上，重复步骤（2）和（3），直到接近目标为止；

（5）接近目标后，应用比例导引方法击中目标。