周宏宇，王小剛，趙亞麗，崔乃剛

(1. 哈爾濱工業(yè)大學(xué)航天學(xué)院，哈爾濱 150001；2. 北京航天晨信科技有限公司，北京 102308)

0 引 言

現(xiàn)代航天活動(dòng)正朝著多元化、頻繁化的方向發(fā)展；在這種背景下，傳統(tǒng)的一次性航天運(yùn)載工具面臨著發(fā)射成本高、安全性低、準(zhǔn)備時(shí)間長(zhǎng)等問(wèn)題，制約了航天事業(yè)的發(fā)展。因此，世界各航天強(qiáng)國(guó)始終在研究更先進(jìn)、更便捷的航天運(yùn)載器，其中最具代表性的便是基于吸氣式組合動(dòng)力的水平起降可重復(fù)使用運(yùn)載器[1-3]。20世紀(jì)末，日本率先提出了一種將變循環(huán)渦扇發(fā)動(dòng)機(jī)和沖壓發(fā)動(dòng)機(jī)組合的方案，證明了其用于天地往返運(yùn)輸?shù)目尚行浴?016年，美國(guó)公布了基于“協(xié)同吸氣式火箭發(fā)動(dòng)機(jī)”的兩級(jí)入軌方案；同年，中國(guó)披露了正在研制的組合動(dòng)力運(yùn)載器項(xiàng)目，計(jì)劃在2030年實(shí)現(xiàn)亞軌道應(yīng)用；在2017全球航天探索大會(huì)上，中國(guó)宣布正在研發(fā)水平起降天地往返運(yùn)載器，現(xiàn)已完成多項(xiàng)關(guān)鍵技術(shù)的地面試驗(yàn)。世界各國(guó)的諸多舉動(dòng)表明，基于吸氣式組合動(dòng)力的水平起飛兩級(jí)入軌方案已成為發(fā)展可重復(fù)使用運(yùn)載器的重要方向[4-6]。

設(shè)計(jì)上升段軌跡是可重復(fù)使用運(yùn)載器執(zhí)行任務(wù)的首要步驟。但對(duì)于采用組合動(dòng)力完成水平起飛、爬升加速的運(yùn)載器來(lái)說(shuō)，其動(dòng)力、氣動(dòng)、軌跡和指標(biāo)間有著很強(qiáng)的相互制約關(guān)系，并且不同動(dòng)力模態(tài)對(duì)飛行狀態(tài)和環(huán)境約束有著不同的需求[7]。因此，上升段軌跡優(yōu)化面臨著約束條件嚴(yán)苛、動(dòng)力切換頻繁、性能指標(biāo)與設(shè)計(jì)變量間映射關(guān)系復(fù)雜等問(wèn)題，求解上升段軌跡優(yōu)化問(wèn)題面臨著極大的挑戰(zhàn)。多年來(lái)，許多學(xué)者針對(duì)水平起飛運(yùn)載器的上升段軌跡優(yōu)化問(wèn)題展開(kāi)了研究。早期主要是借助簡(jiǎn)化數(shù)學(xué)模型解決這類(lèi)優(yōu)化問(wèn)題，例如Corban等[8]利用近似線性模型簡(jiǎn)化了高度和速度的微分方程，降低了軌跡設(shè)計(jì)的難度。隨著最優(yōu)控制理論的發(fā)展，開(kāi)始考慮更符合客觀規(guī)律的數(shù)學(xué)模型。Lu等[9]指出，吸氣式動(dòng)力運(yùn)載器的爬升軌跡優(yōu)化問(wèn)題可轉(zhuǎn)化為典型兩點(diǎn)邊值問(wèn)題求解。Murillo等[10]將優(yōu)化問(wèn)題轉(zhuǎn)化為一個(gè)多項(xiàng)式方程求根問(wèn)題，進(jìn)而提出了一種上升段軌跡快速優(yōu)化方法。Pan等[11]通過(guò)簡(jiǎn)化協(xié)態(tài)方程降低了間接法的求解難度，同時(shí)借助多重打靶法實(shí)現(xiàn)了初值的快速生成。陶超[12]用Gauss偽譜法求解了高超聲速飛行器最優(yōu)上升段軌跡。龔春林等[13]針對(duì)吸氣式動(dòng)力工作模態(tài)多、沖壓模態(tài)約束條件復(fù)雜的問(wèn)題，利用Radau偽譜法建立了基于混合積分變量的軌跡優(yōu)化模型，有效處理了運(yùn)載器多段飛行約束條件。竇立謙等[14]將Gauss偽譜法和(Sequence quadratic program)SQP相結(jié)合，求解了考慮密度、聲速和發(fā)動(dòng)機(jī)推力變化及地球引力等因素的燃料最優(yōu)爬升軌跡。李惠峰等[15]基于有限差分方法和改進(jìn)牛頓法，設(shè)計(jì)了一種以參考面積為同倫參數(shù)的迭代方法，同時(shí)以近似解析解為初值，實(shí)現(xiàn)了吸氣式運(yùn)載器上升軌跡的快速優(yōu)化。

現(xiàn)有文獻(xiàn)有效解決了組合動(dòng)力運(yùn)載器的上升段軌跡優(yōu)化問(wèn)題，但這些研究通常只關(guān)注臨近空間中的高超聲速爬升/巡航段，鮮有關(guān)注從起飛至進(jìn)入稀薄大氣的完整爬升軌跡。當(dāng)針對(duì)跨大空域/速域的完整上升段軌跡進(jìn)行設(shè)計(jì)時(shí)，必須考慮整合多種動(dòng)力模態(tài)的多組合動(dòng)力系統(tǒng)(如渦輪基沖壓發(fā)動(dòng)機(jī)(TBCC)+火箭基沖壓發(fā)動(dòng)機(jī)(RBCC))，涉及更多的軌跡分段、約束條件、耦合因素和優(yōu)化變量，優(yōu)化問(wèn)題的非線性將更強(qiáng)；而現(xiàn)有文獻(xiàn)難以在滿(mǎn)足諸多約束條件下獲得最佳的多組合動(dòng)力使用方式和上升段軌跡。

為解決上述問(wèn)題，針對(duì)上升段中動(dòng)力、狀態(tài)和指標(biāo)間的復(fù)雜耦合關(guān)系以及動(dòng)力系統(tǒng)對(duì)飛行環(huán)境提出的復(fù)雜約束條件，設(shè)計(jì)了一種全新的考慮多重動(dòng)力切換的攻角剖面，并將上升段軌跡優(yōu)化問(wèn)題轉(zhuǎn)化為多參數(shù)尋優(yōu)問(wèn)題。基于該剖面進(jìn)行推導(dǎo)，得到了攻角及其變化率的極值與優(yōu)化參數(shù)間的解析映射關(guān)系，大幅簡(jiǎn)化了優(yōu)化算法需要處理的約束條件，同時(shí)自然滿(mǎn)足軌跡的平滑性要求。然后，提出了一種改進(jìn)的PSO算法用于求解最優(yōu)軌跡：以粒子收斂性分析為基礎(chǔ)，以粒子空間散布狀態(tài)和收斂時(shí)間為參考，以燃料最優(yōu)為指標(biāo)，引入DDPG(Deep deterministic policy gradient)強(qiáng)化學(xué)習(xí)方法對(duì)影響PSO計(jì)算效率的關(guān)鍵參數(shù)進(jìn)行動(dòng)態(tài)調(diào)整；隨著PSO迭代信息的積累和粒子群向最優(yōu)解的逐漸聚集，DDPG能夠智能、準(zhǔn)確、快速地輸出PSO算法控制參數(shù)，進(jìn)而提高算法計(jì)算性能。仿真結(jié)果表明，本文提出的優(yōu)化方法能夠有效解決運(yùn)載器上升段軌跡優(yōu)化問(wèn)題；相對(duì)其它PSO算法，結(jié)合強(qiáng)化學(xué)習(xí)的改進(jìn)PSO算法在計(jì)算效率上有明顯提升。

1 組合動(dòng)力運(yùn)載器運(yùn)動(dòng)數(shù)學(xué)模型

1.1 總體參數(shù)模型

借鑒日本“HYPR”計(jì)劃中提出的“將變循環(huán)渦扇發(fā)動(dòng)機(jī)和沖壓發(fā)動(dòng)機(jī)組合”方案，本文研究的運(yùn)載器先后在TBCC和RBCC的作用下飛行。通過(guò)分析TBCC和RBCC的固有特點(diǎn)和工作模態(tài)，可根據(jù)熱力學(xué)原理建立組合動(dòng)力發(fā)動(dòng)機(jī)性能計(jì)算模型。

1)TBCC性能計(jì)算

TBCC動(dòng)力包含渦輪和沖壓兩種模態(tài)。采用渦輪風(fēng)扇發(fā)動(dòng)機(jī)作為T(mén)BCC的渦輪部分，同時(shí)采用“設(shè)計(jì)點(diǎn)分段定比熱氣動(dòng)熱力計(jì)算”方式，建立渦輪發(fā)動(dòng)機(jī)性能計(jì)算模型；結(jié)合氣流路徑和部件間能量守恒關(guān)系，可得渦輪發(fā)動(dòng)機(jī)的推力和比沖為[16]：

(1)

式中：Wa,tb為進(jìn)入渦輪發(fā)動(dòng)機(jī)的空氣質(zhì)量流量，F(xiàn)s,tb為推力/流量比，msfs為發(fā)動(dòng)機(jī)耗油率，g為重力加速度。Wa,tb=ρc0Stb，其中ρ為環(huán)境大氣密度，Stb為渦輪發(fā)動(dòng)機(jī)進(jìn)氣道入口面積，c0為進(jìn)氣道入口氣流速度；Fs,tb和msfs的計(jì)算方式見(jiàn)文獻(xiàn)[16]。

考慮壓縮、加熱、絕熱膨脹和放熱四個(gè)過(guò)程，TBCC沖壓發(fā)動(dòng)機(jī)的推力和比沖分別為[16]：

(2)

式中：ηf為沖壓發(fā)動(dòng)機(jī)總效率，Hu為燃油熱值，fsc為油氣比，Wa,sc為進(jìn)入沖壓發(fā)動(dòng)機(jī)的空氣質(zhì)量流量。類(lèi)似的，有Wa,sc=ρc0Ssc，其中Ssc為沖壓發(fā)動(dòng)機(jī)進(jìn)氣道入口面積。

2)RBCC性能計(jì)算

設(shè)計(jì)RBCC發(fā)動(dòng)機(jī)的性能為馬赫數(shù)的函數(shù)。一般來(lái)說(shuō)，當(dāng)切換至RBCC時(shí)馬赫數(shù)已超過(guò)3.0，因此僅關(guān)注RBCC亞燃沖壓及以后的工作模態(tài)。

亞燃/超燃沖壓模態(tài)下的推力和比沖分別為：

(3)

(4)

式中：m0為運(yùn)載器起始質(zhì)量，Ma為馬赫數(shù)。

純火箭模態(tài)下的推力和比沖分別為[16]：

(5)

綜上，動(dòng)力系統(tǒng)的待優(yōu)化參數(shù)如表1所示[16]。

3)氣動(dòng)模型

認(rèn)為氣動(dòng)系數(shù)受攻角和馬赫數(shù)的影響，通過(guò)相關(guān)氣動(dòng)參數(shù)擬合，同時(shí)參考類(lèi)似高超聲速飛行器的氣動(dòng)特性，給出如下運(yùn)載器的氣動(dòng)參數(shù)模型：

(6)

式中：cD為阻力系數(shù)，cL為升力系數(shù)，α為攻角。

表1 組合動(dòng)力發(fā)動(dòng)機(jī)待優(yōu)化參數(shù)Table 1 Optimization parameters for combined-cycle-basedengine

1.2 運(yùn)載器上升段運(yùn)動(dòng)模型

縱向剖面內(nèi)的運(yùn)載器上升段運(yùn)動(dòng)模型為[14]：

(7)

式中：V為飛行速度，γ為飛行路徑角，r為運(yùn)載器質(zhì)心到地心的距離，m為運(yùn)載器質(zhì)量，Tm為推力，Isp為比沖，D為氣動(dòng)阻力，L為氣動(dòng)升力，g為重力加速度。

氣動(dòng)力的計(jì)算方式如下：

(8)

式中：q=ρV2/2為飛行動(dòng)壓，Sm為特征面積。

1.3 上升段軌跡優(yōu)化模型

首先，設(shè)置終端約束條件如下：

(9)

考慮結(jié)構(gòu)強(qiáng)度和防熱等因素，取過(guò)程約束為：

(10)

式中：qmax為允許的最大動(dòng)壓，Nm為法向過(guò)載，Nm,max為允許的最大法向過(guò)載。

另外，不同于采用火箭發(fā)動(dòng)機(jī)，吸氣式發(fā)動(dòng)機(jī)對(duì)飛行攻角有著很?chē)?yán)格的約束需求，并且組合動(dòng)力中不同模態(tài)對(duì)飛行攻角的約束也各不相同。因此，針對(duì)各動(dòng)力模態(tài)，分別建立不同的攻角約束：

(11)

式中：αmax,1和αmin,1分別為T(mén)BCC渦輪工作狀態(tài)下的攻角上/下限，αmax,2和αmin,2分別為T(mén)BCC沖壓工作狀態(tài)下的攻角上/下限，αmax,3和αmin,3分別為RBCC亞/超燃沖壓狀態(tài)下的攻角上/下限，αmax,4和αmin,4分別為RBCC火箭模態(tài)下的攻角上/下限；Maq3為RBCC進(jìn)入火箭模態(tài)時(shí)的馬赫數(shù)，|dα/dt|max為允許的最大攻角變化率。

最后，將上升段軌跡優(yōu)化問(wèn)題的性能指標(biāo)取為燃料最省，即minJ=(m0-mf)/m0。

2 上升段攻角剖面設(shè)計(jì)

通過(guò)設(shè)計(jì)攻角—速度剖面完成上升段軌跡設(shè)計(jì)；根據(jù)動(dòng)力形式(TBCC渦輪、TBCC沖壓和RBCC)將攻角剖面分為三段，即將攻角設(shè)計(jì)為速度的分段多項(xiàng)式函數(shù)：

α(V)=

(12)

式中：V0為初始速度，Vq1和Vq2分別對(duì)應(yīng)Maq1和Maq2；ai(i=0,1,…,n1)，bi(i=0,1,…,n2)和ci(i=0,1,…,n3)為待定系數(shù)，共(n1+n2+n3+3)個(gè)。

取n1=n2=n3=3，即用三次多項(xiàng)式描述各動(dòng)力段的攻角，因此共12個(gè)未知系數(shù)。相應(yīng)的攻角剖面示意圖見(jiàn)如圖1，其中V1,0=V0，V2,0=Vq1，V3,0=Vq2，V4,0=VE=Vf。

為求解這12個(gè)未知系數(shù)，需要建立12個(gè)方程。首先，為保證攻角及其變化率在Vq1和Vq2處的連續(xù)性，需滿(mǎn)足以下條件：

(13)

(14)

(15)

(16)

式中：V1,1=V1,0+(Vq1-V1,0)/3，V1,2=Vq1-(Vq1-V1,0)/3，V2,1=(Vq1+Vq2)/2，V3,1=(Vf+Vq2)/2。

圖1 攻角-速度剖面Fig.1 Ascent profile of angle-of-attack (AOA) and velocity

因此，上升段剖面的設(shè)計(jì)變量為：

U1=[α1,1,α1,2,α2,0,α2,1,α3,0,α3,1,α3,3]T

(17)

同時(shí)，發(fā)動(dòng)機(jī)待優(yōu)化參數(shù)為：

U2=[πCL,πCH,T3,Maq1,Maq2]T

(18)

由優(yōu)化約束模型可知，上升段對(duì)攻角提出的約束最為復(fù)雜，若能將攻角約束在優(yōu)化計(jì)算前解決，則能夠大大減少優(yōu)化算法的計(jì)算量。

對(duì)于三階多項(xiàng)式，以區(qū)間[V0,Vq1]為例，攻角的表達(dá)式為：

(19)

對(duì)式(19)求導(dǎo)可得：

(20)

因此：

(21)

(22)

由式(21)可得：

(23)

(24)

(25)

式中：εh是相對(duì)hE而言很小的量。

3 上升段攻角剖面優(yōu)化算法

3.1 PSO算法收斂性分析

PSO算法中每個(gè)粒子都代表優(yōu)化問(wèn)題的一個(gè)潛在解，對(duì)應(yīng)的性能指標(biāo)稱(chēng)為“適應(yīng)度”。設(shè)Mpso個(gè)粒子在空間中搜索Dpso個(gè)參數(shù)，則第i(i=1,2,…,Mpso)個(gè)粒子在Dpso維搜索空間中的位置矢量和速度矢量分別為：

(26)

記第i個(gè)粒子迄今找到的最佳位置為Pi，Pi對(duì)應(yīng)的適應(yīng)度為pbest；同時(shí)記所有粒子搜索到的最佳位置為Pg，Pg對(duì)應(yīng)的適應(yīng)度為gbest。

(27)

PSO算法中，整個(gè)群體按照如下方式進(jìn)化[17]：

(28)

(29)

定義φ1=c1r1+c2r2，φ2=c1r1pi,j+c2r2pg,j，由PSO算法進(jìn)化公式可得：

(30)

設(shè)pi,j和pg,j為常值，分別記為p和g。用符號(hào)E表示某變量的期望，則根據(jù)式(30)可得：

(31)

(32)

當(dāng)φ1-w-1>0時(shí)，粒子在第j維空間中將漸進(jìn)收斂于(p+g)/2。設(shè)算法在相對(duì)誤差小于2%時(shí)收斂，記收斂時(shí)間為ts，則有：

1)當(dāng)0<ξn<1時(shí)，式(30)是一個(gè)欠阻尼系統(tǒng)，收斂時(shí)間為：

(33)

2)由于φ1具有隨機(jī)性，故ξn=1一般不會(huì)出現(xiàn)。

3)當(dāng)ξn>1時(shí)，式(30)是一個(gè)過(guò)阻尼系統(tǒng)，收斂時(shí)間為：

(34)

PSO算法中通常取c1=c2=2.0，因此?E(φ1)=2，其中?E(·)表示某變量的期望；取φ1=2，結(jié)合式(33)～(34)并用w取代ξn和ωn，得到收斂時(shí)間與慣性權(quán)重的關(guān)系如圖2所示。在此指出，粒子運(yùn)動(dòng)是以迭代次數(shù)來(lái)描述的，故圖2中顯示的并非實(shí)際時(shí)間，僅為收斂性分析提供參考。

圖2 收斂時(shí)間與慣性權(quán)重的關(guān)系Fig.2 The relation between convergence time and weight

3.2 基于強(qiáng)化學(xué)習(xí)的改進(jìn)PSO算法

由圖2和式(33)可知，粒子的收斂過(guò)程主要由頻率和阻尼決定，而頻率和阻尼又由慣性權(quán)重、隨機(jī)數(shù)r1和r2以及學(xué)習(xí)因子決定。因此，只需根據(jù)r1和r2動(dòng)態(tài)調(diào)整w，c1和c2，即可實(shí)現(xiàn)對(duì)PSO收斂過(guò)程進(jìn)行控制。然而，以線性權(quán)重、隨機(jī)權(quán)重、指數(shù)權(quán)重和固定權(quán)重為代表的諸多現(xiàn)有研究成果表明，依靠慣性權(quán)重對(duì)PSO算法的搜索過(guò)程進(jìn)行準(zhǔn)確控制是很難實(shí)現(xiàn)的，這不僅需要設(shè)計(jì)者對(duì)研究問(wèn)題的深入了解，還依賴(lài)于大量的仿真實(shí)驗(yàn)分析[18]。

設(shè)Ulearn=[w,c1,c2]T，為了獲得最佳Ulearn，借助強(qiáng)化學(xué)習(xí)方法[19]，通過(guò)構(gòu)建策略網(wǎng)絡(luò)對(duì)PSO進(jìn)化過(guò)程的智能控制，即在每一代進(jìn)化中輸出最優(yōu)Ulearn值，進(jìn)而實(shí)現(xiàn)對(duì)PSO算法收斂過(guò)程的動(dòng)態(tài)智能控制，提高PSO算法性能。

強(qiáng)化學(xué)習(xí)的基本過(guò)程如下：

1)環(huán)境給智能體一個(gè)觀測(cè)值。

2)智能體收到觀測(cè)值后做出一個(gè)動(dòng)作。

3)環(huán)境收到動(dòng)作后做出一系列反應(yīng)，如對(duì)智能體的動(dòng)作給予一個(gè)獎(jiǎng)勵(lì)或給出新的觀測(cè)值。

在PSO算法中，觀測(cè)值定義為表征粒子空間散布狀態(tài)的向量；動(dòng)作空間為3維，分別對(duì)應(yīng)對(duì)w，c1和c2的更新，且動(dòng)作為連續(xù)值。

(35)

式中：n1，n2和n3為事先設(shè)定的超參數(shù)。

2)獎(jiǎng)勵(lì)(記為r)：獎(jiǎng)勵(lì)為一個(gè)反饋信號(hào)，表明當(dāng)前動(dòng)作的優(yōu)劣；取gbest作為獎(jiǎng)勵(lì)值。

3)狀態(tài)(記為s)：環(huán)境狀態(tài)轉(zhuǎn)移根據(jù)式(30)完成，每次進(jìn)化后重新計(jì)算粒子的空間散布狀態(tài)信息。

基于DDPG強(qiáng)化學(xué)習(xí)算法調(diào)整Ulearn。DDPG框架如下：

1)值函數(shù)網(wǎng)絡(luò)Qnet：表征某環(huán)境狀態(tài)下某一動(dòng)作對(duì)應(yīng)的價(jià)值函數(shù)，輸入為[s,d]。

2)策略函數(shù)網(wǎng)絡(luò)Pnet：輸入為s，輸出為d，此網(wǎng)絡(luò)即Ulearn優(yōu)化求解器。

3)片段：PSO從開(kāi)始到終止的全過(guò)程。

基于DDPG的改進(jìn)PSO算法設(shè)計(jì)流程如下：

1)使用權(quán)重θQ和θP隨機(jī)初始化值函數(shù)網(wǎng)絡(luò)Qnet和策略函數(shù)網(wǎng)絡(luò)Pnet。

3)對(duì)于PSO算法全過(guò)程：

(2)隨機(jī)初始化粒子位置和Ulearn，獲取初始狀態(tài)空間的觀測(cè)量。

(3)對(duì)于PSO算法中的每一步迭代：

①根據(jù)當(dāng)前策略網(wǎng)絡(luò)選擇動(dòng)作d=P(s|θP)+。

②執(zhí)行動(dòng)作d，得到獎(jiǎng)勵(lì)r和新的狀態(tài)s。

③將狀態(tài)組(s,d,r,s)保存到經(jīng)驗(yàn)回放池。

④從經(jīng)驗(yàn)回放池中NR次隨機(jī)取樣狀態(tài)組。

⑤計(jì)算取樣值：

(36)

式中：i=1,2,…,NR，kQ為預(yù)先設(shè)定的系數(shù)。

⑥通過(guò)最小化損失函數(shù)來(lái)更新值函數(shù)網(wǎng)絡(luò)，得到新的θQ：

(37)

⑦使用采樣策略梯度更新策略網(wǎng)絡(luò)：

(38)

⑧定義常系數(shù)τθ，更新目標(biāo)網(wǎng)絡(luò)：

(39)

圖3 強(qiáng)化學(xué)習(xí)在PSO中的應(yīng)用Fig.3 The application of reinforcement learning into PSO

經(jīng)過(guò)DDPG大量離線訓(xùn)練后，可得到包含問(wèn)題模型信息的策略網(wǎng)絡(luò)，將其作為PSO算法的Ulearn求解器，即可實(shí)現(xiàn)基于DDPG的改進(jìn)PSO算法，如圖3所示。PSO算法終止的標(biāo)準(zhǔn)是：達(dá)到最大迭代次數(shù)或gbest的值在過(guò)去3代中變化不超過(guò)5%；另外，約束條件借助罰函數(shù)法處理[21]。

4 數(shù)值仿真及結(jié)果分析

首先，設(shè)置上升段過(guò)程約束為：qmax=400 kPa，Nm,max=6，αmin,1=-5°，αmax,1=15°，αmin,2=-8°，αmax,2=8°，αmin,3=-10°，αmax,3=10°，αmin,4=-10°，αmax,4=10°，|dα/dt|max=8°/s。同時(shí)，設(shè)運(yùn)載器起飛質(zhì)量為60000 kg，其中燃料質(zhì)量為50000 kg；特征面積為15 m2，TBCC渦輪和沖壓發(fā)動(dòng)機(jī)的進(jìn)氣道入口面積分別為7.0 m2和6.0 m2；采用液氧/液氫燃料，熱值為120.8 MJ/kg。另外，設(shè)置PSO算法采用40個(gè)粒子，最大迭代次數(shù)為20。

建立四個(gè)仿真場(chǎng)景從而驗(yàn)證算法的適應(yīng)性和有效性，相應(yīng)的上升段初始運(yùn)動(dòng)狀態(tài)和終端狀態(tài)約束如表2所示。最終，PSO算法在迭代13次后終止，參數(shù)優(yōu)化結(jié)果見(jiàn)表3，相應(yīng)仿真結(jié)果見(jiàn)圖4～9。

表2 初始狀態(tài)和終端約束Table 2 The initial states and the terminal constraints

表3 參數(shù)優(yōu)化結(jié)果Table 3 The results of parameter optimization

圖4 速度隨時(shí)間變化情況Fig.4 Time histories of velocity

圖5 高度隨時(shí)間變化情況Fig.5 Time histories of altitude

圖6 飛行路徑角隨時(shí)間變化情況Fig.6 Time histories of flight path angle

圖7 過(guò)程約束隨時(shí)間變化情況Fig.7 Time histories of path constraints

圖8 攻角及其變化率隨時(shí)間變化情況Fig.8 Time histories of AOA and its change rate

圖9 攻角及其變化率隨時(shí)間變化情況(場(chǎng)景1)Fig.9 Time histories of AOA and its change rate (Case1)

從圖4～6可以看出，在所有情況下優(yōu)化算法均能滿(mǎn)足上升段終端約束。四種場(chǎng)景下的上升段飛行時(shí)間分別為71.8 s、86.1 s、74.0 s和103.3 s。四種場(chǎng)景下的上升段終端飛行路徑角分別為14.7°、13.4°、22.1°和14.3°，均滿(mǎn)足表2中給出的終端約束條件。圖7表明，得到上升段軌跡能滿(mǎn)足動(dòng)壓和法向過(guò)載約束，其中四種場(chǎng)景下的最大動(dòng)壓分別為363 kPa、301 kPa、275 kPa和158 kPa，最大法向過(guò)載(絕對(duì)值)分別為5.0、4.2、4.9和5.4。圖8表明，上升段攻角曲線十分平滑，攻角變化率始終保持在限定區(qū)間內(nèi)。圖9以場(chǎng)景1為例給出了上升段攻角及其變化率的仿真曲線；結(jié)果表明，攻角在各動(dòng)力模態(tài)段中均能始終保持在指定的約束范圍內(nèi)，同時(shí)根據(jù)變化率滿(mǎn)足約束條件。

最后，將基本PSO算法(BPSO)和基于高斯擾動(dòng)的PSO算法(GPSO)[22]與本文提出的基于強(qiáng)化學(xué)習(xí)的PSO算法QPSO對(duì)比，驗(yàn)證QPSO在計(jì)算效率上的提升。設(shè)三種PSO算法采用相同的約束處理方式，粒子數(shù)均為40，最大迭代次數(shù)均為20；另外，考慮到PSO算法本身具有隨機(jī)性，對(duì)三種PSO算法均進(jìn)行50次獨(dú)立仿真。以表2中的場(chǎng)景1為例，得到仿真結(jié)果見(jiàn)表4。仿真結(jié)果表明，QPSO搜索到的最優(yōu)結(jié)果和最差結(jié)果均優(yōu)于BPSO和GPSO。從尋優(yōu)能力來(lái)看，QPSO獲得的指標(biāo)較BPSO而言平均高出11.1%，較GPSO而言平均高出6.5%；從計(jì)算效率來(lái)看，QPSO較BPSO而言平均減少35.3%，較GPSO而言平均減少15.4%。另外，從算法輸出穩(wěn)定性來(lái)看(標(biāo)準(zhǔn)差)，QPSO同樣優(yōu)于BPSO和GPSO。綜上，本文提出的基于收斂性分析和強(qiáng)化學(xué)習(xí)的改進(jìn)策略有效、顯著地提升了PSO算法的綜合搜索性能。

表4 不同PSO算法對(duì)比結(jié)果(50次獨(dú)立仿真)

*平均值向上取整

5 結(jié) 論

面向未來(lái)新一代組合動(dòng)力水平起飛可重復(fù)使用運(yùn)載器，針對(duì)跨大空/速域的完整上升段軌跡優(yōu)化問(wèn)題展開(kāi)深入研究，獲得如下成果：

1)考慮多動(dòng)力模態(tài)切換問(wèn)題，設(shè)計(jì)了一種新的攻角剖面，通過(guò)數(shù)學(xué)推導(dǎo)解析求解了攻角約束的滿(mǎn)足條件，減少了優(yōu)化算法需要處理的約束數(shù)量，降低了優(yōu)化問(wèn)題的求解難度。

2)提出了一種改進(jìn)PSO算法。在收斂性分析的基礎(chǔ)上，分析了影響PSO算法尋優(yōu)過(guò)程的關(guān)鍵參數(shù)以及這些參數(shù)PSO尋優(yōu)能力間的映射關(guān)系。在此基礎(chǔ)上，引入強(qiáng)化學(xué)習(xí)機(jī)制對(duì)PSO尋優(yōu)過(guò)程進(jìn)行智能控制，從本質(zhì)上提升了PSO求解效率。

3)通過(guò)多場(chǎng)景數(shù)學(xué)仿真驗(yàn)證了方法的正確性和有效性，同時(shí)通過(guò)對(duì)比分析驗(yàn)證了改進(jìn)PSO算法在尋優(yōu)能力和效率上的提升。