何兆偉，葛 冰，徐利杰，李平岐

(1. 北京宇航系統(tǒng)工程研究所，北京 100076；2. 中國長征火箭有限公司，北京 100070)

0 引言

智慧火箭研制的目標(biāo)是實(shí)現(xiàn)設(shè)計(jì)和試驗(yàn)效率的提升，具體表現(xiàn)為：1周完成總體構(gòu)型優(yōu)化、1個(gè)月完成總體方案設(shè)計(jì)且一次成功率100%[1]。而在火箭總體設(shè)計(jì)中，彈道設(shè)計(jì)起著極其重要的作用，火箭總體方案、設(shè)計(jì)參數(shù)、運(yùn)載性能與飛行方案都是根據(jù)彈道設(shè)計(jì)的結(jié)果而確定的[2]。運(yùn)載火箭彈道設(shè)計(jì)是一種典型的非線性規(guī)劃問題[3]，設(shè)計(jì)過程除了要符合火箭自身構(gòu)型特點(diǎn)和有效載荷入軌要求外，還要滿足發(fā)射場、落區(qū)、飛行安全、地面跟蹤測量與控制等約束條件。算法對初值的設(shè)置較為敏感，常出現(xiàn)運(yùn)算發(fā)散的情況[4]，特別是在發(fā)射任務(wù)類型變化時(shí)，現(xiàn)有工程設(shè)計(jì)方法還較為依賴人員的自身經(jīng)驗(yàn)，需要進(jìn)行部分手動調(diào)整和優(yōu)化。

本文從解決上述問題出發(fā)，對不同初值給定方法[5-8]進(jìn)行了大量仿真研究。最終通過引入多屬性權(quán)衡空間探索及全局優(yōu)化思想，結(jié)合彈道設(shè)計(jì)需求建立了權(quán)衡空間中方案的評價(jià)指標(biāo)，給出了初值搜索方法并據(jù)此重塑了工程彈道設(shè)計(jì)流程。同時(shí)，以較為復(fù)雜的三級半且末級兩次啟動的液體運(yùn)載火箭為對象進(jìn)行了仿真驗(yàn)證。

1 彈道設(shè)計(jì)模型

1.1 彈道計(jì)算

彈道計(jì)算采用旋轉(zhuǎn)橢球地球模型，考慮攝動項(xiàng)影響，所用動力學(xué)模型和飛行程序設(shè)計(jì)方法參見文獻(xiàn)[9]。

工程彈道計(jì)算中使用“彈道級”的概念以提升設(shè)計(jì)配置的通用性和適用范圍，每一個(gè)飛行段即為一個(gè)彈道級，包括結(jié)構(gòu)級(含質(zhì)量特性、加注諸元、動力參數(shù)等)、導(dǎo)引、關(guān)機(jī)、控制、氣動、飛行程序、分離和通用屬性等對象。每個(gè)彈道級按圖1所示流程進(jìn)行計(jì)算并積分，全部彈道級積分完成即完成了一條彈道的計(jì)算。

圖1 彈道級計(jì)算流程Fig.1 Calculation process of trajectory stage

1.2 迭代設(shè)計(jì)

彈道迭代設(shè)計(jì)的思路為通過對控制矢量X=[x1,x2,…,xn]的攝動來計(jì)算對入軌參數(shù)，即目標(biāo)矢量Y=[y1,y2,…,yn]的影響，通過Y與入軌要求的偏差解算X的增量，直至偏差符合精度要求[3,10]。針對半長軸、近地點(diǎn)高度、遠(yuǎn)地點(diǎn)高度、入軌點(diǎn)矢徑等長度目標(biāo)變量，通常精度要求取1×10-1m；針對軌道傾角、近地點(diǎn)幅角、升交點(diǎn)經(jīng)度、彈道傾角、分離時(shí)刻攻角等角度目標(biāo)變量，通常精度要求取1×10-4°。

數(shù)學(xué)模型為構(gòu)建雅克比矩陣

(1)

通過高斯消元法求解X的增量

(2)

設(shè)計(jì)流程如圖2所示。

圖2 彈道迭代設(shè)計(jì)流程Fig.2 Design process of trajectory

針對低軌圓軌道，Y的各分量一般選取半長軸a、入軌點(diǎn)地心矢徑r、入軌點(diǎn)速度v、入軌點(diǎn)彈道傾角θ等；針對高軌地球同步轉(zhuǎn)移軌道，Y的各分量一般選取近地點(diǎn)高度hp、遠(yuǎn)地點(diǎn)高度ha、近地點(diǎn)幅角ω、軌道傾角i等。而X的各分量一般在考慮火箭構(gòu)型特點(diǎn)、設(shè)計(jì)約束及對Y敏感度的基礎(chǔ)上，選取射向A0、不同彈道級的俯仰程序角φi、偏航程序角ψi、關(guān)機(jī)時(shí)間ti等參數(shù)。

X和Y的選取和對應(yīng)關(guān)系，以及其對設(shè)計(jì)收斂性的影響，也體現(xiàn)出了任務(wù)類型變化對火箭彈道設(shè)計(jì)過程的影響，以及人為經(jīng)驗(yàn)在設(shè)計(jì)過程中起到的作用。

2 初值確定方法

2.1 問題描述

如圖2所示，彈道設(shè)計(jì)流程最初要選定控制矢量X的分量并給定初值，參數(shù)選擇適當(dāng)則算法收斂，得到滿足入軌精度要求的彈道設(shè)計(jì)結(jié)果。該算法及流程的優(yōu)點(diǎn)在于可快速得到精確解，缺點(diǎn)在于對初值敏感，若初值距真值較遠(yuǎn)、不在收斂域內(nèi)，則計(jì)算發(fā)散。因此，面向試樣或應(yīng)用飛行階段的運(yùn)載火箭轉(zhuǎn)換常規(guī)任務(wù)類型、承接特殊軌道任務(wù)等工程實(shí)際場景，或開展新型運(yùn)載火箭構(gòu)型論證時(shí)，往往無法通過直接將目標(biāo)矢量Y更改至有效載荷入軌要求的方式自動得到設(shè)計(jì)結(jié)果，而需要人為手動調(diào)整，結(jié)合設(shè)計(jì)人員的經(jīng)驗(yàn)對初值進(jìn)行探索直至其進(jìn)入收斂域內(nèi)。方法的自動化、智能化程度低，設(shè)計(jì)效率有待提升。

2.2 方法與流程

結(jié)合權(quán)衡空間探索與全局優(yōu)化的思想確定彈道設(shè)計(jì)的初值。

根據(jù)多屬性權(quán)衡空間探索方法[11]，控制矢量X的不同取值可達(dá)到不同的效用和成本。這就構(gòu)成了一個(gè)權(quán)衡空間，代表了決策者所有的設(shè)計(jì)選擇，空間中每個(gè)點(diǎn)都是一個(gè)備選方案，探索權(quán)衡空間就是對備選方案進(jìn)行權(quán)衡評估尋找最優(yōu)解[12]。為兼顧局部收斂及全局搜索能力，采用改進(jìn)粒子群算法(PSO算法)[13]對權(quán)衡空間進(jìn)行全域搜索。即將每個(gè)控制矢量Xi(i=1,2,…,m)抽象為沒有質(zhì)量和體積的微粒，并延伸到n維空間，粒子Xi在n維空間里的位置表示為一個(gè)矢量[xi,1,xi,2,…,xi,n]，每個(gè)粒子還有一個(gè)速度矢量[vi,1,vi,2,…,vi,n]決定了它們的搜索方向和距離。所有粒子都有一個(gè)由被優(yōu)化函數(shù)決定的評價(jià)指標(biāo)值f，粒子知道自身到目前為止發(fā)現(xiàn)的最好位置(pbest)和現(xiàn)在的位置，即自身經(jīng)驗(yàn)；也知道整體粒子群體發(fā)現(xiàn)的最好位置(gbest)，即同伴經(jīng)驗(yàn)。粒子通過自身經(jīng)驗(yàn)和同伴經(jīng)驗(yàn)來決定下一步的搜索運(yùn)動，以達(dá)到全局探索和局部開采的平衡。

取評價(jià)指標(biāo)

(3)

權(quán)衡空間搜索的步驟如下：

1)隨機(jī)生成權(quán)衡空間中m個(gè)控制矢量Xi的位置和速度；

2)將當(dāng)前每個(gè)Xi的位置和評價(jià)值fi存在Xi的pbest中，將所有pbest中fi最小的Xi的位置和評價(jià)值存在gbest中；

3)用式(4)更新每個(gè)Xi的速度和位置；

vi,j(t+1)=wvi,j(t)+c1r1[pi,j-xi,j(t)]+
c2r2[pg,j-xi,j(t)]
xi,j(t+1)=xi,j(t)+vi,j(t+1),j=1,2,…,n

(4)

4)將每個(gè)Xi的評價(jià)值fi與其經(jīng)歷過的最小值比較，若較小則更新其pbest，比較當(dāng)前所有pbest和gbest，若較小則更新gbest；

5)若滿足停止條件(迭代輪次達(dá)到k，或gbest中的f值小于等于預(yù)設(shè)的截至條件f*)，則輸出結(jié)果，否則返回3)步。

將初值確定方法嵌入設(shè)計(jì)流程，重塑的工程彈道設(shè)計(jì)流程見圖3?？梢钥闯?，初值確定方法與現(xiàn)行成熟的工程設(shè)計(jì)算法有機(jī)銜接(利用現(xiàn)有流程中的彈道計(jì)算進(jìn)行評價(jià)值計(jì)算，得到收斂域內(nèi)的初值后交還給迭代設(shè)計(jì))，程序復(fù)雜度低，使得該方法完全符合工程實(shí)際應(yīng)用需要。

圖3 增加初值計(jì)算后的彈道設(shè)計(jì)流程Fig.3 Design process of trajectory with initial values calculation

3 仿真計(jì)算效果

3.1 仿真條件

以較為復(fù)雜的三級半、末級兩次啟動、常規(guī)發(fā)射任務(wù)為高軌GTO載荷的液體運(yùn)載火箭為對象，完成新低軌任務(wù)彈道設(shè)計(jì)，在滿足各項(xiàng)約束的同時(shí)，終端入軌條件如表1所示。

表1 終端條件

采用分解策略[4]對三級半運(yùn)載火箭進(jìn)行彈道設(shè)計(jì)，求解策略及控制矢量、目標(biāo)矢量的選取見表2。

表2 求解策略

3.2 仿真結(jié)果分析

火箭常規(guī)GTO任務(wù)控制矢量X的設(shè)計(jì)結(jié)果為[-15.124 1°, - 85.452 0°, 2.434 6°, 177.044 0 s]，以此作為初值，應(yīng)用工程現(xiàn)行設(shè)計(jì)流程(圖2)的收斂性如圖4所示?？梢钥闯?，目標(biāo)變量在迭代約10輪以后開始發(fā)散。

圖4 工程現(xiàn)行設(shè)計(jì)流程的收斂性曲線Fig.4 Convergence curve using original design process

應(yīng)用2.2節(jié)權(quán)衡空間搜索方法進(jìn)行初值計(jì)算，取粒子數(shù)m=30，迭代輪次k=500，評價(jià)值截止條件f*≤1×10- 8，權(quán)重調(diào)節(jié)系數(shù)ξj均取1，慣性權(quán)重w的計(jì)算中wmax=0.9和wmin=0.4，學(xué)習(xí)因子c1和c2均取1.5。

經(jīng)仿真計(jì)算，得到評價(jià)值f的變化曲線如圖5所示。最小值fmin=1×10- 6，對應(yīng)gbest中控制矢量X的位置為[19.697 1°, -97.719 3°, 25.983 0°, 44.996 7 s]，以此為初值進(jìn)行彈道迭代設(shè)計(jì)的收斂性如圖6所示。由圖6可以看出，迭代約30輪左右收斂，目標(biāo)矢量Y最終收斂至[500.000 0 km, 178.000 0°, 27.500 0°, 560.000 0 km]，滿足終端條件精度要求，X的最終設(shè)計(jì)結(jié)果(真初值)為[19.675 4°, -97.814 3°, 25.719 4°, 44.894 1 s]。由圖7可以看出，X逐步接近真初值的收斂過程。

圖5 評價(jià)值f曲線Fig.5 Curve of evaluation value

圖6 新流程下的目標(biāo)矢量收斂性曲線Fig.6 Convergence curve of target vector using new design process

圖7 新流程下的控制矢量收斂性曲線Fig.7 Convergence curve of control vector using new design process

4 結(jié)論

本文著重解決多約束運(yùn)載火箭彈道設(shè)計(jì)的初值敏感問題。在梳理了當(dāng)前工程彈道計(jì)算及迭代設(shè)計(jì)詳細(xì)流程的基礎(chǔ)上，將多屬性權(quán)衡空間探索和全局優(yōu)化思想有機(jī)結(jié)合，提出了彈道設(shè)計(jì)初值確定方法，并重塑了彈道設(shè)計(jì)流程。仿真結(jié)果表明：

1)按新流程設(shè)計(jì)的彈道滿足全部過程約束條件及最終入軌精度要求。

2)針對不同終端約束條件，可在不改變控制和目標(biāo)變量選取，且不改變設(shè)計(jì)初值的情況下，自動收斂得到設(shè)計(jì)結(jié)果。

3)方法解決了初值敏感性問題，具有快速、全局的特點(diǎn)，并與現(xiàn)行成熟的工程設(shè)計(jì)算法有機(jī)銜接，復(fù)雜度低，符合工程實(shí)際應(yīng)用需要。

應(yīng)注意的是，在火箭應(yīng)用發(fā)射常規(guī)任務(wù)彈道設(shè)計(jì)時(shí)，依然是傳統(tǒng)設(shè)計(jì)方法計(jì)算效率更高，本文給出的方法及流程主要面向新構(gòu)型、新任務(wù)類型下的設(shè)計(jì)工作。同時(shí)，方法針對如評價(jià)值f的截至條件等參數(shù)的取值，目前尚未形成統(tǒng)一的標(biāo)準(zhǔn)，要結(jié)合具體的任務(wù)問題來設(shè)置，因此，還需繼續(xù)深入探索，不斷提高方法的計(jì)算效率。