杜 昕，劉會(huì)龍

（中國空氣動(dòng)力研究與發(fā)展中心高速空氣動(dòng)力所，綿陽621000）

探月飛船跳躍式再入軌跡可達(dá)域分析

杜 昕，劉會(huì)龍

（中國空氣動(dòng)力研究與發(fā)展中心高速空氣動(dòng)力所，綿陽621000）

針對(duì)探月飛船跳躍式再入軌跡分段多、段與段相互耦合、可達(dá)域求解與分析較直接再入更加困難的問題，給出了跳躍式再入軌跡可達(dá)域的數(shù)學(xué)描述，在此基礎(chǔ)上將可達(dá)域求解問題拆分為兩類最優(yōu)控制問題，并建立了相應(yīng)的優(yōu)化模型。采用基于高斯偽譜法的兩步優(yōu)化策略進(jìn)行求解，得到了跳躍式再入軌跡可達(dá)域邊界。最后分析了初始條件（再入角、再入方位角）對(duì)可達(dá)域的影響。仿真結(jié)果表明兩步優(yōu)化策略能兼顧計(jì)算精度和計(jì)算效率。

探月飛船；跳躍式再入；可達(dá)域；高斯偽譜法；兩步優(yōu)化策略

1 引言

探月飛船返回地球時(shí)，將采用跳躍式再入以降低過載、熱流等［1］。再入軌跡可達(dá)域是工程上比較關(guān)心的一個(gè)問題，其分析結(jié)論可為任務(wù)方案總體設(shè)計(jì)提供參考。跳躍式再入能量更高，可達(dá)域范圍更大，但軌跡特性與直接再入差別較大，需要對(duì)其可達(dá)域問題進(jìn)行專門研究。

再入軌跡可達(dá)域問題實(shí)際上是多個(gè)最優(yōu)控制問題的組合，需要在沿初始速度方向的一系列縱程點(diǎn)上求取最大的橫程絕對(duì)值。目前，計(jì)算再入軌跡可達(dá)域主要有四種方法：第一種方法采用平衡滑翔假設(shè)來降低問題的復(fù)雜度，且不考慮路徑約束，通過坐標(biāo)變換，將問題轉(zhuǎn)化為一系列更為簡單的終端縱程自由、橫程最大的最優(yōu)控制問題，并采用間接法求解［2?3］；第二種方法是采用直接優(yōu)化算法（如偽譜法）進(jìn)行求解［4?6］；第三種方法基于航天飛機(jī)再入制導(dǎo)算法，其核心思想就是飛行器的縱程與飛行過程中的阻力加速度成反比［7?8］；第四種方法運(yùn)用擬平衡滑翔條件，將可達(dá)域問題轉(zhuǎn)化為一系列更為簡單的單變量尋根問題，并采用近似最優(yōu)控制律求解［9］。

第一種方法由Vinh在上世紀(jì)70年代提出。他對(duì)問題作了很多簡化，忽略了科氏加速度，并且不考慮路徑約束［2］。實(shí)際上，科氏加速度會(huì)造成可達(dá)域的非對(duì)稱性，是一個(gè)必須要考慮的因素。而不考慮任何路徑約束也與實(shí)際問題出入較大。第二種方法的最大缺點(diǎn)就是求解速度較慢，這也是直接法的缺點(diǎn)。第三種方法同時(shí)考慮了科氏加速度和路徑約束，能夠較為快速地計(jì)算可達(dá)域。第四種方法引入了擬平衡滑翔條件［10］，在考慮路徑約束的條件下能夠快速求得最優(yōu)解。總的說來，第三種方法求解快速、精度較高且使用范圍廣，是目前為止比較出色的計(jì)算落點(diǎn)可達(dá)域的方法，航天飛機(jī)的故障飛行管理（Shuttle Abort Flight Manager，SAFM）模塊正是采用了這種方法［11］。

第三種和第四種方法不適用于小升阻比的飛船再入，而第一種方法對(duì)于跳躍式再入精度較差。本文主要是離線計(jì)算可達(dá)域并對(duì)其進(jìn)行分析，對(duì)于計(jì)算速度的要求并不高，因此采用第二種方法進(jìn)行求解。首先，采用基于高斯偽譜法的兩層優(yōu)化策略求解最優(yōu)控制問題，第一層優(yōu)化算法采用較少的高斯節(jié)點(diǎn)給出一個(gè)精度不高的優(yōu)化結(jié)果，第二層優(yōu)化算法采用較多的高斯節(jié)點(diǎn)給出精度較高的優(yōu)化結(jié)果；其次，按照橫程的符號(hào)將整個(gè)可達(dá)域分為兩邊，同一邊里上一個(gè)邊界點(diǎn)的優(yōu)化結(jié)果作為下一個(gè)邊界點(diǎn)的優(yōu)化初值，逐點(diǎn)推進(jìn)，直至求解完成；最后通過仿真算例驗(yàn)證方法的有效性。

2 再入軌跡可達(dá)域問題描述

2.1三自由度再入運(yùn)動(dòng)方程

考慮地球?yàn)樾D(zhuǎn)圓球的三自由度再入運(yùn)動(dòng)方程如式（1）～（6）所示［12］：

其中：r為地心距，θ為經(jīng)度，?為緯度，V為飛行器速度，γ為飛行路徑角，ψ為速度方位角（正北方向?yàn)榱?，順時(shí)針旋轉(zhuǎn)為正），σ為滾轉(zhuǎn)角（控制變量），ω為地球自轉(zhuǎn)角速度，g為重力加速度，L和D分別為升力加速度和阻力加速度，定義如式（7）～（8）所示：

其中：ρ為大氣密度，m為飛行器質(zhì)量，S為飛行器參考迎風(fēng)面積。大氣密度模型采用美國1976標(biāo)準(zhǔn)大氣模型。

2.2可達(dá)域問題的數(shù)學(xué)描述

其中，σ 0，t( )表示從0時(shí)刻到t時(shí)刻的滾轉(zhuǎn)角剖面，其容許集Ut()可以表示為式（10）：

其中，C20，t[ ]表示在時(shí)間區(qū)間 0，t[ ]上具有連續(xù)一階導(dǎo)數(shù)和二階導(dǎo)數(shù)的連續(xù)函數(shù)域。

定義P?X為滿足約束的狀態(tài)域如式（11）：

當(dāng)高度達(dá)到10 km時(shí)，認(rèn)為達(dá)到開傘條件，因此，定義開傘狀態(tài)域D?X為式（12）：

可行的再入軌跡是指滾轉(zhuǎn)角函數(shù)在其容許集Ut()內(nèi)，滿足路徑約束，且終端狀態(tài)在D內(nèi)的軌跡?？蛇_(dá)域（Reachable Set，RS）是指從初始狀態(tài)x0出發(fā)的可行軌跡所能達(dá)到的終端狀態(tài)的集合，表示為式（13）：

而在再入問題里通常所說的可達(dá)域?qū)嶋H上指的是可達(dá)覆蓋區(qū)（Reachable Footprint，RFP），定義為式（14）：

RS和RFP的定義如圖1所示［13］。在本文后面的討論中，默認(rèn)可達(dá)域?yàn)镽FP，它是由再入點(diǎn)狀態(tài)、再入飛行器氣動(dòng)特性以及路徑約束共同決定，不同的再入點(diǎn)狀態(tài)會(huì)得到完全不同的可達(dá)域。在應(yīng)急或故障情況下，再入點(diǎn)的狀態(tài)不能事先預(yù)測，為了快速評(píng)估可能的著陸點(diǎn)以規(guī)劃應(yīng)急再入軌跡，需要在線快速生成可達(dá)域。

3 可達(dá)域求解策略

計(jì)算再入飛行器的軌跡可達(dá)域?qū)嵸|(zhì)是計(jì)算可達(dá)域的邊界。而求解完整的可達(dá)域邊界需要兩個(gè)步驟：一是優(yōu)化計(jì)算沿初始速度方向縱程的最小值和最大值；二是在縱程最小值和最大值之間選取一系列離散點(diǎn)，在這些點(diǎn)上優(yōu)化計(jì)算橫程的最小值和最大值（對(duì)應(yīng)橫程符號(hào)一正一負(fù)）。相應(yīng)的優(yōu)化模型如下：

1）路徑約束

再入問題的典型路徑約束為動(dòng)壓、過載以及熱流密度約束，如式（15）～（17）所示［13］。

此外，為防止跳出高度過高導(dǎo)致任務(wù)失敗，路徑約束還應(yīng)包括一個(gè)高度約束。這里設(shè)定最大高度不超過300 km，如式（18）所示：

2）終端約束

以高度作為判斷再入段結(jié)束的標(biāo)志，即當(dāng)高度達(dá)到開傘高度10 km時(shí)認(rèn)為再入段結(jié)束，達(dá)到開傘狀態(tài)。因此有終端約束如式（19）：

對(duì)于步驟一中的優(yōu)化模型，終端約束只有式（19）一項(xiàng)，但對(duì)于步驟二中的優(yōu)化模型，終端約束還應(yīng)包括式（20）：

其中，RV0為沿初始速度方向的縱程，其定義如圖2所示，圖中為再入初始點(diǎn)的經(jīng)緯度，為開傘點(diǎn)的經(jīng)緯度，為沿初始速度方向的橫程，Rt為再入初始點(diǎn)與開傘點(diǎn)之間的大圓弧距離。Rspec為選定的各個(gè)離散點(diǎn)處的RV0值。RV0、ZV0的計(jì)算公式如式（21）～（23）所示［14］。

其中，ψ為再入初始點(diǎn)的速度方位角，Ψ為再入初始點(diǎn)到開傘點(diǎn)的視線角。

3）目標(biāo)函數(shù)

對(duì)于步驟一中的優(yōu)化模型，目標(biāo)函數(shù)為式

（24）或者式（25）：

對(duì)于步驟二中的優(yōu)化模型，目標(biāo)函數(shù)為式（26）或者式（27）：

對(duì)于可達(dá)域問題，對(duì)求解速度的要求高于對(duì)求解精度的要求，因此，要選擇計(jì)算效率較高的優(yōu)化算法。本文選擇高斯偽譜法進(jìn)行優(yōu)化計(jì)算，并采用兩步優(yōu)化策略：第一步采用較少的高斯節(jié)點(diǎn)（例如取6個(gè)節(jié)點(diǎn)），計(jì)算得到一個(gè)較為粗略的結(jié)果作為下一步優(yōu)化的初值；第二步取較多節(jié)點(diǎn)（如取20個(gè)節(jié)點(diǎn)），并采用第一步的計(jì)算結(jié)果進(jìn)行插值，得到當(dāng)前多節(jié)點(diǎn)下的設(shè)計(jì)變量初值，計(jì)算得到一個(gè)精度較高的最優(yōu)解。高斯偽譜法最終將最優(yōu)控制問題轉(zhuǎn)化為非線性規(guī)劃問題，本文采用SQP算法求解所有的非線性規(guī)劃問題。

同時(shí)，在求解時(shí)，可達(dá)域邊界同一邊（橫程符號(hào)相同）上相鄰的兩點(diǎn)，上一個(gè)點(diǎn)的計(jì)算結(jié)果作為下一個(gè)點(diǎn)的初值，逐點(diǎn)推進(jìn)，可以提高計(jì)算效率。

4 仿真驗(yàn)證

再入點(diǎn)狀態(tài)為：再入點(diǎn)高度120 km，再入點(diǎn)經(jīng)度0°，再入點(diǎn)緯度0°，再入速度11 km／s，再入角－6°，再入點(diǎn)速度方位角45°。

飛船參數(shù)取Apollo 11號(hào)返回艙參數(shù)［1］，質(zhì)量5500 kg，大底截面直徑 3.9 m，駐點(diǎn)曲率半徑2.5 m，升力系數(shù)0.387 73，阻力系數(shù)1.289 1。

路徑約束如下［1，15］：q≤30 kPa，n≤7g0，Q·≤6 MW／m2。

第一步，先優(yōu)化得到沿初始速度方向最大的縱程和最小的縱程。計(jì)算得到最大縱程約為10 816 km，最小縱程約為1674 km。

第二步，在1674 km至10 816 km之間取一系列點(diǎn)，優(yōu)化得到這些點(diǎn)處的橫程最大值（正向最大和負(fù)向最大），這樣就求得了可達(dá)域的邊界。最后將這些邊界點(diǎn)連接起來，就得到了可達(dá)域。

優(yōu)化結(jié)果如圖3～8所示。其中圖3給出了可達(dá)域邊界軌跡的地面航跡，圖4給出了可達(dá)域，圖5給出了可達(dá)域邊界軌跡的高度曲線，圖6給出了可達(dá)域邊界軌跡的過載曲線，圖7給出了再入角對(duì)可達(dá)域的影響，圖8給出了再入航向角對(duì)可達(dá)域的影響。

從圖3～4中可以看到，跳躍式再入的可達(dá)域是一個(gè)長條狀的區(qū)域，若不考慮地球自轉(zhuǎn)，那么這個(gè)區(qū)域應(yīng)該是對(duì)稱的。另外，可達(dá)域邊界應(yīng)是光滑邊界，不光滑的點(diǎn)是由于優(yōu)化算法在這一點(diǎn)處計(jì)算的結(jié)果不太理想所導(dǎo)致，另外一個(gè)原因是由于離散點(diǎn)的個(gè)數(shù)取的不夠多。

計(jì)算程序采用MATLAB語言編碼，在安裝Pentium 4（2.5 GHz）處理器的微機(jī)上運(yùn)行。優(yōu)化得到一條邊界軌跡的時(shí)間平均約30 s，得到整個(gè)可達(dá)域的總時(shí)間約12 min。當(dāng)然，取點(diǎn)的密集度越大，計(jì)算總時(shí)間會(huì)相應(yīng)增長。

從圖5中可以看到，最短縱程對(duì)應(yīng)的兩條軌跡是直接再入軌跡，軌跡沒有出現(xiàn)跳躍。隨著縱程的增大，軌跡出現(xiàn)跳躍，并且跳躍高度越來越高。而且是兩條軌跡一組，高度曲線較為接近，這是因?yàn)檫@兩條軌跡的縱程一致，只是橫程一正一負(fù)。

從圖6中可以看到，最短縱程對(duì)應(yīng)的兩條軌跡的過載峰值幾乎達(dá)到了約束的最大值7g0，而其他軌跡的過載峰值則均小于7g0。與高度曲線相同的是：過載曲線也是兩條軌跡一組。這是由于高度、過載等均是縱向運(yùn)動(dòng)狀態(tài)決定的量，與橫程運(yùn)動(dòng)關(guān)系不大。

圖7 給出了再入角對(duì)可達(dá)域的影響。從圖中可以看到，再入角的絕對(duì)值越大，即初始再入軌跡越陡，可達(dá)域會(huì)更靠近再入點(diǎn)，但并不明顯，這是因?yàn)榇嬖谶^載峰值的約束，可行的再入軌跡不可能達(dá)到很小的縱程。從圖中還可以看到，初始再入軌跡越陡，可達(dá)域的范圍越窄，新的邊界點(diǎn)全部落在原邊界線以內(nèi)。

圖8 給出了再入航向角對(duì)可達(dá)域的影響。由于本文中航向角的定義是從正北起算，順時(shí)針方向?yàn)檎?。因此，?dāng)航向角增大時(shí)，可達(dá)域整體向南移動(dòng)；當(dāng)航向角減小時(shí)，可達(dá)域整體向北移動(dòng)。

5 結(jié)論

本文采用兩層直接優(yōu)化策略計(jì)算得到跳躍式再入軌跡的可達(dá)域并進(jìn)行了分析，結(jié)論如下：

1）建立了跳躍式再入軌跡可達(dá)域問題的優(yōu)化模型，并采用基于高斯偽譜法的兩層優(yōu)化策略進(jìn)行求解。仿真結(jié)果表明這種優(yōu)化策略能較好地兼顧精度和計(jì)算效率。

2）可達(dá)域邊界理論上為光滑邊界，不光滑的點(diǎn)是由于優(yōu)化算法在這一點(diǎn)處計(jì)算的結(jié)果不太理想以及離散點(diǎn)個(gè)數(shù)較少所導(dǎo)致，離散點(diǎn)越多，可達(dá)域范圍越精確，但計(jì)算時(shí)間越長。

3）分析了再入初始條件對(duì)可達(dá)域的影響。初始再入軌跡越陡，可達(dá)域的范圍越窄，且整體更靠近再入點(diǎn)。當(dāng)航向角增大時(shí)，可達(dá)域整體向南移動(dòng)；當(dāng)航向角減小時(shí)，可達(dá)域整體向北移動(dòng)。

4）本文采用直接法求解跳躍式再入軌跡可達(dá)域。由于直接法的計(jì)算速度難以適應(yīng)在線任務(wù)的要求，因此僅適用于線下的任務(wù)設(shè)計(jì)與分析。跳躍式再入軌跡可達(dá)域的在線生成還需要進(jìn)一步研究。

（References）

［1］ Bairstow S H.Reentry Guidance with Extended Range Capa?bility for Low L／D Spacecraft［D］.Massachusetts：Massachu?setts Institute of Technology，2006.

［2］ Vinh N X.Optimal Trajectories in Atmospheric Flight［M］. New York：Elsevier，1981：222?222.

［3］ Ngo A D，Doman D B.Footprint determination for resusable launch vehicles experiencing control effector failures［R］. AIAA 2002?4775，2002.

［4］ Fahroo F，Doman D B，Ngo A D.Modeling issues in footprint generation for reusable launch vehicles［C］／／2003 IEEE Aer?ospace Conference.Piscataway，NJ，2003.

［5］ Bollino K，Ross I M，Doman D B.Six?degrees?of?freedom trajectory optimization for reusable launch vehicle footprint de?termination［C］／／AAS／AIAA Astrodynamics Specialists Con?ference.Springfield，VA，2007.

［6］ Fahroo F，Doman D B.A direct method for approach and landing trajectory reshaping with failure effect estimation［R］. AIAA?2004?4772，2004.

［7］ Harpold J C，Graves C A.Shuttle entry guidance［R］. NASA?TM?79949，1979.

［8］ Saraf A，Leavitt J A，Mease K D.Landing footprints compu?tation for entry vehicles［C］／／AIAA Guidance，Navigation，and Control Conference and Exhibit.Providence，2004.

［9］ Lu P.Rapid generation of accurate entry landing footprints［J］.Journal of Guidance，Control，and Dynamics，2010，33（3）：756?767.

［10］ Lu P.Asymptotic analysis of quasi?equilibrium glide in lifting entry flight［J］.Journal of Guidance，Control，and Dynam?ics，2006，29（3）：662?670.

［11］ Jackson M C，Sreaube T M，Barrows T J，etal.Onboard de?termination of vehicle glide capability for the Shuttle abort flight manager（SAFM）［C］／／IEEE Aerospace and Electron?ic System Society，Core Technologies for Space Systems Con?ference.Piscataway，NJ，2002.

［12］ Lu P.Predictor?corrector entry guidance for low?lifting vehi?cles［J］.Journal of Guidance，Control，and Dynamics，2008，31（4）：1067?1075.

［13］ Benito J，Mease K D.characterizing the controllable entry states and reachable sites for planetary landing［C］／／6th In?ternational Planetary Probe Workshop.Atlanta，GA，2008.

［14］ 趙漢元.飛行器再入動(dòng)力學(xué)和制導(dǎo)［M］.長沙：國防科技大學(xué)出版社，1997：354?361. Zhao Hanyuan.Flight Vehicle Reentry Dynamics and Guid?ance［M］.Changsha：National University of Defense Tech?nology Press，1997：354?361.（in Chinese）

（責(zé)任編輯：龍晉偉）

Analysis of Reachable Sets of Lunar Module Skip Entry Trajectory

DU Xin，LIU Huilong

（High Speed Aerodynamics Institute of China Aerodynamics Research and Development Center，Mianyang 621000，China）

The multiple phases of the skip entry trajectory are coupled with each other，so the land?ing footprint computation and analysis of the skip entry is much more difficult than that of the direct entry.To solve this problem，the mathematical descriptions of the reachable set and the landing foot?print were given.On this basis，the landing footprint problem was converted into two optimal control problems，and the corresponding optimal models were established.Then，a two?step optimization strategy based on Gauss pseudospectral method（GPM）was used to obtain the landing footprint. The effects of initial condition（entry flight path angle and heading angle）on landing footprint were analyzed.The simulation results indicate that the two?step optimization strategy has good perform?ance on both computational accuracy and efficiency

lunar module；skip entry；reachable sets；Gauss pseudospectral method；two?step opti?mization strategy

V412

：A

：1674?5825（2017）02?0163?05

2015?12?14；

2017?02?17

國家自然科學(xué)基金（11372345）

杜昕，男，博士，工程師，研究方向?yàn)轱w行器軌跡規(guī)劃與制導(dǎo)控制。E?mail：232783682＠qq.com