采用滾動偽譜優化的組合動力飛行器上升段制導方法

唐湘佶,李兆亭,張洪波

(國防科技大學 空天科學學院,湖南 長沙 410073)

組合動力飛行器是未來空天飛行器發展的重要方向,包括吸氣式臨近空間高超聲速飛行器、水平起降空天往返飛行器[1]等。本文針對可水平起降的渦輪基組合動力(TBCC)飛行器,飛行任務場景設定為水平起飛,加速爬升后在臨近空間進行高超聲速巡航。該場景適用于高超聲速飛機、兩級入軌空天飛行器一子級[2]等飛行器。

上述飛行場景的任務需求和組合動力飛行器上升段的飛行特性要求制導系統具有較好的自適應能力。上升過程飛行跨度大,具有寬速域多模態氣動/推進強耦合特性[2];動力系統對控制及其動態過程的約束嚴格,性能受狀態與環境影響大;不同速度區間內動力形式不同,系統特性多變,非線性強;飛行器模型與環境存在不確定性,飛行過程中誤差累積效應顯著。復雜多變的系統和飛行環境加上兩者間的強耦合關系導致制導難度大,需要制導算法具有較好的精度以及自適應性。

上升制導中較為常見的標稱軌跡跟蹤制導對于較大偏差和擾動的糾偏能力有限,對于復雜多變的系統模型,制導律設計難度大,且難以滿足對自適應能力的需求。提高制導自適應能力的技術途徑之一是基于在線解算最優控制的制導技術,這一方法在工程中已被應用,例如運載火箭上升段的迭代制導[3]等。隨著計算機技術的發展以及數值計算方法的深入應用,一系列基于軌跡優化技術的制導新方法被提出[4],其中最優控制問題在線求解方法包括間接法[5-6]、直接法[4,7]和混合法[8]。

直接法中的偽譜法具有快速求解最優控制問題的潛力,可用于最優控制實時計算,實現最優反饋控制[9]。國內外對偽譜法在最優反饋控制問題及制導方面的應用進行了深入研究。ROSS等[10]將基于偽譜法的實時最優反饋控制應用于衛星姿態控制。BOLLINO等[11]研究了基于偽譜法的X-33飛行器再入段最優非線性反饋制導。張友安等[12]研究了偽譜法在高超聲速飛行器再入制導中的應用,理論上分析了閉環控制系統的有界穩定。閆循良等[13]將改進的偽譜反饋控制應用于運載火箭的遠程變軌制導,通過狀態量縮減降低了制導問題規模,保證了制導算法實時性。許東歡等[14]針對空天飛行器研究了基于偽譜法的閉環制導算法可行性。張志國等[4]研究了偽譜法在大氣層外運載火箭上升段制導律中的應用,與迭代制導的對比結果驗證了該方法的制導精度。上述研究表明了經典偽譜法可用于制導,但針對的問題約束相對較少,控制指令計算難度相對較低,部分研究中未探討方法對誤差的適應性及計算實時性。

本文針對渦輪基組合動力飛行器的上升過程,研究了采用滾動偽譜優化的制導律。應用團隊自研偽譜法求解器SPTOS,利用分段偽譜法在線求解最優控制。分段偽譜法是對經典偽譜法的改進,采用分段擬合提高收斂速率,該方法尤其適用于變量變化不夠平滑、有間斷點等情形[15]。制導過程以雙層滾動優化框架為基礎,引入離散網格的同步滾動更新改進求解效率,制導周期根據實時狀態偏差自適應調整。介紹了仿真及制導過程中所用動力學模型和相關約束,及采用滾動偽譜優化的制導方法,展示了仿真結果并進行了討論。

1 飛行器動力學模型及約束

針對制導過程,主要研究飛行器的質心運動,而加速爬升的運動過程主要在縱平面內。基于上述分析,建立飛行器質心運動方程及相關約束條件。

1.1 質心運動方程

僅考慮縱平面內的運動,將狀態量取為x=(HRvθm)T,分別表示高度,航程,速度,速度傾角以及飛行器質量。控制量取為u=(φtα)T,分別表示節流系數與攻角。采用圓球假設并忽略地球自轉,建立如下飛行器質心運動方程:

(1)

式中:re為地球半徑;ρ為大氣密度;FP,FL,FD分別為發動機推力、升力與阻力;KS為空氣/燃料流量轉換系數;Sc為進氣口截面積;CA(Ma,α)為等效進氣面積系數,是馬赫數和攻角的函數,通過插值計算。

氣動力的計算公式如下:

(2)

式中:Sref為飛行器參考面積;CL(Ma,α),CD(Ma,α)分別為升力系數和阻力系數,與馬赫數和攻角有關,通過插值計算。

本文考慮的組合動力系統包括3種工作模態:在Ma<2時是渦輪動力,Ma=2～6為沖壓模態,Ma>6時是超燃沖壓動力。推力可統一描述成如下形式:

(3)

式中:g0為水平面重力加速度;Isp(Ma,φt)為發動機比沖,與馬赫數和節流系數有關,通過插值計算。上述飛行器模型參數及插值數據參考文獻[16]。

組合動力系統的特性如圖1所示。主要性能參數受狀態和控制的影響大,由圖1(a)可見,隨著速度的增加,比沖整體減小;不同動力形式下比沖受速度的影響不同;節流系數提高,比沖增加。如圖1(b)所示,在低速段,進氣效率受攻角影響小;跨聲速段的進氣效率有所下降;在沖壓動力段,氣體壓縮過程使進氣效率不斷提升,由于需要借助機體實現來流的減速增壓,高速推進時的進氣效率受攻角影響明顯。

圖1 組合動力系統特性

1.2 約束條件

鑒于飛行任務要求及飛行器自身特性,實際飛行過程中存在諸多約束條件,本文考慮3類典型約束:狀態量約束、控制量約束及過程約束。

①狀態量約束。狀態量在初始時刻t0及末時刻tf的約束形式為

H(t0)=H0,v(t0)=v0,θ(t0)=θ0,m(t0)=m0
H(tf)=Hf,v(tf)=vf,θ(tf)∈[θmin,θmax]

②控制量約束。節流系數、攻角及其變化率均存在上下界,約束形式為

式中:KH為熱流密度常數。

上述符號中,下標min表示下界,max表示上界;下標0表示初始時刻的值,f表示終端時刻的值。

2 采用滾動偽譜優化的閉路制導方法

本文制導方法的核心為數值最優制導律。在雙層滾動優化框架下,通過外層滾動優化使剩余時間內制導指令隨狀態反饋更新,實現閉路制導過程,并通過制導周期的在線調整進一步提高方法的自適應能力;內層利用分段偽譜法實時求解最優控制,在內層迭代過程中引入網格更新策略以提高計算實時性。

2.1 基于雙層滾動優化的制導策略

制導策略如圖2所示。

圖2 基于雙層滾動優化的制導策略

為提高初次指令解算效率,將分段偽譜法離線規劃的最優軌跡作為首次指令更新時的參考。制導過程中單次指令解算的具體步驟如下:

①根據實時狀態反饋x(ti)建立多約束最優控制問題Pi(x(ti)),該問題的初始狀態與時刻即為狀態反饋x(ti)及相應的采樣時刻ti;

②通過分段偽譜法將連續問題離散化為非線性規劃(NLP)問題,離散時的初始離散網格根據前次解算的離散網格Mi-1進行更新;

③對前次優化結果Ui-1,Xi-1進行插值,獲得新的離散點處的迭代初值,然后求解NLP問題;

⑤當出現求解結果不收斂、求解時間過長等情況時,在相應時刻不進行制導指令更新,沿用前次計算結果,直至成功在線求解后再次更新。

制導周期Tg采用自適應調整方法,由實時狀態偏差Δx在線控制。自適應周期制導過程如圖3所示,圖中,tc為在線解算時間。制導過程中按一定采樣周期實時獲取飛行狀態,并計算與參考狀態的偏差,設定偏差容限ΔxLT,在Δx(ti)≥ΔxLT的一系列采樣時刻ti進行制導指令的在線解算與更新。實際制導周期Tg受偏差累積特性及在線解算時間tc影響,為可變周期,當解算時間非常短時該部分影響可忽略。

圖3 自適應周期制導過程

本文以高度偏差ΔH作為周期自適應調整依據,為適應飛行過程中不同的飛行狀態與不同階段的制導需求,建立偏差容限函數,形式為

式中:Hd為誤差容限切換高度。

制導周期的自適應調整能使指令更新頻率與實際狀態偏差累積特性相匹配。偏差容限的作用一方面在于避免誤差累積過大導致制導修正的難度過高;另一方面,當狀態偏差較小時,更新前后的控制量差別相對較小,指令解算速度更快。因此設定合適的偏差容限能夠較好地平衡計算效率與制導精度。

2.2 分段偽譜法求解最優控制問題

多約束非線性最優控制問題具有通用的標準形式,包括性能指標和各類約束。本文以燃料最省作為優化目標,等價于使終端質量最大,問題形式如下:

如圖4所示,分段偽譜法的基本原理是將連續問題直接離散化,利用數值優化算法求解;之后評估離散解對原問題解的逼近程度,若達不到精度要求,則進行離散網格的細化,通過構造序列NLP問題迭代求解;否則直接輸出結果。

圖4 分段偽譜法求解流程

本文選擇Legendre-Gauss-Radau(LGR)點作為離散點,其分布由Legendre正交多項式的根決定。LGR點的定義域為[-1,1),便于描述各分段點處連續條件。LGR點與各區間內時間具有如下關系:

離散化的關鍵是動力學微分約束的轉換,其余約束可直接應用于離散變量集{X,U}。動力學約束的轉換關鍵在于對微分算子的有限逼近。各區間歸一化后,分別用拉格朗日插值多項式近似原問題的狀態量,LGR點數目Nk確定時,插值基函數形式固定。計算公式如下:

將決策變量取為{X,U}及終端時刻tf,則原連續問題轉換成如下NLP問題:

式中:l=1,2,…,Nk;k=1,2,…,K。

本文運用開源求解器IPOPT求解上述大規模NLP問題,其基本原理為內點法,給定解的初始參考值,經過迭代可收斂至NLP問題最優解。序列NLP問題構造過程中的數值解評估方法及離散網格細化方法參考文獻[15]。

2.3 離散網格滾動更新

離散網格關系到NLP問題的規模及問題特性,是求解效率的主要影響因素之一。與分段偽譜法的網格細化不同,滾動更新針對初始離散網格。建立合適的初始離散網格能降低NLP問題難度,還能避免構造序列NLP問題,提高求解效率。

離散網格M的滾動更新過程為對3個決定參數(分段數K、各區間的占比hk與離散點數Nk)的迭代計算。根據文獻[17]中的收斂理論,對于分段偽譜法,當解的逼近程度一定時,區間長度與段內離散點數近似具有等比關系,區間越短,達到一定逼近程度所需要的插值多項式階數越低,離散點個數越少。在理想情況下,滾動優化過程中每次在線求解的結果是對上一次計算結果的截取。為相應調整NLP問題的規模,離散網格也應進行截取。

圖5 離散網格更新示意圖

在繼承離散網格Mi-1的基礎上,根據當前時刻t′0=t0+Δt,對初始離散網格的3個參數進行調整,方法如下:

①以前次成功求解時的離散網格Mi-1為參考,判斷當前時刻t′0所處區間順序z;

3 仿真分析

仿真中初始狀態與期望終端狀態如表1所示,各變量含義如前所述。

表1 初始狀態與期望終端狀態

表2 變量范圍約束

采樣周期取0.1 s,偏差容限函數參數ΔHmax1=30 m,ΔHmax2=10 m,Hd=43 km。

仿真硬件條件為Intel Core i7-4712MQ、2.3 GHz的PC。采用4階定步長Runge-Kutta方法進行數值積分,步長取0.05 s。

取大氣密度偏差為10%,升力系數偏差為10%,阻力系數偏差為10%,推力偏差為2%。組合極值拉偏仿真的結果如圖6和圖7所示,控制量及其變化率均在約束范圍內,終端速度傾角為0.07°,也滿足約束。

圖6 速度和高度隨時間變化曲線

圖7 攻角和節流系數隨時間變化曲線

由圖可見,加速到期望速度的時間相比于高度爬升過程更短,這是由于隨高度增加,大氣密度減小,耗油率增加,因此為節省燃料,優先完成加速過程,這體現出在誤差條件下,迭代更新后的制導指令仍保持最優控制的原則。雖然指令的整體動態變化顯著,但一個周期內指令是根據離散解進行插值得到,在小范圍內能夠保證控制量相對平緩,易于控制系統實現。

動壓與熱流密度均在上限約束范圍內,計算結果如圖8所示。仿真結果表明:制導算法具有可行性;采用滾動偽譜優化的制導方法能夠在滿足基本約束的條件下完成上升制導任務。

圖8 熱流密度與動壓隨時間變化曲線

指令更新結果與相應參考狀態的對比結果如圖9和圖10所示。以一個制導周期為例,為突出偏差影響,在高度-速度域內進行參考狀態量的對比。

圖9 參考狀態對比

由圖9可見,制導過程中,剩余時間內制導指令隨狀態反饋自適應修正,修正后的控制量對應的參考狀態相應改變。

指令迭代更新對比如圖10所示。外層滾動優化輸入狀態反饋輸出修正后的控制,狀態反饋隱含了誤差信息,為保持最省燃料爬升,剩余時間內攻角指令不斷調整。

以2個制導周期內的攻角指令為例,在引入狀態反饋后,剩余時間內控制發生改變;僅截取一個制導周期內的指令作用于飛行器,該周期內的指令變化較為規律。

自適應制導周期如圖11所示,在線求解次數為34次,每次求解均收斂至最優解。平均周期約為4 s,最小周期為0.6 s,此時高度超過Hd,指令更新頻率高。

圖11 自適應制導周期

在線計算耗時如圖12所示。由圖12可看出,平均在線計算時間約為270 ms,大部分時間在200 ms左右,滿足一定的實時性要求。在滾動優化的內層迭代過程中,離散網格同步滾動更新,改善了NLP問題的特性與規模;基于前次優化結果的迭代使優化算法收斂更快,提高了在線求解實時性。

圖12 在線計算時間

本文僅從計算方法的角度尋求在線快速求解。在實際應用中,結合硬件設備優化計算程序,可進一步提高求解效率。

為進一步驗證本文方法對誤差和飛行器特性的適應性,進行了500次的蒙特卡洛打靶仿真。各誤差項的3σ值如表3所示。

表3 蒙特卡洛打靶仿真誤差項

終端高度誤差與速度誤差的統計結果如表4所示。由表可見,狀態偏差均值為負,在最省燃料爬升過程中推進效率接近最優,在誤差條件下容易出現飛行器爬升能力趨近飽和,控制修正能力有限,使加速爬升的末端高度狀態低于期望值。

表4 打靶仿真統計結果

隨著周期自適應調整,任務末段更新頻率的提高保證了一定的終端制導精度。終端高度與速度精度整體均在可接受范圍內,在不同誤差條件下,制導方法均能夠適應復雜多變的組合動力系統特性。

4 結論

本文研究了采用滾動偽譜優化的閉路制導方法。在雙層滾動優化框架下,根據實時狀態反饋建立多約束最優控制問題,加入離散網格同步更新,運用分段偽譜法進行快速計算,并在線自適應調整制導周期,從而實現制導指令的迭代更新。對制導過程進行數值仿真,主要結論如下:

①制導方法的精度滿足任務需求,在線求解過程始終滿足優化指標,具有最優反饋性質。根據狀態反饋進行指令的迭代更新,能夠抑制模型參數不確定和外界干擾引起的誤差影響,使飛行器達到期望的飛行狀態,方法具有一定的魯棒性。

②引入離散解迭代以及離散網格滾動更新,提高了分段偽譜法在線求解多約束最優控制問題的效率,使該方法滿足一定實時性要求。

③制導周期的自適應調整能夠較好地平衡制導效率與制導精度,提高了制導方法自適應性。