基于迭代制導(dǎo)的運(yùn)載火箭入軌級(jí)彈道設(shè)計(jì)及仿真

高家智，陳小前

（1.國(guó)防科技大學(xué)空天科學(xué)學(xué)院，長(zhǎng)沙，410073；2.太原衛(wèi)星發(fā)射中心，太原，030027）

0 引 言

迭代制導(dǎo)制導(dǎo)精度高、任務(wù)適應(yīng)性強(qiáng)、箭上飛行軟件簡(jiǎn)單、對(duì)地面諸元準(zhǔn)備要求相對(duì)較低，在運(yùn)載火箭上應(yīng)用越來(lái)越廣泛[1～3]。文獻(xiàn)[4]針對(duì)主發(fā)動(dòng)機(jī)推力為常值的情況，通過(guò)優(yōu)化開關(guān)機(jī)時(shí)間來(lái)滿足終端位置約束，引入權(quán)重因子提高迭代求解的精度；文獻(xiàn)[5]針對(duì)多約束入軌條件，提出一種將迭代制導(dǎo)與數(shù)值積分相結(jié)合的軌跡預(yù)測(cè)制導(dǎo)方法，增加的推進(jìn)劑的消耗在可接受范圍之內(nèi)；文獻(xiàn)[6]在地心慣性系中建立航天器模型，以目標(biāo)軌道要素為終端約束條件，得到一種簡(jiǎn)單有效的控件變軌迭代制導(dǎo)算法；針對(duì)飛行過(guò)程中發(fā)動(dòng)機(jī)推力下降的重大故障，文獻(xiàn)[7]基于迭代制導(dǎo)算法對(duì)彈道進(jìn)行了重構(gòu)，在完成故障識(shí)別和能力預(yù)測(cè)后，利用迭代制導(dǎo)將有效載荷送入預(yù)先設(shè)計(jì)好的救援軌道上去。

Boeing 公司針對(duì)阿波羅計(jì)劃使用的土星5 運(yùn)載火箭迭代制導(dǎo)方法的不足，開發(fā)了線性正切制導(dǎo)律（Linear Tangent Guidance，LTG），并在飛馬火箭上升段入軌控制和離軌段返回控制中得到應(yīng)用[8]。本文方法利用了其閉路制導(dǎo)律和預(yù)估公式，但制導(dǎo)計(jì)算的思路與文獻(xiàn)提供的相反，引入了制導(dǎo)迭代變量，采用了根據(jù)初值求解兩點(diǎn)邊值問(wèn)題的Newton-Raphson 方法。最后利用數(shù)值仿真驗(yàn)證了本文所設(shè)計(jì)方法的有效性。

1 最優(yōu)控制問(wèn)題及其制導(dǎo)律

在不考慮稀薄大氣的情況下，火箭在真空段內(nèi)的運(yùn)動(dòng)方程可表示如下[9]：

終端約束條件為

式中V，R分別為火箭在地心發(fā)射慣性坐標(biāo)系中的速度和位置；G為重力矢量，G=-ω2R=-3R，其中，fM= 3.986 × 105km3s2為引力常數(shù)；ω為地球自轉(zhuǎn)角速度；為秒流量；a為加速度的模，a=，其中，P為發(fā)動(dòng)機(jī)推力，m為質(zhì)量；Ue為有效排氣速度；u為推力方向的單位向量。

定義哈密頓函數(shù)H和增廣函數(shù)F如下：

這里的λ，α和β是拉格朗日乘子。

定義如下變量：

式中φ，ψ分別為俯仰角和偏航角。根據(jù)最優(yōu)控制的必要條件有如下歐拉方程和橫截條件：

由歐拉方程可解出：

式中uR=unit(R)為火箭質(zhì)心地心矢徑的單位向量。

由橫截條件可得出：

相對(duì)于地球半徑而言，火箭因推進(jìn)引起高度的變化量是小量，并且λ與uR的夾角在90°附近，因此在式（8）中可忽略高階向量，得：

該方程的解為

式中t為飛行時(shí)間；λK為定常單位向量；為定常速率向量；K為定常時(shí)間，一般情況下K≈Tgo2，Tgo為剩余飛行時(shí)間。如果是常值加速度，則K=Tgo2。

可以證明，對(duì)于圓軌道入軌時(shí)，下式成立：

即兩矢量正交。

由于在迭代過(guò)程中，不斷刷新K、λK和，為簡(jiǎn)便起見(jiàn)，省去下標(biāo)K，即有：

考察式（13），存在如下特性：

b）如果采用地球平面假設(shè)，因?yàn)棣?0，則：

在本文的制導(dǎo)研究中，仍考慮地球是球體，且使用上述有關(guān)假設(shè)：

b）ω=λ˙；

c）忽略式（8）中的 3ω2(λ·uR)uR；

那么，制導(dǎo)律為

簡(jiǎn)化為L(zhǎng)TG 制導(dǎo)律：

2 迭代制導(dǎo)的實(shí)現(xiàn)

2.1 端點(diǎn)狀態(tài)預(yù)估

運(yùn)載火箭的迭代制導(dǎo)屬于復(fù)雜的兩點(diǎn)邊值問(wèn)題，關(guān)機(jī)時(shí)刻火箭位置和速度的預(yù)估精度關(guān)系到迭代制導(dǎo)的精度、穩(wěn)定性和魯棒性。數(shù)值積分是一種精確的解決辦法，但卻是無(wú)法在箭上實(shí)時(shí)實(shí)現(xiàn)的辦法，所以推導(dǎo)簡(jiǎn)化的解析代數(shù)公式就是一項(xiàng)根本任務(wù)。

在制導(dǎo)律式（15）的控制作用下，火箭在真空段的運(yùn)動(dòng)方程為

假設(shè)剩余飛行時(shí)間Tgo是已知的，且存在如下封閉形式的積分：

式中F1，F(xiàn)2，F(xiàn)3為待定系數(shù)。

積分式（18）有：

式中VD，RD分別為末端速度和位置矢量。

改寫式（19）、（20），得到：

再積分式（17），得到：

式中Vgo，Rgo分別為待增速度和位置矢量。

因此，如果已經(jīng)知道λ和，就可以根據(jù)積分公式LT、ST、QT得出Vgo和Rgo，從而在當(dāng)前位置R和速度V下，由制導(dǎo)律式（15）進(jìn)行制導(dǎo)的端點(diǎn)的預(yù)估公式為

式中VP，RP分別為預(yù)測(cè)的末端速度和位置矢量。

2.2 縱平面迭代制導(dǎo)

運(yùn)載火箭的歐拉角采用3-2-1 型表示的發(fā)慣系動(dòng)力學(xué)方程組，可知只有側(cè)平面的運(yùn)動(dòng)影響縱平面的運(yùn)動(dòng)，但縱平面的運(yùn)動(dòng)不會(huì)影響側(cè)平面的運(yùn)動(dòng)，特別當(dāng)偏航角很小時(shí)飛行控制可以解耦。因此可以先進(jìn)行縱平面的迭代控制，再進(jìn)行側(cè)平面的迭代控制。

縱平面內(nèi)的迭代控制采用標(biāo)準(zhǔn)入軌參數(shù)作為邊界條件。這里的邊界條件是采用極坐標(biāo)邊界條件而非直角指標(biāo)邊界條件，因?yàn)榘腴L(zhǎng)軸和偏心率與極坐標(biāo)直接相關(guān)，所以在縱平面內(nèi)的運(yùn)動(dòng)控制邊界條件是滿足標(biāo)準(zhǔn)入軌點(diǎn)對(duì)應(yīng)的地心矢徑、速度大小和當(dāng)?shù)厮俣葍A角，這樣就保證了橢圓軌道的半長(zhǎng)軸和偏心率，加上側(cè)平面的標(biāo)準(zhǔn)軌道跟蹤導(dǎo)引，控制了zK和，理論上也是控制軌道根數(shù)的5 個(gè)變量。升交點(diǎn)赤經(jīng)Ω主要由發(fā)射窗口的起飛時(shí)刻、火箭飛行時(shí)間和方位角決定，近地點(diǎn)幅角雖然是放開的，但跟蹤Z通道的標(biāo)準(zhǔn)軌道這種約束本身也對(duì)第6 個(gè)軌道根數(shù)進(jìn)行了隱式控制，這種方法一定程度約束了近地點(diǎn)幅角偏差的放大。

這里采用極坐標(biāo)邊界條件的另一個(gè)目的是提高迭代制導(dǎo)的魯棒穩(wěn)定性和收斂性，減小迭代制導(dǎo)結(jié)尾階段的過(guò)早發(fā)散效應(yīng)，三自由度制導(dǎo)仿真計(jì)算表明選擇這種邊界條件與視加速度反饋措施一起可以有效提高迭代制導(dǎo)的精度，使跳出迭代的剩余常值飛行時(shí)間小于0.5～1.0 s，幾乎迭代至最后，提高了制導(dǎo)精度。

縱平面內(nèi)的迭代控制采用 LTG 制導(dǎo)律方法，預(yù)測(cè)縱平面內(nèi)關(guān)機(jī)時(shí)刻發(fā)慣系內(nèi)的直角坐標(biāo)參數(shù)Xp、、Yp和，然后進(jìn)一步預(yù)測(cè)關(guān)機(jī)時(shí)刻的位置矢徑Rp、速度大小Vp和當(dāng)?shù)厮俣葍A角Θp，即：

迭代公式采用牛頓迭代法，即：

式中α,β為與俯仰角有關(guān)的制導(dǎo)變量。

因?yàn)椴捎昧藰O坐標(biāo)邊界條件，沒(méi)有顯式解，校正計(jì)算需要采用3 變量的牛頓快速迭代法，涉及3×3 雅可比逆矩陣的實(shí)時(shí)計(jì)算，3×3 雅可比逆矩陣不需要數(shù)值計(jì)算，可以有解析表達(dá)式，并且因?yàn)樯弦徊接?jì)算得到了的初值，下一步牛頓迭代計(jì)算只需一次就很精確，完全滿足制導(dǎo)收斂精度。

由迭代解α,β,Tgo可以計(jì)算LTG 制導(dǎo)律的乘子向量λ和得到：

發(fā)慣系縱平面推力方向的單位控制向量：

記u=[u(1),u(2)]T，故當(dāng)前時(shí)刻滿足邊界條件的俯仰姿態(tài)角應(yīng)為

2.3 側(cè)平面迭代制導(dǎo)

發(fā)射慣性系中的姿態(tài)角為俯仰角和偏航角為φ和ψ，初始質(zhì)量m0，秒流量，軸向推力F在發(fā)射慣性系中推力分量表示為

動(dòng)力學(xué)方程為

若偏航角為常值規(guī)律，火箭推力通常假設(shè)是常值，則推力的Z通道分量為常值，那么側(cè)向通道在關(guān)機(jī)點(diǎn)的狀態(tài)可以解析預(yù)測(cè)，即：

利用縱平面預(yù)測(cè)得到的縱向關(guān)機(jī)點(diǎn)狀態(tài)和關(guān)機(jī)時(shí)間以及上述側(cè)向通道關(guān)機(jī)點(diǎn)狀態(tài)公式預(yù)測(cè)，可以進(jìn)行軌道傾角i的解析預(yù)測(cè)計(jì)算，方法為：在縱向通道迭代計(jì)算完畢，預(yù)估縱向關(guān)機(jī)速度和位置，再根據(jù)上述常值偏航角和常值推力計(jì)算關(guān)機(jī)時(shí)的發(fā)慣系側(cè)向速度和側(cè)向位置，然后轉(zhuǎn)換到發(fā)射地慣系計(jì)算軌道傾角，則i=i(ψ)，根據(jù)軌道傾角期望值ic，容易單變量迭代求解偏航角指令。

那么對(duì)偏航角進(jìn)行迭代滿足所要求的軌道傾角屬于單變量迭代，很容易快速求解。

3 仿真結(jié)果及分析

仿真起始條件為：發(fā)慣系下起滑點(diǎn)及第3 級(jí)關(guān)機(jī)點(diǎn)參數(shù)X=593 125.71 m，Y=184 422.32 m，Z=-77 388 m，VX=6676.607 m/s，Vy=648.986 m/s，Vz=-286.533 m/s，標(biāo)準(zhǔn)點(diǎn)火時(shí)刻532.3 s，推力28 000 N，初始質(zhì)量1206.3 kg，秒耗量10 kg/s，標(biāo)準(zhǔn)滑行時(shí)間303.7 s[10]。

目標(biāo)衛(wèi)星軌道參數(shù)為：半長(zhǎng)軸a=6860.788 km，偏心率e=0，軌道傾角i=97.3 °。

制導(dǎo)精度要求：半長(zhǎng)軸偏差小于1 km，偏心率偏差小于0.0001，軌道傾角偏差小于0.01°。

經(jīng)過(guò)仿真計(jì)算，發(fā)慣系位移及速度變化情況見(jiàn)圖1所示。

圖1 發(fā)慣系位移及速度變化情況Fig.1 Position and Velocity Changes in Launch Inertial System

俯仰角、偏航角以及其它參數(shù)變化曲線如圖2 所示，滑行軌道、迭代制導(dǎo)軌道與運(yùn)行軌道的放大圖如圖3 所示。輸出結(jié)果如下：半長(zhǎng)軸6860.768 km，偏心率0.000 012，軌道傾角97.3072°，制導(dǎo)時(shí)間49.35 s，升交點(diǎn)角距163.42°?？梢?jiàn)長(zhǎng)半軸誤差20 m,偏心率偏差0.000 012，軌道傾角偏差0.0072°，均滿足制導(dǎo)精度要求。

圖2 主要飛行變量變化情況Fig.2 Major Flight Variable Variations

圖3 滑行軌道、迭代制導(dǎo)軌道及運(yùn)行軌道放大圖Fig.3 Enlarged View of Dragging Track,Iterated Guidance Track,and Orbit

4 結(jié)束語(yǔ)

本文對(duì)多級(jí)運(yùn)載火箭入軌段的迭代制導(dǎo)方法進(jìn)行了較為深入研究，從仿真結(jié)果來(lái)看本文所提出的入軌級(jí)迭代制導(dǎo)方法制導(dǎo)精度高，且計(jì)算簡(jiǎn)單、易于實(shí)現(xiàn)。