藺君, 何英姿, 黃盤興
(1.北京控制工程研究所, 北京 100190; 2.空間智能控制技術重點實驗室, 北京 100190)
高超聲速飛行器再入包含飛行器控制、制導與規劃等方面,是一個高度復雜的綜合性問題。再入軌跡規劃是分析飛行器再入能力的重要手段,是現代飛行器設計中的重要內容。特別是針對現代隨控布局飛行器,軌跡規劃影響著總體、氣動布局、制導控制、動力和結構等多個分系統的設計。
再入軌跡規劃包含直接法和間接法。直接法通過將連續的最優控制問題離散化、參數化,轉化為包含優化指標及參數約束的參數優化問題,利用優化方法得到最優解。間接法則基于極值原理和變分法將最優控制變量表示成狀態變量和協態變量的函數,將最優軌跡優化問題轉化為哈密頓邊值問題(HBVP),通過求解該HBVP,獲得相應的最優控制變量和最優軌跡[1]。
間接法求解最優軌跡問題時,協態方程與橫截條件等求解過程復雜繁瑣,且對初值估計精度要求較高。高超聲速飛行器的非線性特性及再入過程的多種約束條件,進一步增加了間接法求解過程的復雜度[2]。無實際物理意義的協態變量加大了初值估計的難度;在含有路徑約束時,必須將路徑約束轉化為終端約束才能解算;含隱式約束條件時,相應的Legendre乘子很難消去,進一步加大了計算量。間接法的這些局限性限制了其應用和發展。直接法根據參數化方法的不同,可分為僅離散控制變量的直接打靶法和同時離散控制變量與狀態變量的配點法。單純離散化和參數化控制輸入的直接打靶法對初值較為敏感。直接配點法則同時離散化和參數化控制和狀態,為了克服等距離配點易引起“龍格”現象,正交配點法(也稱為偽譜法)選取全局插值多項式近似狀態量和控制量[3]。
高斯偽譜法(GPM)是Benson在Legendre偽譜法基礎上提出的一種改進方法[4]。GPM將狀態和控制時間歷程在一系列高斯點上離散化,然后用這些離散的狀態與控制分別構造Lagrange插值多項式去逼近真實的狀態與控制,再通過對狀態量求導來代替動力學微分方程,將連續系統最優控制問題轉化為受一系列代數約束的參數優化問題。GPM可保證轉化得到的非線性規劃問題卡羅需- 庫恩- 塔克(KKT)條件均與其1階最優性條件兩點邊值問題的離散形式完全等價。Huntington對Benson的工作做了進一步擴展,證明了在動力學約束和路徑約束同時存在的條件下,用GPM離散得到的非線性規劃問題KKT條件仍與HBVP等價[5]。
基于偽譜法的軌跡優化,已在貨運飛船返回艙返回軌道設計[6]、火箭上升段軌跡優化[7-8]及在線軌跡優化和再入制導[9]、月球定點著陸[10]、制導炸彈最大射程優化[11]、高超聲速飛行器再入軌跡優化[3]以及帶有禁飛區[12-14]、航路點約束[15]、總吸熱量約束[16]的再入軌跡優化等領域得到應用。為了提高軌跡優化的時效性,火箭返回著陸軌跡優化結合凸優化理論及方法,對推力進行無損松弛,快速得到再入優化軌跡[17]。
然而,針對高超聲速飛行器再入軌跡優化這一問題,現階段研究較多的是針對無動力再入展開。為了提高再入過程飛行器的任務靈活性和機動能力,發展帶推力高超聲速飛行器,可在再入過程中進行多次點火、熄火,從而改變飛行器再入軌跡和落點。本文針對再入過程中多次點火這一需求,展開帶推力高超聲速飛行器再入軌跡優化研究。
考慮如下帶推力高超聲速飛行器再入動力學模型[18]:
(1)
式中:r為地心距;θ為經度;φ為緯度;v為飛行器飛行速度;γ為航跡角;ψ為航向角;T為發動機推力;α為攻角;FD為阻力,FD=ρv2CDSref2,ρ為大氣密度,CD為阻力系數,Sref為飛行器參考面積;m為飛行器質量;g為重力加速度;Ω為地球自轉角速度;FL為升力,FL=ρv2CLSref2,可通過國際大氣數據表查得,CL為升力系數;σ為傾側角;ce為發動機燃料質量流量。
為了保證飛行器結構安全,飛行器再入過程需要滿足如下動壓、過載、熱流密度等強約束條件:
動壓:
Q=ρv2/2≤Qmax,
(2)
過載:
(3)
熱流密度:
(4)

為成功保證飛行任務,飛行器接近任務終端時,需要滿足終端約束條件,通常包括如下約束:
終端高度約束:
|hf-h|≤Δh,
(5)
式中:hf為終端高度;h為高度;Δh為高度允許偏差。
終端速度約束:
|vf-v|≤Δv,
(6)
式中:vf為終端速度;Δv為速度允許偏差。
落點位置約束:
|θf-θ|≤Δθ,|φf-φ|≤Δφ,
(7)
式中:θf為終端經度;Δθ為經度允許偏差;φf為終端緯度;Δφ為緯度允許偏差。
飛行器再入過程中需要考慮控制量約束,即
|α|≤αmax,|σ|≤σmax,
(8)
式中:αmax為最大攻角;σmax為最大傾側角。
帶推力再入飛行器還應滿足如下發動機工作總時長約束:
(9)
式中:N為發動機點火次數;κν為第ν次點火時長;tz由發動機比沖及推進劑總量決定。
考慮如下一般性非線性系統:
(10)
式中:t為時間;x(t)為狀態;u(t)為控制輸入;t0為初始時間;tf為終端時間。
再入過程應滿足的邊界條件為
Υ(x(t0),t0,x(tf),tf)=0,
(11)
路徑約束為
Γ(x(t),u(t),t)≤0.
(12)
對一般性的Bolza型最優控制性能指標,形如(13)式,

(13)
式中:Φ(x(t0),t0,x(tf),tf)表示初始和終端狀態性能指標;G(x(t),u(t),t)表示積分性能指標。

(14)
(11)式可改寫為
Υ(x(-1),t0,x(1),tf)=0,
(15)
(12)式可改寫為
Γ(x(τ),u(τ),τ)≤0,
(16)
(13)式可改寫為

(17)
高斯偽譜法通過選取M階Legendre-Gauss(LG)點,并以τ0τ(0)=-1為初始節點,在LG點對系統方程進行離散[19]。
(14)式在選取的LG節點可改寫為代數約束形式,即

(18)
式中:Dki為微分矩陣,由Legendre多項式遞推關系及其導數計算公式迭代計算;X(τk為狀態變量x在LG節點τk處的值;U(τk為輸入量u在τk處的值;f(X(τk),U(τk);t0,tf)表示(10)式在LG節點τk處的值。
(17)式由高斯積分公式,可近似為

(19)
式中:τf為終端歸一化時間;ωk為LG點權值;G(X(τk),U(τk),τk;t0,tf)表示積分性能指標G(x(t),u(t),t)在LG節點τk處的值。
則(14)式在終端時刻可轉化為代數約束,即

(20)
同理,(15)式在LG節點離散化后,近似為
Υ(X(τ0),t0,X(τf),tf)=0,
(21)
于是(16)式變為
Γ(X(τk),U(τk),τk;t0,tf)≤0,k=1,2,…,M.
(22)
最終,(14)式~(17)式的解,可通過求解含有多種約束條件(18)式~(22)式的非線性規劃來確定。
帶推力高超聲速飛行器發動機開關機時間對再入軌跡影響顯著。發動機的不連續點火將再入軌跡分割成多個階段,再入軌跡在優化過程中需要考慮分段點處系統狀態和控制輸入的連接問題。
假設1每次點火時刻為tν,點火推進時長為κν,ν= 1, 2,…,N.
忽略發動機開關機的動態特性,帶推力高超聲速飛行器再入過程中發動機點火N次,則利用高斯偽譜法進行再入軌跡優化問題被分為2N+1段。若每段包含LG節點數均為K,則分段高斯偽譜法總的配點數為(2N+1)×(K+1)+1. 假設系統動力學方程的維度為n,系統控制輸入維度為q,則分段高斯偽譜法總的優化參數量為((2N+1)×(K+1)+1)×(n+q)。

(23)
從而總的性能指標為
(24)
在各個時間區間,路徑約束和邊界約束為
(25)

(26)

(14)式在離散高斯點轉化為如下多段代數方程:
(27)
在各個時間區間內,終端狀態約束為
(28)
為了保證各段區間的連續性,系統狀態和時間應滿足如下條件:
(29)
在各個時間節點,系統控制輸入及飛行器質量還應滿足如下連續性約束:
(30)
帶推力高超聲速飛行器采用分段GPM進行再入軌跡優化時,飛行器的控制輸入選取為攻角和傾側角。
在各個時間節點,系統控制輸入及飛行器質量還應滿足連續性約束。即對于采用火箭發動機在再入段進行推進以增加飛行器機動能力和飛行航程時,推進劑耗盡后,如果本級火箭推進器被拋掉,則此時飛行器總質量應滿足:
(31)

由于飛行器重量及結構限制,再入飛行器多級點火增加飛行器航程及機動能力的方式不能采用過多級的推進器。本文假設飛行器可最多進行2次點火助推,以兼顧性能提升和結構限制。
通過如上數值離散化,帶推力高超聲速飛行器再入軌跡優化問題即可以轉化為分段高斯偽譜法求解目標函數(23)式和(24)式、路徑約束(25)式、邊界條件約束(26)式、系統動態過程(27)式、終端狀態約束(28)式及分段高斯偽譜法連續性約束(29)式、(30)式或(31)式的最優控制問題。
以美國波音公司CAV-H高超聲速飛行器為例,對再入飛行器無動力再入與帶推力再入進行對比。飛行器質量m=907.2 kg,參考面積Sref=0.483 9 m2,Kn=1.65×10-8.
由CAV-H飛行器氣動力系數,將其進行近似化處理,阻力系數和升力系數[20-21]可描述為
CD=0.024 67+0.000 714 3α2+0.325 2e-0.279Ma,
(32)
CL=-0.234 2+0.051 36α+0.294 3e-0.100 7Ma,
(33)
式中:Ma為馬赫數。
由(32)式和(33)式給出的飛行器阻力系數和升力系數,可得到CAV-H飛行器最大升阻比攻角剖面,如圖1所示。

圖1 最大升阻比攻角剖面Fig.1 Angle of attack profile for max lift-to-drag ratio
假設CAV-H飛行器再入的初始條件和終端條件如表1所示。

表1 再入初始和終端條件

選取目標函數為橫向航程最大,在選取初始航向角為90°時,最大橫向航程可等效為最大緯度,即J=maxφf,φf為終端時刻緯度。在MATLAB軟件仿真環境下,利用GPM進行軌跡規劃,無動力再入過程優化變量為攻角和傾側角。采用文獻[20]中的方法,得到無動力再入軌跡。
仿真結果表明,再入軌跡的最大橫向航程可達48.136 7°,約為5 358.4 km,縱向航程可達50.690 5°,約為5 642.7 km. 高超聲速飛行器再入軌跡,包括高度、速度、航跡傾斜角、航跡方位角、終端條件均可滿足,且再入攻角和傾側角均位于控制輸入邊界以內;動壓、過載、熱流密度均小于允許的最大值。無動力再入時,高超聲速飛行器再入軌跡較平滑。
在無動力再入軌跡優化基礎上,分別在不同時刻施加發動機推力,并考慮再入過程路徑約束和邊界約束,研究發動機點火時刻以及點火助推時長對再入飛行器估計的影響。選取如圖2所示A、B、C、D共4個時刻,分別對應再入初始時刻A、再入過程中度過最大熱流密度時刻后,飛行器高度再次達到最高點時刻C、點火時刻A和點火時刻C對應中間高度時刻B,此時飛行器處于接近最大下降速度以及將近飛行器接近半程航程時刻D.
3.3.1 再入點火時刻分析
假設2發動機空質量為18 kg,推力為10 kN,發動機比沖Isp為260 s,助推時間為20 s.
發動機點火后,一次性充分燃燒推進劑,完成助推。由ce=T/(Ispg0)(g0為海平面重力加速度)可知,質量流量為3.925 2 kg/s. 由假設2可知,推進劑總質量為78.5 kg.

圖2 帶推力再入點火時刻Fig.2 Ignition timing of powered reentry
帶推力飛行器再入初始質量為1 003.7 kg,采取與無動力再入相同的初始條件(見表1)及約束條件,在給定A、B、C、D4個時刻分別進行點火助推。帶推力飛行器達到終端時刻質量為飛行器質量與發動機空質量之和,為925.2 kg. 帶推力再入與無動力再入終端狀態如表2所示。
圖3所示為點火時刻助推軌跡。由圖3可見,發動機點火助推可滿足再入軌跡對終端狀態的約束包括速度和高度要求。飛行器在再入初期維持大攻角,飛行器氣動力達到一定臨界值后,飛行器轉入最大升阻比攻角飛行器,以獲得最大的飛行航程和側向航程。由圖3(a)可見,發動機在時刻A和時刻B點火后,飛行器再入軌跡變化不大,基本與再入無動力軌跡類似,但飛行時間、橫程、縱程較無動力明顯增大。在點火時刻C和時刻D點火后,飛行器軌跡出現小幅波動。飛行器在到達時刻C之前,氣動力已經產生明顯作用,通過點火增加飛行器速度后,為了滿足飛行器終端約束保證速度和高度要求,再入軌跡需要進行多次調整,導致軌跡出現一定波動。

表2 不同點火時刻助推及無動力再入終端狀態對比
由圖3仿真結果可知,帶推力高超聲速飛行器可大幅提高飛行器再入射程。發動機點火時間越早,飛行器射程增加越大,但飛行器飛行時間也越長,導致飛行器總的吸熱量增加。飛行器在時刻A和時刻B點火需要進行攻角和傾側角的雙重調整,從大攻角調整至11°左右的最大升阻比攻角,傾側角調整至0°,需要利用飛行器姿控發動機進行調姿,能量消耗較大。在時刻C和時刻D點火僅需要利用飛行器氣動力進行操控,能量損耗小。
3.3.2 再入點火時長分析
在圖2所示的4個時刻火點,分別施加不同推力大小、不同助推時長的再入點火助推形式,對飛行器最終狀態進行對比,如表3所示。

圖3 不同點火時刻助推軌跡Fig.3 Powered reentry trajectories at different ignition times

表3 不同點火助推時長終端狀態對比
由表3可知:在再入初始時刻A施加不同時長的助推,相當于改變再入初始速度,對飛行器總的飛行時間有影響,但改變較小;在時刻B點火,由于飛行器氣動能力不足,難以保持平衡滑翔姿態,點火助推加劇飛行器下降速度,導致再入飛行時間變短;在時刻C和時刻D,升力已足夠維持飛行器平衡滑翔,點火后,飛行時間均有小幅增長。由此也可推斷,在時刻C之后點火助推,飛行器飛行時間較無動力時,均會增加,且助推時間越長,總飛行時間也相應增加。在4個時刻點火,當總沖一定時,飛行器終端狀態變化不大。即大推力短時間助推和小推力長時間助推,在選取最大橫向航程的目標函數下,對飛行器終端改變基本一致。
為了進一步說明發動機推力在不同點火時刻的作用,將發動機推力分為兩次進行點火助推,分別在時刻C和時刻D進行增程助推設計。
假設3再入過程中發動機可進行兩次點火,即p=2,且點火時長均為10 s.
帶推力高超聲速飛行器的再入初始狀態和終端狀態與無動力再入時相同(見表1),目標函數為橫向跨度最大。根據假設3,完整的再入軌跡分為5段,飛行器質量為

再入過程中發動機點火時,攻角維持在最大升阻比攻角,傾側角維持在0°. 采用分段高斯偽譜法和數值積分對帶推力高超聲速飛行器進行再入軌跡優化后,得到再入軌跡如圖4所示。
由地面經緯圖4(a)顯示,帶推力再入軌跡的最大橫向航程可達51.754 8°,約為5 761.2 km,縱向航程可達54.054 0°,約為6 017.1 km,飛行器射程為8 330.5 km. 高超聲速飛行器再入軌跡如圖4(b)~圖4(e)所示,包括高度、速度、航跡傾斜角、航跡方位角,且均滿足終端約束條件。再入攻角和傾側角如圖4(f)~圖4(g)所示,再入攻角和傾側角均位于控制輸入邊界以內,在分段點處對發動機進行兩次點火,分別產生推力,使飛行器速度增加,再入攻角和傾側角對應當前環境的最大升阻比攻角和0°傾側角。動壓、過載、熱流密度如圖4(h)~圖4(j)均小于允許的最大值。由圖4(k)帶推力再入三維軌跡可見,高超聲速飛行器可平穩滑翔飛行。帶推力高超聲速飛行器質量變化如圖4(l)所示。
圖5所示為高超聲速飛行器帶推力與無動力再入軌跡對比。由圖5可見:高超聲速飛行器帶推力較無動力再入時,經度增大約為3.363 5°,緯度增大3.618 1°,飛行距離增加約為549.907 7 km;再入過程中由于發動機兩次點火以及傾側角保持在0°,再入軌跡較無動力時,呈現一定的波動。
圖6所示為高超聲速飛行器不同點火時刻助推與無動力再入攻角對比。由圖6可見,在點火時刻C和D,攻角基本維持在最大升阻比攻角,因此將點火后助推段的攻角剖面選擇為最大升阻比攻角[22]。傾側角維持在0°,以保持飛行器升力僅存在于飛行航跡方向,從而盡量拉升飛行器,以取得最大的航程。

圖4 分段GPM帶推力再入軌跡優化Fig.4 Powered reentry trajectory optimization using MGPM

圖5 帶推力與無動力三維再入軌跡對比Fig.5 Comparison of 3D powered and unpowered reentry trajectories

圖6 不同點火時刻助推與無動力再入攻角對比Fig.6 Comparison of angles of attack at different ignition times during powered and unpowered reentries
利用高斯偽譜法得到的飛行器再入攻角和傾側角指令后,飛行器利用數值積分方法可得到在標稱軌跡下的飛行器飛行軌跡。表4所示為跟蹤標稱軌跡終端狀態對比。由表4可見,帶推力高超聲速飛行器可較好地實現對標稱軌跡的跟蹤,飛行器終端狀態誤差很小。

表4 跟蹤標稱軌跡終端狀態對比
本文給出基于分段GPM和數值積分對帶推力高超聲速飛行器非連續點火再入軌跡進行優化。優化求解到的軌跡滿足分段GPM連續性約束及過載、動壓、熱流密度約束。分段GPM方法將再入軌跡在發動機點火時刻進行分段,發動機關機時利用分段GPM進行軌跡優化,發動機開機后利用數值積分,計算再入軌跡。
數值算例結果顯示,分段GPM和數值積分可以有效生成滿足連續性約束及常規約束的再入軌跡,在總沖一定時,發動機點火時刻對再入軌跡影響明顯。高超聲速飛行器兩次點火產生的速度增量可有效增加飛行器的縱向航程和橫向航程,并有效提高飛行器再入機動能力和靈活性;分段GPM規劃得到的控制輸入可作為參考標稱軌跡指令,輔助制導系統設計。