多滑翔飛行器時間協同軌跡快速規劃

劉哲,陸浩然,鄭偉,*,聞國光,王奕迪,周祥

1. 國防科技大學 空天科學學院,長沙 410073

2. 北京宇航系統工程研究所,北京 100076

3. 北京交通大學 理學院,北京 100093

高超聲速滑翔飛行器是一種在臨近空間飛行,最大速度可達馬赫數20以上的升力式飛行器,具有速度快、側向運動能力強、精度高、強突防的特點。然而隨著反導武器防御系統的不斷完善和發展,單個飛行器的生存能力下降,傳統的多彈同時飽和攻擊方式突防能力有限。因此,高超聲速飛行器集群中的個體可通過精確協同控制再入總時間及終端狀態的方式,多角度、多批次抵近目標,從而實現“探測-打擊-評估”一體化協同再入,有效提升集群的綜合效能。

相關學者從2005年起開始對多彈協同算法的研究。文獻[1]最早設計了固定總時間約束下的單彈最優制導律(ITCG),實現了多彈同時到達目標點；文獻[2-3]分別提出了“領彈-從彈”及雙層協同架構,后續的多彈協同末制導方法多基于這兩種架構,協同方式主要集中在時間協同。近年來,文獻[4]首先論證了高超聲速飛行器與常規導彈協同的可行性,而后對于高超集群協同軌跡規劃與制導的研究逐漸興起。

相比于多彈、多無人機協同,滑翔飛行器集群為欠驅動系統,控制能力有限,再入環境復雜,多彈協同方法無法直接應用,因此多采用“標準軌跡法[5]”“預測校正法[6-9]”或優化方法實現多滑翔飛行器時間協同。文獻[5-6]通過調節航向角誤差走廊寬度的方式改變側向運動幅度,間接控制再入總時間；文獻[6]在文獻[5]的基礎上,采用神經網絡擬合了航向角誤差走廊邊界與到達時間的關系,動態改變走廊邊界進行時間調節。文獻[7-8]基于預測校正方法,直接迭代調整控制量剖面以實現時間控制；文獻[9]利用模型預測靜態規劃(MPSP)算法替代了預測校正的迭代過程；文獻[10-11]將“領-從”模式及ITCG方法應用到了多高超聲速飛行器的末制導段。上述方法存在的問題在于: ① 難以獲取再入總時間的可調范圍并確定合理的協同時間; ② 標準軌跡方法中航向角誤差走廊的增大會導致終端誤差的增大; ③ 預 測校正方法一般只能滿足終端位置與協同時間需求,很難同時滿足終端角度約束,迭代求解的收斂性受初值影響較大。

為滿足終端位置與角度約束,偽譜法等優化方法將問題離散化后求解,但計算用時較長,且同時滿足協同時間約束的文獻相對較少。MPSP算法通過解整個時間段的協態向量來更新控制變量解,不需要直接對非線性動力學模型采用線性化假設,計算效率、數值精度高,但大多需要對控制量進行預設計以獲取初始猜想解[12-13]。

近年來,序列凸化算法以其快速收斂性得到了廣泛應用[14-25]。文獻[14]利用序列二階錐規劃(SOCP)求解了再入飛行器的最速下降軌跡,但未涉及給定協同時間約束下的規劃方法；文獻[15]直接將協同時間加入約束中,采用平面相對運動學模型求解,未考慮實際的動力學環境；文獻[16-17]分別通過減弱模型非線性,避免凸化處理的方式提高算法收斂性。此外,相關學者通過在求解過程中調整信賴域半徑的方式加快收斂速度,主要包括將信賴域作為“硬約束”與“軟約束”兩類[18-22]。文獻[18-20]將信賴域作為硬性約束加入到求解模型中,并根據各次迭代過程的線性化誤差等指標調整信賴域半徑,但求解給定協同時間下的規劃問題時,常出現存在多條最優軌跡的情況,導致解會在兩條最優軌跡間不斷振蕩,序列迭代算法難以收斂；文獻[21-22]分別將信賴域及線性化誤差加入罰函數中,從而在保證線性化誤差盡可能小的同時加快收斂速度,但此類方法在處理本文問題時會使得罰函數中懲罰項過多,各項的系數大小難以確定,不利于實際應用。同時,不合適的信賴域懲罰項更新策略會導致算法收斂性下降,甚至搜索不到最優解。

針對以上問題,本文首先給出了“探測-打擊-評估”一體化模式下集群再入的協同規劃方案與數學模型,進一步將協同軌跡規劃問題轉化為3個子問題進行求解,即：飛行器再入總時間可調范圍的快速獲取、期望協同時間的確定、給定協同時間約束下的軌跡規劃。將飛行路徑角初值剖面定義為二次曲線,并將終端經緯度偏差、總時間偏差以及飛行路徑角相對于預設剖面的偏離程度加入罰函數,提高了問題的可解性。提出了一種罰函數與信賴域系數自適應策略,在求解過程中不斷調整目標函數中時間項所占權重,并根據罰函數值的變化收縮信賴域,從而在飛行路徑角平滑剖面附近逐步搜索得到最優解,提高了序列迭代算法求解目標函數中含時間項問題的收斂性。實現了不同再入任務下多滑翔飛行器時間協同軌跡的快速規劃,并與偽譜法及MPSP算法求解結果進行了對比。

1 集群時間協同再入問題描述

1.1 集群時間協同再入策略

為實現“探測-打擊-評估”一體化協同再入,采用兩種時間協同方案： ① 不同起始點的多個集群確定協同再入時間,同時到達目標點,實現多方向飽和攻擊； ② 單個集群內部各飛行器通過采取不同軌跡形式控制再入時間,從而使同一集群內部成員多批次到達目標點。其中前序飛行器負責探測,并將態勢信息向后傳遞,后序飛行器負責攻擊及效果評估,實現時域協同。圖1 給出了“探測-打擊-評估”一體化集群協同再入模式的示意圖。

圖1 “探測-打擊-評估”一體化協同再入

兩種協同再入方式首先均需快速確定集群內部各飛行器的再入總時間范圍,確定協同到達時間后對各子飛行器軌跡進行精確的時間控制。可見,高超聲速集群滑翔段時間協同對軌跡解算的快速性有較高要求,是一種新的再入軌跡規劃問題,即在滿足式(1)、式(6)、式(12)、式(13)的約束條件下,優化得到再入總時間范圍,進而進行固定期望協同時間約束下的軌跡規劃。

1.2 飛行器數學模型

對于高超聲速滑翔飛行器再入集群中的第i個成員,其以能量為自變量的無量綱三自由度動力學模型為

(1)

(2)

式中：ρ為大氣密度；Sref為參考面積；m為飛行器質量；CL和CD分別為升力系數、阻力系數,其值由當前時刻的攻角α及馬赫數Ma二維插值得到。模型的自變量為無量綱化后的能量：

e=1/r-V2/2

(3)

則在地心距為r處,飛行器的速度為

(4)

再入過程中控制量為攻角α和傾側角σ,攻角剖面一般設置為固定形式,則再入軌跡可由傾側角剖面σ(e)唯一確定。

1.3 約束條件

1.3.1 控制量約束

在進行滑翔段軌跡規劃時,采用“固定攻角剖面+可控制傾側角剖面”的控制方式。根據工程經驗,攻角一般定義為線性剖面的形式：

(5)

式中：αmax、αmax L/D分別為滑翔飛行器所允許的最大攻角及最大升阻比所對應的攻角；V1、V2為攻角剖面轉折點對應的速度值。本文限制飛行過程中升力方向一直向上,則傾側角σ的允許變化范圍可以表示為

σmin≤|σ|≤σmax

(6)

為了避免控制量振蕩對于軌跡求解收斂性的影響,采用文獻[14]中的控制量選取方式,轉化后的新控制量可分解為向量形式u=[u1,u2]：

u1=cosσ,u2=sinσ

(7)

則僅需對控制量的分量進行約束：

(u1)min≤u1≤(u1)max

(8)

由于分量u1為標量,因此式(8)滿足凸約束條件。此外控制量還需滿足：

(9)

1.3.2 過程約束

(10)

式中：KQ=9.43×10-5,標準大氣密度ρ可展開為

ρ=ρ0e-h/hs

(11)

式中:hs為基準高度。

h≥hQ(e),h≥hq(e),h≥hn(e)

(12)

1.3.3 終端約束

高超聲速飛行器滑翔段結束后,為使飛行器能夠充分利用自身機動能力順利到達預定目標點,需要對滑翔段結束點的位置、速度等狀態量進行約束,則當飛行器到達終端能量狀態ef時,終端狀態量需滿足：

(13)

2 考慮協同時間約束的改進序列凸化方法

由于目前相對缺少對多滑翔飛行器協同軌跡規劃算法的研究,協同時間約束難以滿足,因此本文針對實現“探測-打擊-評估”一體化協同再入的任務需求,提出了協同到達時間約束下的軌跡規劃方法。首先,根據兩種集群時間協同再入形式,分別提出了協同時間的確定策略。而后利用序列凸化算法快速解算協同軌跡。圖2給出了時間協同軌跡規劃的具體流程。

圖2 時間協同軌跡規劃流程

2.1 時間協同軌跡規劃問題

2.1.1 協同時間的確定策略

首先,通過時間協同決策模塊綜合飛行器及任務目標點的狀態信息,利用序列凸化算法快速優化得到各飛行器滑翔段的總時間可調范圍,進而由上層時間協調模塊根據任務場景確定協同到達時間或到達時間序列。接下來,將協同時間下發至各飛行器,各飛行器分別進行固定期望協同時間約束下的軌跡規劃問題求解。不同任務場景下的協同時間確定策略為

1) 多個飛行器同起點同時開始滑翔,序列多波次到達同一目標點。

步驟1利用序列凸化算法快速解算該起點到目標點的滑翔最長時間tmax以及最短時間tmin。

步驟2按集群內部飛行器個數確定到達時間間隔Δt,生成到達時間序列ti：

(14)

式中：n為集群中飛行器總個數；i為飛行器編號,i=1,2,…,n-1。

步驟3將協調時間指令ti下發至各個飛行器,各子飛行器分別進行固定滑翔段總時間約束下的軌跡規劃。

2) 多個飛行器不同起始點同時開始滑翔,同一時間到達同一目標點。

步驟1分別優化得到各個起點到終端目標點的滑翔最長時間tmax,i以及最短時間tmin,i。

步驟2多飛行器的協同到達時間t可依據式(15)確定,以使得多飛行器同時到達目標點。

t∈[tmin,1,tmax,1]∩…∩[tmin,i,tmax,i]

(15)

在本文算例中,不同起點的飛行器滑翔段協同時間t按式(16)計算：

(16)

步驟3將協同時間指令下發至各子飛行器,各子飛行器進行固定滑翔段總時間約束下的軌跡協同規劃。

2.1.2 協同時間約束下的軌跡優化模型

針對以上任務場景需求,將協同軌跡規劃問題轉化為2個子優化問題,即P1：再入滑翔段飛行最長、最短總時間優化問題；P2：固定總時間約束下的軌跡規劃問題。結合各類約束條件,滑翔段軌跡優化問題求解模型可以被初步描述為

(17)

式中：tf為滑翔段終端時間；δtf為終端時間相對于期望協同到達時間的偏差；R0為地球半徑。

2.2 考慮時間約束的預設剖面設計

將再入過程分為初始下降段及滑翔段,初始下降段傾側角設置為常值。定義剩余航程s為飛行器當前位置到目標點連線在地球表面投影的大圓弧長度,其表達式為

s=R0arccos[sinφsinφf+cosφcosφfcos(λf-λ)]

(18)

對時間求導可得

(19)

式中：Δψ為航向角誤差。由e=1/r-V2/2可得

(20)

則無量綱化剩余航程s相對于能量e的變化率為

(21)

由式(21)可知,飛行器在縱向平面內航程s的變化主要由飛行路徑角γ、地心距r和阻力加速度D所決定。由于以能量為離散自變量進行軌跡設計,各個離散點處e的是固定的,且r近似為1,因此速度v和阻力加速度D的變化形式相對固定,s的大小主要取決于飛行路徑角γ的剖面形式[23]。仿真實驗表明,極值總時間或滿足期望時間約束對應的軌跡有時并不唯一,為提高算法收斂性,并且使滑翔段軌跡規劃結果盡量平滑,將飛行路徑角γ的初值剖面定義為二次曲線的形式：

(22)

式中：γ0、γf為初始點及目標點處的飛行路徑角,則當中間能量狀態處的飛行路徑角γmid確定后,即可確定飛行路徑角全程剖面,積分式(21)得到總航程s。迭代求取γmid,當式(23)成立時,可認為飛行路徑角初值剖面滿足終端航程要求。

(23)

2.3 時間協同軌跡優化問題的凸化處理

為應用SOCP算法,需對目標函數中的時間項進行凸化處理。此外,將P0中的連續非線性動力學方程線性化,轉化為離散型有限維最優控制問題進行求解。

2.3.1 時間項的凸化

將飛行時間轉化為模型狀態量的線性表達式,再入段末端時刻tf可表示為

(24)

(25)

將式(2)及的表達式代入式(25)可得

(26)

式中：ρr為大氣密度ρ相對于地心距r的導數;Tr為T相對于地心距r的導數。

2.3.2 動力學約束的凸化

按式(7)重新選取控制量后,系統動力學模型可整理為線性的凸約束形式：

(27)

式中:w(x,e)為與地球旋轉有關的項,

(28)

(29)

(30)

2.3.3 問題的離散化

軌跡規劃問題為連續型無窮維非線性規劃問題,求解時需要對其進行離散化。本文以能量為自變量,采用梯形離散化方法,在滑翔段能量變化區間[e0,ef]內取均勻分布的N+1個離散點{e0,e1,e2,…,eN},相鄰兩個離散點的距離Δe滿足:

Δe=(ef-e0)/N

(31)

則各個離散點對應的能量ei為

ei=e0+i·Δei=0,1,…,N

(32)

離散后的狀態量和控制量可分別表示為(x0,x1,…,xN)和{u0,u1,…,uN}。第i個離散點處控制量為ui=[cosσi,sinσi],第i個離散點處附近對應的動力學約束可以被線性化為

(Aixi+Βiui+bi)]=0i=1,2,…,N

(33)

(34)

式中：ri為第i個離散點處的地心距;Tri,τi為線性化系數,由第k-1次序列迭代的解代入得到。

2.4 固定時間約束下的序列凸化算法

若僅以式(34)作為目標函數,則會導致控制量求解結果不滿足約束條件[14],因此引入復合型目標函數,將終端經緯度偏差及時間項加入罰函數,提高問題求解的可行性。此外,含時間約束的優化問題可能存在多個解,致使序列迭代難以收斂,通過罰函數及信賴域系數的自適應調整避免該問題,提高算法的收斂性。

2.4.1 終端約束的處理

考慮到終端約束過于嚴格會導致沒有可行解,適當放寬終端經緯度的限制,在終端位置約束不能滿足的情況下接受最接近目標點的解。定義誤差變量λe,φe,則每一次迭代求解結果中終端經緯度均滿足：

|λ(ef)-λf|≤λe, |φ(ef)-φf|≤φe

(35)

進一步將λe,φe加入目標函數,使求解得到的軌跡終點與目標點經緯度偏差最小:

J=pλe+qλφ

(36)

式中:p、q為自適應系數,取值由2.4.3節給出。

2.4.2 目標函數中時間項的處理

針對問題P1,直接將線性化后的時間項加入到目標函數中,即

(37)

針對問題P2,在求解模型中加入期望總飛行時間約束：

(38)

式中：Tri、τi第i個離散點處的時間項系數；Tc為協同時間；ΔT為當前迭代過程中的時間誤差,將ΔT加入到目標函數中,即

J=mΔT

(39)

式中:m為自適應調整的系數,其取值策略由2.4.3 節給出。當ΔT趨于0時,序列凸化求解得到的軌跡將滿足協同時間的要求。

2.4.3 罰函數系數及信賴域自適應調整策略

在實際求解過程中,滑翔段最長、最短總時間及固定總時間等指標對應的軌跡形式并不唯一,存在多個解對應的目標函數值相同的情況,從而導致在序列迭代過程的末段,求解得到的結果在2條不同的軌跡形式間不斷振蕩,難以收斂。因此本文采用一種罰函數系數及信賴域自適應調整策略,即將飛行路徑角相對于預設平滑剖面的偏離程度加入罰函數,而后不斷降低其在罰函數中的權值系數,同時收縮信賴域,以保證求解過程能夠收斂,并且得到的軌跡較為平滑。飛行路徑角剖面搜索范圍的變化趨勢如圖3所示。

圖3 飛行路徑角剖面搜索范圍變化示意圖

引入式(22)對應的飛行路徑角預設剖面作為求解過程中的軟約束,即在第i個離散點處飛行路徑角γi需滿足：

|γi-γ0,i|≤Δγi,i=1,2,…,N

(40)

在序列迭代的開始階段,取較大的n值,此時將在圖2中的“Ⅰ”區域內搜索可行飛行路徑角剖面,求解的主要指標為軌跡的平滑性。隨著序列迭代的進行,逐步降低n值,飛行路徑角預設剖面對軌跡求解的影響逐漸減弱,搜索區域將擴大至區域“Ⅱ”,有利于在平滑的初始軌跡附近搜索得到最優解。則問題P1、P2的目標函數可表示為

(41)

式中:p、q、m、n為自適應調整的系數,第k+1次迭代求解時各系數表達式為

(42)

由于動力學方程及各類約束的線性化過程基于小擾動假設,因此只有待求解變量在原參考點附近取值,即相鄰兩次求解結果滿足式(43)所示的信賴域約束時,離散后的求解模型才是對原非線性問題的有效近似。若信賴域半徑過小,則求解過程中搜索的范圍有限,可能收斂不到最優解；若信賴域半徑過大,雖然易于找到最優解,但對于存在多個最優解的情況,序列凸化算法難以收斂,因此采用一種信賴域自適應收縮策略,定義第次迭代求解時的信賴域約束為

|x(k+1)-x(k)|≤η(k)δ

(43)

式中:x=[r,λ,φ,γ,ψ]T;δ為各狀態量的初始信賴域[δr,δλ,δφ,δγ,δψ]T;η(k)為隨序列迭代次數改變的系數,其表達式為

(44)

式中:ω=1-|(J(k-1)-J(k-2))/J(k-2)|;ε為預先設定的閾值,0<ε<1/2;η(k)的選取由收斂系數ξ及ω共同決定。系數ω反映了代價函數值J的變化,當ω小于閾值ε時,認為代價函數J的值下降不明顯,解已經趨于最優解附近,為加快收斂速度,避免因存在多個最優解導致求解結果在多條軌跡間振蕩的情況,引入收斂系數ξ,滿足0<ξ<1。ξ值過大或過小會導致軌跡精度較差或收斂較慢,因此本文取ξ=0.5；當ω大于閾值ε時,令ξ=1。

在迭代初始階段,信賴域半徑δ可選取為一個足夠大的數值,以保證能夠在較大的初始解空間范圍內搜索得到近似最優解。在后續迭代過程中,隨著信賴域半徑的減小,最優解的搜索空間逐漸減小,易于序列凸化算法的收斂,相比于固定信賴域半徑方法的求解結果更接近全局最優解。

2.5 序列迭代求解過程

在對問題凸化和離散化之后,多高超聲速飛行器滑翔段時間協同彈道規劃問題可被分為P1、P2兩個序列二階錐規劃問題,兩個子問題的求解模型可表示為

(45)

2) 固定總時間約束下的軌跡規劃問題

(46)

本文所采用的序列二階錐規劃算法流程如圖4所示,具體的求解步驟如下。

圖4 序列凸化算法流程圖

步驟1令迭代次數k=0,迭代求取飛行路徑角初值剖面,獲取軌跡的初始解x(0),并將x(0)定義為參考軌跡x(k)。

步驟2對動力學方程及各類約束進行凸化、離散化,將參考軌跡代入計算動力學約束中的系數矩陣A(x(k))、B(x(k))、b(x(1))等,得到SOCP形式的問題求解模型。

步驟3利用SOCP算法得到最優解x(k+1),更新參考軌跡。

步驟4判斷|x(k+1)-x(k)|<ε是否滿足,其中ε為預先設定的常值向量,若收斂條件不滿足,則更新罰函數系數m、n、p、q及信賴域半徑η(k-1)δ,轉步驟2,若收斂條件滿足,則轉步驟5。

步驟5得到優化問題的最優解為x(k+1),序列迭代結束。

3 仿真結果與分析

為驗證所提滑翔段時間協同軌跡規劃算法的有效性,以3個滑翔飛行器構成的集群為例,基于2種協同任務場景進行了數值仿真實驗。首先,針對問題P1,分別運用改進序列凸化算法與基于非線性規劃求解器GPOPS的h-p自適應Radau偽譜法,求解了最長、最短再入時間對應的軌跡。而后在確定協同時間的基礎上求解了問題P2,分析了含時間約束的優化問題求解過程中,本文算法與傳統序列凸化算法的收斂效果。最后對比了序列凸化算法與偽譜法的軌跡規劃結果。

3.1 多飛行器不同起始點協同再入

表1給出了各飛行器及目標點的初始狀態信息,表2給出了不同算法對于飛行極值時間的求解結果。

表1 各飛行器及目標點初始狀態

表2 序列凸化(SOCP)與偽譜法(RPM)結果對比

由表2可知,各飛行器的再入總時間可調范圍分別為144.5 s,269.6 s,239.6 s。由各飛行器與目標點的相對位置可以看出,初始航向角對飛行器的再入總時間極值有著較大影響,飛行器2、3相比于飛行器1的可調范圍顯著增大,當初始速度方向指向目標時,飛行器的再入總時間將有更大的可調范圍。SOCP算法得到的總時間可調范圍與偽譜法基本一致,偽譜法求得的最短時間稍小于SOCP算法,總時間可調范圍稍大,但SOCP算法的計算時間要明顯優于偽譜法。根據時間協同策略可確定3個飛行器的協同再入總時間為1 430 s,利用本文提出的固定期望總時間約束下的軌跡規劃方法求解再入協同軌跡,飛行器1～3的三維協同軌跡與其各自的最速下降軌跡地面投影對比如圖5和圖6 所示。表3 給出了各飛行器軌跡數值積分驗證結果的終端狀態精度情況。

圖5 協同軌跡三維圖

圖6 協同軌跡地面投影

表3 協同軌跡積分驗證結果

由圖5 可知,各飛行器在到達時間相同的基礎上,能夠實現終端飛行路徑角及航向角的一致。圖6 中各飛行器為實現最短時間再入,速度方向基本保持在起始點與目標點連線的縱向平面內,軌跡的地面投影均較為平直。為實現再入時間的一致,各飛行器在初始階段偏離了最速下降方向,通過側向機動的方式延緩到達時間,從而實現時間協同。根據表3中的積分驗證結果,終端最大高度誤差為140 m,最大再入總時間誤差為6.28 s,終端飛行路徑角與航向角誤差均小于1°,滿足精度要求。

3.2 多飛行器同起點序列式再入

假定3個滑翔飛行器從同起點再入,在滑翔段起始點分離,要求序列到達同一目標點,驗證了算法的有效性。表4 給出了仿真初始條件,圖7～圖9給出了各飛行器軌跡對比及熱流、動壓、過載約束滿足情況,圖10以問題P1的求解過程為例,展示了本文方法對于序列凸化算法收斂性的改進情況。

圖7 序列到達軌跡高度-能量剖面對比

圖8 序列到達軌跡地面投影對比

圖9 滑翔段動壓、過載、熱流變化情況

表4 仿真初始條件

圖7和圖8表明最長時間到達軌跡呈周期性

進一步將本文變系數方法與傳統序列凸化方法、固定罰函數及信賴域收縮系數方法的迭代求解過程進行對比。圖10給出了本文方法對于求解最短時間再入問題收斂過程的影響情況。圖中,縱軸為各方法在第次迭代過程中求解得到的滑翔段總時間相對于最優解的絕對誤差。可以看出,本文變系數方法最先收斂至最優解；傳統序列凸化方法求解結果不斷振蕩,說明存在多組解對應的目標函數值相同,使得算法難以收斂；罰函數

圖10 不同方法最短再入時間求解結果的絕對誤差隨迭代次數變化情況

及信賴域系數均固定的方法在開始階段目標函數值收斂趨勢與變系數方法一致,但在第3次迭代后偏離最優解,說明目標函數中的時間項受到了其他指標影響,未在序列優化后期起主要作用。

由表5可得,本文方法與偽譜法優化得到的到達時間基本一致,終端狀態量精度均滿足要求。傳統序列凸化方法未收斂至最優解,導致終端存在500 m高度偏差,最短再入總時間同樣出現較大誤差。接下來驗證了本文方法在處理給定期望總時間約束下的軌跡求解的有效性,并與偽譜法、MPSP方法進行對比。

表5 不同方法終端精度對比

選取1 500 s作為期望總時間,分別進行給定總時間約束條件下的軌跡解算,圖11和圖12給出了3種方法的縱向剖面及地面軌跡投影的結果對比。其中,序列凸化方法總用時4.24 s,MPSP方法用時10.82 s,偽譜法用時105.02 s。3種方法求解結果都能夠滿足期望總時間約束,終端精度均處于合理范圍。由規劃結果可知,同一滑翔段總時間對應的軌跡形式并不唯一,三者優化得到的地面軌跡基本一致,但序列凸化方法計算得到的縱向軌跡相比于偽譜法及MPSP方法上下跳躍幅度明顯減小,軌跡更為平滑,且解算用時更短。此外,相比于序列凸化及偽譜法,MPSP方法的收斂性對于控制量初值選取較為敏感。

圖11 高度-能量剖面對比

圖12 地面投影軌跡對比

3.3 算法適應性測試

進一步,分別在初始點狀態大范圍變化及小范圍攝動2種情況下,對算法的收斂性及快速性進行仿真測試。首先,選取期望再入總時間為1 400 s,目標終端狀態為

xf=[hf,vf,λf,φf]=[20 km,2 km/s,155°,50°]

選取多個不同的起始點狀態,以驗證算法的收斂性和快速性。各次仿真的初始狀態及算法各次求解用時情況如表6所示。圖13和圖14給出了各次求解得到的軌跡積分驗證結果。

圖13 三維軌跡

表6 仿真初始條件

各次仿真優化求解的平均時間為3.64 s,最大總時間偏差為7 s ,最大終端高度偏差為516 m,最大終端速度偏差為62 m/s,最大終端經度偏差為0.63°,最大終端緯度偏差為0.46°。

最后,通過蒙特卡洛方法測試所提固定總時間約束下軌跡規劃算法對于初始經緯度變化的適應性。目標終端狀態為

xf=[20 km,2 km/s,155°,

50°,-0.5°±1°,90°±1°]

初始點狀態的取值上限定義為

x0max=[50 100 m,6 050 m/s,90°,45°,0°,90°]

初始點狀態的取值下限定義為

x0min=[49 900 m,5 980 m/s,88°,50°,-1°,90°]

在上述初始點范圍內隨機取值并利用序列凸化算法求解滿足再入總時間約束的可行軌跡。共進行50次仿真,各次求解得到的軌跡通過數值積分驗證后的結果如圖15～圖17所示。

圖15 蒙特卡洛仿真結果

圖16 終端高度分布

圖17 再入總時間分布

由仿真結果統計可得,序列凸化算法各次仿真平均CPU用時為3.92 s；終端高度平均偏差為136.05 m,最大偏差為214.59 m；再入總時間平均偏差為2.37 s,最大時間偏差為8.04 s；平均經緯度偏差為0.12°,終端航向角及飛行路徑角均滿足誤差要求。

4 結 論

1) 針對多高超聲速滑翔飛行器集群協同再入問題,設計了“探測-打擊-評估”一體化模式下,集群內部時間協同的軌跡規劃方案。

2) 給出了集群協同時間的確定策略,解決了給定期望再入協同時間約束下的多飛行器軌跡快速規劃問題。

3) 提出了一種罰函數與信賴域自適應調整策略,能夠避免存在多個最優解情況下的軌跡振蕩問題,提高了序列凸化算法的收斂性。

4) 仿真結果表明,所提改進序列凸化算法能夠適應不同的時間協同再入任務,收斂性受初值影響相對MPSP算法較小,求解結果相比于偽譜法軌跡更加平滑,計算速度優于偽譜法。