含滑行時(shí)間約束的真空段彈道設(shè)計(jì)研究

胡冬生，劉 楠，童科偉，李 爍，張 皓

(1. 中國(guó)運(yùn)載火箭技術(shù)研究院研究發(fā)展中心，北京 100076；2. 中國(guó)科學(xué)院空間應(yīng)用工程與技術(shù)中心，北京 100094)

0 引言

在進(jìn)行運(yùn)載火箭彈道設(shè)計(jì)時(shí)，為了提升火箭執(zhí)行太陽同步軌道或近地軌道發(fā)射任務(wù)的性能，通常需要末級(jí)具備兩次點(diǎn)火能力，設(shè)置無動(dòng)力滑行段，以減小速度損失，如我國(guó)的長(zhǎng)征八號(hào)、長(zhǎng)征四號(hào)丙、長(zhǎng)征六號(hào)火箭等。未來大量低軌星座任務(wù)的發(fā)射也將依賴于這類彈道設(shè)計(jì)。但受到動(dòng)力系統(tǒng)工程研制、推進(jìn)劑管理及測(cè)控資源等方面因素的制約，末級(jí)滑行時(shí)間往往受限，無法通過霍曼變軌方式來最大程度發(fā)揮運(yùn)載火箭的能力[1]。這就帶來了含滑行時(shí)間約束的末級(jí)多次點(diǎn)火彈道設(shè)計(jì)問題，這類問題一般需要火箭滑行時(shí)間取約束區(qū)間內(nèi)的最大極限值，以更好地利用重力轉(zhuǎn)彎，因而問題轉(zhuǎn)化為固定滑行時(shí)長(zhǎng)的彈道設(shè)計(jì)。

國(guó)內(nèi)在進(jìn)行帶滑行段飛行任務(wù)的彈道設(shè)計(jì)時(shí)，通常采用參數(shù)化的工程方法，通過對(duì)程序角、關(guān)機(jī)時(shí)間的迭代和優(yōu)化來獲得滿足入軌約束的最大運(yùn)載能力[2]。彈道分段增多，設(shè)計(jì)變量個(gè)數(shù)相應(yīng)增加，使得彈道設(shè)計(jì)復(fù)雜，影響了設(shè)計(jì)效率。在將最優(yōu)控制理論應(yīng)用于含滑行段的彈道優(yōu)化方面，崔乃剛等[3]研究了大氣層外二三級(jí)間有滑行段的最優(yōu)軌跡生成方法，需要對(duì)滑行時(shí)間進(jìn)行在線優(yōu)化，無法解決本文提出的問題。

國(guó)外Gath等[4]、Lu等[5]、Jezewski[6]主要基于主矢量理論對(duì)多次點(diǎn)火飛行器的真空段軌跡和制導(dǎo)設(shè)計(jì)進(jìn)行了研究，但均未對(duì)滑行時(shí)間進(jìn)行約束，而是需要優(yōu)化滑行時(shí)間，不適用于我國(guó)運(yùn)載火箭設(shè)計(jì)的具體情況。Pan等[7]進(jìn)一步研究推導(dǎo)了含有滑行時(shí)間約束的切換條件(Switching Condition)，但未給出應(yīng)用于火箭末級(jí)多次點(diǎn)火彈道設(shè)計(jì)的方法和說明。

本文在前人研究的基礎(chǔ)上，將含固定滑行時(shí)長(zhǎng)的切換條件引入真空段含滑行時(shí)間約束的末級(jí)兩次點(diǎn)火彈道設(shè)計(jì)問題中，將彈道優(yōu)化轉(zhuǎn)換成對(duì)兩點(diǎn)邊值問題的求解和對(duì)運(yùn)動(dòng)方程的積分，同時(shí)對(duì)比分析了兩點(diǎn)邊值問題迭代求解過程中幾種滑行段彈道預(yù)報(bào)方法的特點(diǎn)和效果，以及該方法與傳統(tǒng)彈道設(shè)計(jì)方法的優(yōu)劣，以期為該類彈道設(shè)計(jì)問題提供新的思路和方法。

1 含滑行時(shí)間約束的彈道設(shè)計(jì)建模

1.1 火箭末級(jí)飛行過程

火箭末級(jí)一般在大氣層外飛行，通過在兩次工作之間加入滑行段來實(shí)現(xiàn)重力轉(zhuǎn)彎，改變火箭飛行方向，從而節(jié)省推進(jìn)劑。對(duì)于兩次點(diǎn)火的火箭來說，要先后經(jīng)歷末級(jí)一次動(dòng)力飛行段、滑行段、末級(jí)二次動(dòng)力飛行段，最終實(shí)現(xiàn)入軌，飛行過程見圖1。其中τl1表示末級(jí)一次點(diǎn)火時(shí)間，由大氣層內(nèi)飛行時(shí)間或前一級(jí)飛行時(shí)間確定；τkl1表示末級(jí)一次關(guān)機(jī)時(shí)間，為設(shè)計(jì)變量；τl2表示末級(jí)二次點(diǎn)火時(shí)間，當(dāng)滑行時(shí)間固定時(shí)，由τkl1確定；τkl2表示末級(jí)二次關(guān)機(jī)時(shí)間，為設(shè)計(jì)變量。

圖1 火箭末級(jí)飛行過程Fig.1 Flight process of last stage

1.2 運(yùn)動(dòng)方程

火箭在真空段飛行時(shí)，所受的力主要是地球引力和推力，其在發(fā)射慣性坐標(biāo)系的運(yùn)動(dòng)方程為

(1)

式中，r為火箭相對(duì)地心的位置矢量，V為速度矢量；g(r)為引力加速度；T為常值推力；Ib為發(fā)動(dòng)機(jī)推力方向的單位矢量，假定推力方向與火箭縱軸方向重合；m(t)為火箭當(dāng)前質(zhì)量。

采用線性引力場(chǎng)假設(shè)，對(duì)運(yùn)動(dòng)方程進(jìn)行歸一化處理[8]，得

(2)

1.3 邊界條件及性能指標(biāo)

(3)

真空段最優(yōu)飛行彈道為推進(jìn)劑消耗最少的飛行軌跡，在常值推力和滑行時(shí)間固定的情況下，即要求飛行時(shí)間最短，其性能指標(biāo)設(shè)置如下

minJ=τf

(4)

2 含滑行時(shí)間約束的彈道優(yōu)化思路

2.1 最優(yōu)化條件

應(yīng)用最優(yōu)控制理論，哈密頓函數(shù)為

(5)

式中，λr，λV為協(xié)態(tài)變量，ν為拉格朗日乘子。

根據(jù)龐特里亞金極小值原理，最優(yōu)縱軸方向滿足

(6)

協(xié)態(tài)變量滿足

(7)

協(xié)態(tài)變量和狀態(tài)變量的解析解可表示為[9]

(8)

(9)

式中，Ω(τ)，Γ(τ)分別為狀態(tài)轉(zhuǎn)移矩陣、轉(zhuǎn)換矩陣，Ic(τ)，Is(τ)為推力積分，表達(dá)式見文獻(xiàn)[9]。當(dāng)火箭處于滑行段時(shí)，僅受地球引力作用，Ic(τ)≡0，Is(τ)≡0，可以用二體運(yùn)動(dòng)解等彈道預(yù)報(bào)方法來代替式(9)，以獲得更為精確的狀態(tài)變量值。

由橫截條件和哈密頓函數(shù)終值條件，最終可得[9]

(10)

(11)

對(duì)于滑行時(shí)間固定的火箭末級(jí)來說，公式(12)的切換條件成立[7]

|λV(τl2)|-|λV(τkl1)|=0

(12)

2.2 程序角計(jì)算

在彈道設(shè)計(jì)中，一般設(shè)計(jì)偏航程序角ψcx=0，軌道傾角約束通過對(duì)發(fā)射方位角的迭代來滿足。由此可推導(dǎo)出λr(τ)和λV(τ)的第3個(gè)分量均為0，俯仰程序角公式變?yōu)?/p>

(13)

2.3 最優(yōu)彈道求解方法

當(dāng)火箭處于滑行段時(shí)，滑行段結(jié)束時(shí)刻的位置、速度僅與一次關(guān)機(jī)點(diǎn)的位置、速度相關(guān)，因此不進(jìn)行迭代，僅進(jìn)行固定滑行時(shí)間內(nèi)的彈道積分。

當(dāng)火箭處于末級(jí)二次動(dòng)力飛行段時(shí)，通過迭代λr0=[λr0xλr0y0]T，λV0=[λV0xλV0y0]T和τkl25個(gè)變量來求解滿足式(3)、(10)、(11)中5個(gè)終端約束和最優(yōu)化條件的當(dāng)前協(xié)態(tài)變量、最優(yōu)程序角及關(guān)機(jī)時(shí)間，再代入原始的運(yùn)動(dòng)方程進(jìn)行彈道積分，進(jìn)而計(jì)算出滿足入軌約束要求和最優(yōu)性的末級(jí)二次動(dòng)力飛行段彈道。該飛行段彈道積分以半長(zhǎng)軸約束為終止條件。至此完成了真空段彈道的設(shè)計(jì)計(jì)算。

具體設(shè)計(jì)流程見圖2。圖中，a*表示軌道半長(zhǎng)軸約束，a表示當(dāng)前時(shí)刻對(duì)應(yīng)的半長(zhǎng)軸；ε1和ε2為小量，根據(jù)設(shè)計(jì)精度要求進(jìn)行設(shè)定。在各飛行段的迭代計(jì)算中，本文采用Levenberg-Marquardt算法進(jìn)行求解，具有對(duì)迭代初值不敏感、求解速度快等優(yōu)點(diǎn)。

圖2 含滑行時(shí)間約束的真空段彈道設(shè)計(jì)流程Fig.2 Design process of vacuum trajectory with coasting time constraint

3 滑行段彈道預(yù)報(bào)方法的影響

由第2章中兩點(diǎn)邊值問題的迭代求解過程可以看出，在末級(jí)一次動(dòng)力飛行段針對(duì)每組迭代初值的快速解算過程均需要用到滑行段彈道解析解，本章就幾種預(yù)報(bào)方法對(duì)滑行段彈道解算和彈道設(shè)計(jì)的影響進(jìn)行分析。除了式(9)所示的解析解外，另外兩種常用的滑行段運(yùn)動(dòng)解析解分別為基于二體運(yùn)動(dòng)的Goodyear方法[10]和基于中間軌道理論、含J2項(xiàng)攝動(dòng)的Vinti方法[11-12]。

3.1 滑行段不同彈道預(yù)報(bào)方法的精度

若滑行段開始時(shí)刻火箭在發(fā)射慣性系的位置為[73 997.5，6 538 780.4，1 303.3]m，速度為[4 080.187，580.715，59.195]m·s-1，則分別利用式(9)解析解、Goodyear方法、含J2項(xiàng)Vinti方法和含J2項(xiàng)數(shù)值積分計(jì)算出滑行600 s結(jié)束時(shí)的位置、速度，如表1所示。

表1 幾種彈道預(yù)報(bào)方法的計(jì)算結(jié)果

由表1可見，式(9)解析解與數(shù)值積分的偏差較大，Goodyear結(jié)果偏差較小，Vinti結(jié)果最為接近數(shù)值積分結(jié)果，位置偏差為10 m量級(jí)，速度偏差小于0.01 m·s-1。在實(shí)際彈道設(shè)計(jì)中，可以考慮采用Goodyear或Vinti方法來進(jìn)行滑行段彈道解算。

3.2 滑行段不同彈道預(yù)報(bào)方法對(duì)彈道設(shè)計(jì)的影響

以文獻(xiàn)[3]中某火箭相關(guān)參數(shù)為算例，輸入?yún)?shù)見表2。末級(jí)滑行時(shí)間為600 s，目標(biāo)軌道為500 km圓軌道。

表2 真空段彈道設(shè)計(jì)輸入?yún)?shù)

分別采用基于Goodyear預(yù)報(bào)的彈道設(shè)計(jì)方法、基于Vinti預(yù)報(bào)的彈道設(shè)計(jì)方法開展仿真分析，并與傳統(tǒng)彈道設(shè)計(jì)結(jié)果進(jìn)行對(duì)比。計(jì)算結(jié)果見表3及圖3～5。

表3 不同彈道設(shè)計(jì)方法的計(jì)算結(jié)果

圖3 火箭末級(jí)飛行地心距Fig.3 Geocentric distance of last stage flight

圖4 火箭末級(jí)飛行速度Fig.4 Velocity of last stage flight

圖5 火箭末級(jí)飛行俯仰程序角Fig.5 Pitch program angle of last stage flight

在該算例中，傳統(tǒng)彈道設(shè)計(jì)方法涉及點(diǎn)火時(shí)刻程序角、動(dòng)力飛行段程序角斜率、關(guān)機(jī)時(shí)間等6個(gè)設(shè)計(jì)變量，約束參數(shù)僅有入軌點(diǎn)地心距、速度和彈道傾角3個(gè)，遠(yuǎn)少于設(shè)計(jì)變量個(gè)數(shù)，因此需要進(jìn)行大量、充分的優(yōu)化運(yùn)算；而本文提出的方法僅需針對(duì)兩點(diǎn)邊值問題進(jìn)行不斷的迭代計(jì)算和彈道積分，思路相對(duì)簡(jiǎn)潔。仿真結(jié)果表明，3種方法均可以得到滿足入軌約束的飛行彈道，關(guān)機(jī)時(shí)間、入軌質(zhì)量較為接近，其中基于Vinti的彈道設(shè)計(jì)結(jié)果與傳統(tǒng)方法設(shè)計(jì)結(jié)果更為吻合，包括一次、二次關(guān)機(jī)時(shí)間等，最終的入軌質(zhì)量?jī)H相差8.7 kg。而基于Goodyear的彈道設(shè)計(jì)方法采用了線性引力場(chǎng)假設(shè)和滑行段二體運(yùn)動(dòng)解析解，迭代模型精度有所損失，進(jìn)而導(dǎo)致了在偏心率的入軌約束上精度稍差，入軌質(zhì)量較其他2種方法小大約180 kg，二次動(dòng)力飛行段程序角斜率也更大，達(dá)到約8(°)/s， 將給火箭姿態(tài)控制帶來不利影響。

4 多級(jí)火箭含滑行時(shí)間約束的彈道設(shè)計(jì)

4.1 多級(jí)火箭含滑行時(shí)間約束的彈道設(shè)計(jì)思路

對(duì)于末級(jí)可以兩次啟動(dòng)的三級(jí)火箭而言，通常二級(jí)和三級(jí)均在大氣層外飛行，本文提出的彈道設(shè)計(jì)方法仍然適用，只需在2.3節(jié)的設(shè)計(jì)流程中增加二級(jí)飛行段。由于二級(jí)飛行時(shí)間固定，因此迭代參數(shù)不變，通過設(shè)置迭代初值λr0=[λr0xλr0y0]T，λV0=[λV0xλV0y0]T,τkl1，τkl2和式(8)、(9)依次快速解算出二級(jí)飛行段、三級(jí)一次動(dòng)力飛行段、滑行段、三級(jí)二次動(dòng)力飛行段的狀態(tài)變量和協(xié)態(tài)變量，其余過程與2.3節(jié)相同。

本節(jié)以某常溫三級(jí)火箭為例進(jìn)行仿真計(jì)算[13]。

4.2 彈道設(shè)計(jì)算例仿真

火箭二級(jí)和三級(jí)相關(guān)參數(shù)取值見表4。三級(jí)滑行時(shí)間為550 s，目標(biāo)軌道為600 km太陽同步軌道。

分別采用基于Goodyear預(yù)報(bào)的彈道設(shè)計(jì)方法、基于Vinti預(yù)報(bào)的彈道設(shè)計(jì)方法、傳統(tǒng)彈道設(shè)計(jì)方法開展仿真分析。計(jì)算結(jié)果見表5及圖6～8。

表4 三級(jí)火箭真空段彈道設(shè)計(jì)輸入?yún)?shù)

表5 三級(jí)火箭不同彈道設(shè)計(jì)方法的計(jì)算結(jié)果

圖6 火箭二、三級(jí)飛行地心距Fig.6 Geocentric distance of second and third stage flight

圖7 火箭二、三級(jí)飛行速度Fig.7 Velocity of second and third stage flight

圖8 火箭二、三級(jí)飛行俯仰程序角Fig.8 Pitch program angle of second and third stage flight

在該算例中，采用傳統(tǒng)彈道設(shè)計(jì)方法時(shí)由于增加了二級(jí)飛行段，設(shè)計(jì)變量增加到8個(gè)，而約束參數(shù)仍然只有入軌點(diǎn)地心距、速度和彈道傾角3個(gè)，這就進(jìn)一步加大了優(yōu)化難度；而本文提出的方法相比3.2節(jié)算例并沒有增加迭代變量的個(gè)數(shù)，求解過程幾乎完全一致，僅增加了二級(jí)飛行段的迭代和彈道積分。

仿真結(jié)果表明，基于Vinti的彈道設(shè)計(jì)結(jié)果與傳統(tǒng)方法設(shè)計(jì)結(jié)果更為吻合，最終的入軌質(zhì)量?jī)H相差0.8 kg；而基于Goodyear的彈道設(shè)計(jì)結(jié)果在入軌質(zhì)量上較其他兩種方法減小近70 kg，二次動(dòng)力飛行段程序角斜率也更大，達(dá)到3.09(°)/s，與3.2節(jié)算例出現(xiàn)的情況一致。

5 兩種彈道設(shè)計(jì)思路對(duì)比

傳統(tǒng)彈道設(shè)計(jì)方法本質(zhì)上是將程序角、關(guān)機(jī)時(shí)間等設(shè)計(jì)變量代入完整、準(zhǔn)確的運(yùn)動(dòng)方程進(jìn)行數(shù)值積分，進(jìn)而開展迭代和優(yōu)化，該方法中的程序角采用了工程上常用的線性化處理，在程序角設(shè)計(jì)規(guī)律上有時(shí)會(huì)損失一定的最優(yōu)性。此外，彈道設(shè)計(jì)的最優(yōu)性主要靠?jī)?yōu)化算法來保證，優(yōu)化算法的魯棒性、尋優(yōu)能力都將影響最終的彈道設(shè)計(jì)結(jié)果。尤其是在含無動(dòng)力滑行段的情況下，設(shè)計(jì)變量個(gè)數(shù)遠(yuǎn)遠(yuǎn)多于約束個(gè)數(shù)，滑行段的存在也進(jìn)一步增強(qiáng)了問題的非線性，使得彈道優(yōu)化難度加大，對(duì)程序角初值的設(shè)置更為敏感，甚至陷入局部最優(yōu)解。

本文提出的基于線性引力場(chǎng)和滑行段切換條件的彈道設(shè)計(jì)方法相對(duì)簡(jiǎn)潔很多，其基本的迭代變量?jī)H有初始的4個(gè)協(xié)態(tài)初值和2個(gè)關(guān)機(jī)時(shí)間，彈道設(shè)計(jì)的最優(yōu)性完全依靠橫截條件、哈密頓函數(shù)終值條件和切換條件來保證。但同時(shí)，該方法在迭代過程中使用了簡(jiǎn)化的運(yùn)動(dòng)方程解析解，且最優(yōu)性條件也是基于線性引力場(chǎng)推導(dǎo)出來的，這也在一定程度上損失了最優(yōu)性。此外，通過仿真可以發(fā)現(xiàn)，滑行段彈道預(yù)報(bào)模型越精確，則迭代產(chǎn)生的程序角越接近最優(yōu)解，最終得到的入軌質(zhì)量就越大，也越容易滿足對(duì)入軌軌道的約束精度，同時(shí)也不會(huì)產(chǎn)生較大的程序角斜率。因此，在含滑行段的彈道設(shè)計(jì)中，應(yīng)當(dāng)優(yōu)先采用Vinti預(yù)報(bào)方法開展滑行段彈道解算。

通過對(duì)比和仿真可知，兩種彈道設(shè)計(jì)方法各有特點(diǎn)，通?？梢垣@得幾乎相同的設(shè)計(jì)結(jié)果。當(dāng)應(yīng)用于級(jí)數(shù)較多、含滑行段的火箭彈道設(shè)計(jì)問題時(shí)，基于線性引力場(chǎng)和切換條件的彈道設(shè)計(jì)方法能夠在一定程度上簡(jiǎn)化問題的求解，保證得到問題的近似最優(yōu)解，且在小推重比情況下更具有性能優(yōu)勢(shì)[8]；此外，使用的Vinti滑行段彈道預(yù)報(bào)方法也可以應(yīng)用在傳統(tǒng)彈道設(shè)計(jì)的滑行段計(jì)算中，從而避免滑行段彈道積分，大大提高彈道設(shè)計(jì)的效率。

6 結(jié)論

本文在線性引力場(chǎng)假設(shè)的基礎(chǔ)上，將含固定滑行時(shí)長(zhǎng)的切換條件引入彈道設(shè)計(jì)中，從而將含滑行時(shí)間約束的彈道優(yōu)化問題轉(zhuǎn)換成對(duì)兩點(diǎn)邊值問題的求解和對(duì)運(yùn)動(dòng)方程的積分。通過多個(gè)算例仿真得出以下結(jié)論：

1)基于線性引力場(chǎng)和滑行段切換條件的彈道設(shè)計(jì)方法思路簡(jiǎn)潔，能夠適用于對(duì)滑行時(shí)間有約束、末級(jí)兩次點(diǎn)火的兩級(jí)或多級(jí)火箭，設(shè)計(jì)結(jié)果與傳統(tǒng)設(shè)計(jì)方法幾乎一致；

2)滑行段彈道預(yù)報(bào)精度對(duì)最終的彈道優(yōu)化結(jié)果具有一定的影響，模型精度越高，優(yōu)化結(jié)果越好；

3)本文方法和思路應(yīng)用范圍較廣，能夠在各種推重比情況下獲得最優(yōu)彈道，在總體方案論證或初步設(shè)計(jì)階段可以替代傳統(tǒng)設(shè)計(jì)方法，以保證獲得更好的運(yùn)載性能。