馮子鑫 薛文超 張 冉 齊洪勝
近些年來(lái),由于可以顯著降低太空探索的成本,關(guān)于可回收火箭的研究受到廣泛的關(guān)注.20 世紀(jì)50~60 年代,美國(guó)在航天飛機(jī)上首次實(shí)現(xiàn)航天運(yùn)載器的部分重復(fù)使用,隨后進(jìn)行大量的可重復(fù)使用運(yùn)載器的研究[1].21 世紀(jì)以來(lái),SpaceX、Blue Origin、Maste Space System 等商業(yè)公司同樣開展大量關(guān)于可回收火箭的研究[2-3].動(dòng)力下降作為可回收火箭中的一個(gè)重要問(wèn)題,要求利用稠密大氣的空氣阻力和火箭發(fā)動(dòng)機(jī)的推力來(lái)實(shí)現(xiàn)火箭的定點(diǎn)軟著陸,是火箭著陸的最后一環(huán).由于環(huán)境存在許多干擾,在動(dòng)力下降段中常常存在大量的不確定性,尤其是未知的風(fēng)場(chǎng),這給動(dòng)力下降段的制導(dǎo)控制帶來(lái)挑戰(zhàn).
目前針對(duì)火箭動(dòng)力下降問(wèn)題的制導(dǎo)方法主要可以分為三類: 解析制導(dǎo)方法,基于軌跡優(yōu)化的制導(dǎo)方法,基于學(xué)習(xí)的制導(dǎo)方法.解析制導(dǎo)方法是指制導(dǎo)指令由火箭狀態(tài)量解析表示的一系列方法,包括多項(xiàng)式制導(dǎo)[4]、重力轉(zhuǎn)彎制導(dǎo)[5]和近似最優(yōu)解析制導(dǎo)[6](零控位置/速度誤差制導(dǎo))等.這些制導(dǎo)方法具有形式簡(jiǎn)單、易于實(shí)現(xiàn)等優(yōu)點(diǎn),已經(jīng)在阿波羅登月、嫦娥登月等實(shí)際場(chǎng)景中得到應(yīng)用.但是為了獲得解析的制導(dǎo)指令,往往使用較為簡(jiǎn)單的問(wèn)題模型,難以處理大氣層內(nèi)的動(dòng)力下降問(wèn)題.基于學(xué)習(xí)的制導(dǎo)方法[7-8]則是近年來(lái)隨著深度學(xué)習(xí)興起而發(fā)展起來(lái)的方法,通過(guò)訓(xùn)練可以顯著提高方法對(duì)初始狀態(tài)和模型不確定性的適應(yīng)性.但是基于學(xué)習(xí)的方法具有難以設(shè)計(jì)獎(jiǎng)勵(lì)函數(shù)、可解釋性差等缺點(diǎn),因此目前還難以在實(shí)際場(chǎng)景中進(jìn)行應(yīng)用.
基于軌跡優(yōu)化的制導(dǎo)方法[9-10]將控制和狀態(tài)作為待優(yōu)化變量,在一定的性能指標(biāo)下通過(guò)求解優(yōu)化問(wèn)題來(lái)獲得制導(dǎo)指令,其主要分為間接法和直接法.間接法是將動(dòng)力下降問(wèn)題建模為最優(yōu)控制問(wèn)題,通過(guò)極大值原理將最優(yōu)控制問(wèn)題轉(zhuǎn)化為兩點(diǎn)邊值問(wèn)題進(jìn)行求解[11-12].但是對(duì)于較為復(fù)雜的問(wèn)題模型來(lái)說(shuō),求解轉(zhuǎn)化后的兩點(diǎn)邊值問(wèn)題同樣較為困難,因此間接法更多只能用在一些較為特殊的、化簡(jiǎn)后的模型.直接法則是指通過(guò)對(duì)原始最優(yōu)控制問(wèn)題進(jìn)行直接描述和尋優(yōu)的方法.凸優(yōu)化是一種在動(dòng)力下降問(wèn)題中得到廣泛應(yīng)用的直接法,基于凸優(yōu)化的軌跡規(guī)劃將原始的非凸問(wèn)題轉(zhuǎn)化或近似轉(zhuǎn)化為凸問(wèn)題,再通過(guò)離散化將無(wú)窮維的連續(xù)系統(tǒng)轉(zhuǎn)化為有限維的離散系統(tǒng),最后通過(guò)成熟的凸優(yōu)化求解算法進(jìn)行快速求解.文獻(xiàn)[13-14]首次在動(dòng)力下降問(wèn)題中提出無(wú)損凸化的概念,通過(guò)引入松弛變量將控制幅值約束轉(zhuǎn)化為凸約束并推廣到更一般的情形,從理論上嚴(yán)格證明轉(zhuǎn)化前后最優(yōu)解的等價(jià)性.文獻(xiàn)[15]則將上述結(jié)果推廣到同時(shí)具有控制和狀態(tài)約束的情況,同樣證明轉(zhuǎn)化前后最優(yōu)解的等價(jià)性.但是上述無(wú)損凸化的方法只對(duì)特定的約束形式成立,難以推廣到動(dòng)力著陸段的一般情形.文獻(xiàn)[16]提出序列線性化的方法,針對(duì)非凸的、一般的非線性不等式,每次迭代使用上一次迭代的軌跡并在其附近進(jìn)行一階近似展開,求解近似的凸優(yōu)化問(wèn)題,使用信賴域技術(shù),通過(guò)多步迭代使近似最優(yōu)解收斂.隨后同樣有一系列工作使用凸優(yōu)化作為軌跡優(yōu)化的工具[17-19].目前,基于凸優(yōu)化的軌跡規(guī)劃方法已經(jīng)被廣泛應(yīng)用于可回收火箭的動(dòng)力著陸問(wèn)題中.但是上述方法針對(duì)的均是確定性問(wèn)題,未考慮不確定性的影響.
對(duì)于軌跡規(guī)劃中不確定性的研究則相對(duì)較少.一種考慮不確定性的方法是基于敏感度矩陣(系統(tǒng)狀態(tài)關(guān)于未知參數(shù)的雅可比矩陣)的方法[20-23].文獻(xiàn)[20]將敏感度矩陣作為增廣狀態(tài)并在原始的最優(yōu)控制問(wèn)題中增加關(guān)于增廣狀態(tài)的初始和終端約束,通過(guò)求解增廣的最優(yōu)控制問(wèn)題來(lái)進(jìn)行魯棒軌跡規(guī)劃.文獻(xiàn)[21]則設(shè)計(jì)一個(gè)基于敏感度矩陣和風(fēng)場(chǎng)觀測(cè)器的臨近最優(yōu)制導(dǎo)律,并且為處理真實(shí)風(fēng)場(chǎng)的不確定性,加入擾動(dòng)觀測(cè)補(bǔ)償器以提升魯棒性.文獻(xiàn)[22]考慮存在參數(shù)不確定性的軌跡跟蹤問(wèn)題,通過(guò)最小化狀態(tài)敏感度和輸入敏感度進(jìn)行魯棒軌跡規(guī)劃.文獻(xiàn)[23]將敏感度矩陣作為增廣狀態(tài),對(duì)不等式約束進(jìn)行一階近似并使用敏感度矩陣計(jì)算一階項(xiàng),將原問(wèn)題轉(zhuǎn)化為一個(gè)易求解的優(yōu)化問(wèn)題.此外,基于多項(xiàng)式混沌展開(Polynomial chaos expansion,PCE)的方法[24]則是軌跡規(guī)劃中處理不確定性的一種典型方法,其主要思想為,將存在不確定性的隨機(jī)過(guò)程展開為一系列正交的隨機(jī)多項(xiàng)式基底函數(shù)的加權(quán)和,從而建立不確定輸入與不確定輸出的映射關(guān)系.文獻(xiàn)[25]則將PCE 與凸優(yōu)化的方法相結(jié)合,在僅考慮初值不確定性的情況下,通過(guò)序列線性化的方法降低求解問(wèn)題的計(jì)算代價(jià).文獻(xiàn)[26]通過(guò)引入譜分解技術(shù)來(lái)自適應(yīng)更新基函數(shù),以處理長(zhǎng)時(shí)間的優(yōu)化場(chǎng)景可能導(dǎo)致PCE 方法發(fā)散的問(wèn)題.文獻(xiàn)[27]通過(guò)特征選擇和序列采樣提出序列增強(qiáng)PCE 方法,相比于原始的PCE 方法具有更高的計(jì)算效率和精度.
可以看出,現(xiàn)有對(duì)存在不確定性的動(dòng)力下降問(wèn)題或者軌跡規(guī)劃的研究要么對(duì)不確定性(或者風(fēng)場(chǎng))的建模較為簡(jiǎn)單[23-24],要么使用隨機(jī)的PCE 方法,通過(guò)期望與方差等統(tǒng)計(jì)特性來(lái)衡量性能指標(biāo).這種方法需要把原始的隨機(jī)問(wèn)題轉(zhuǎn)化為更高維的確定性問(wèn)題,大大增加求解的復(fù)雜度.此外,PCE 方法的最優(yōu)指標(biāo)是在期望和方差等統(tǒng)計(jì)學(xué)意義下的平均最優(yōu),難以應(yīng)對(duì)極端場(chǎng)景.本文則使用一種確定性的魯棒優(yōu)化(Robust optimization,RO)方法[23],該方法轉(zhuǎn)化后的問(wèn)題維度遠(yuǎn)低于PCE 方法,并且考慮在參數(shù)最差時(shí)的最優(yōu)性,可以處理更加極端的情況.本文的主要工作包括: 首先,建立更加貼近實(shí)際場(chǎng)景的風(fēng)場(chǎng)模型并在確定性框架下給出動(dòng)力下降的RO 問(wèn)題;其次,使用一種對(duì)不等式約束采取一階近似的RO 方法,得到一個(gè)可以易于求解的單階段RO 算法;然后,通過(guò)理論分析單階段RO 算法在終端約束范圍較小時(shí)可能無(wú)解的缺陷,進(jìn)一步提出一種改進(jìn)的多階段RO 算法;最后,通過(guò)仿真說(shuō)明多階段RO 算法對(duì)未知不確定風(fēng)場(chǎng)同時(shí)具有較高的最終狀態(tài)精度和較強(qiáng)的魯棒性.
大氣層內(nèi)的可回收火箭的動(dòng)力下降問(wèn)題是指火箭在發(fā)動(dòng)機(jī)推力、重力和空氣阻力的影響下從高空狀態(tài)被引導(dǎo)到指定的著陸點(diǎn).在本文中,x軸、y軸和z軸分別指向東、北和天空.為了簡(jiǎn)便,文中的下標(biāo) (·)x、(·)y和 (·)z分別代表對(duì)應(yīng)量在x軸、y軸和z軸上的分量.考慮如下大氣層內(nèi)火箭的動(dòng)力學(xué)模型
其中,r=[rx,ry,rz]T∈R3,v=[vx,vy,vz]T∈R3分別是位置向量和速度向量,g代表重力加速度,m∈R是火箭質(zhì)量,κ∈R代表火箭恒定的燃料消耗系數(shù),u=[ux,uy,uz]T∈R3是系統(tǒng)的控制輸入,代表火箭引擎提供的推力,其滿足如下幅值約束
其中,D(rz,v,vw)代表空氣阻力,vw(rz)∈R3代表風(fēng)場(chǎng)向量,ρ(rz)∈R是大氣密度,Sref∈R,CD∈R分別代表火箭的參考面積和空氣阻力系數(shù).
設(shè)初始時(shí)刻為t0=0,該問(wèn)題的初始條件為
終端時(shí)刻tf自由,在終端時(shí)刻需要滿足終端條件
rf,vf代表期望的終端位置和速度.性能指標(biāo)為終端時(shí)刻和推力二次型積分之間的加權(quán)和
其中,γ代表加權(quán)權(quán)重.
記X=[rT,vT,m]T,X0=[,,m0]T,那么動(dòng)力學(xué)模型和初值可以寫為如下緊湊形式
標(biāo)稱優(yōu)化(Nominal optimization,NO)問(wèn)題1.給定初始狀態(tài)X0以及目標(biāo)終端狀態(tài)Xf,選擇合適的控制律u,在vw=(rz)的情況下,滿足式(3)、式(5)和式(6)的同時(shí)最小化性能指標(biāo)(4).
NO 問(wèn)題使用標(biāo)稱風(fēng)場(chǎng)信息進(jìn)行計(jì)算,因此可以直接使用標(biāo)準(zhǔn)優(yōu)化算法進(jìn)行求解.基于NO 問(wèn)題可以得出如下NO 算法.
NO 算法.輸入X0,rf,vf,求解NO 問(wèn)題,輸出最優(yōu)控制和相應(yīng)的狀態(tài)軌跡與終端時(shí)刻.
第1.1 節(jié)中的NO 問(wèn)題使用標(biāo)稱風(fēng)場(chǎng),但在實(shí)際場(chǎng)景中,風(fēng)場(chǎng)具有很強(qiáng)的不確定性.為了在模型中考慮不確定性的影響,在風(fēng)場(chǎng)模型中引入不確定參數(shù).
假設(shè)風(fēng)場(chǎng)在y軸方向存在不確定性,不確定性風(fēng)場(chǎng)為是未知參數(shù),d代表未知參數(shù)的維度.具體而言
未知參數(shù)s的取值范圍滿足
此外,由于實(shí)際場(chǎng)景中z軸方向的風(fēng)速可以忽略不計(jì),因此假設(shè)z軸方向不存在風(fēng).對(duì)于x軸方向的風(fēng)場(chǎng),為了簡(jiǎn)單起見,本文假設(shè)x軸方向的風(fēng)場(chǎng)是確定已知的,滿足
其中,k,b均為已知參數(shù).最后,記vw(rz,s)代表上述定義的、包含不確定參數(shù)s的風(fēng)場(chǎng).
在第1.2 節(jié)給出的不確定風(fēng)場(chǎng)vw(rz,s)下,考慮如下魯棒優(yōu)化問(wèn)題:
魯棒優(yōu)化(RO)問(wèn)題2.
其中,Lri,Lvi代表事先確定的終端約束范圍.
注1.RO 問(wèn)題2 在滿足對(duì)任意參數(shù)都不違背約束條件的可行解中尋求最小化優(yōu)化指標(biāo)的解,這使得RO 問(wèn)題2 的最優(yōu)解即使在最差的參數(shù)條件下也能擁有較好的效果.
注2.注意到RO 問(wèn)題2 的終端約束和NO 問(wèn)題1 的終端約束不同,這是因?yàn)樵陲L(fēng)場(chǎng)存在不確定性的情況下,使火箭嚴(yán)格滿足給定的終端條件式(3)是難以實(shí)現(xiàn)的,而要求火箭盡可能落在給定終端狀態(tài)的一定范圍之內(nèi)更加具有可行性.
注3.在式(9)中選擇不同的約束范圍Lri,Lvi,是因?yàn)榈?.2 節(jié)定義的風(fēng)場(chǎng)不確定性出現(xiàn)在y軸方向,所以相比于y軸方向而言,火箭由風(fēng)場(chǎng)不確定性導(dǎo)致的、在x軸和z軸方向的偏差較小,自然也應(yīng)該取不同的約束范圍.
一般而言,直接獲得RO 問(wèn)題2 的最優(yōu)解是十分困難的,特別是當(dāng)系統(tǒng)的動(dòng)態(tài)方程較為復(fù)雜時(shí).接下來(lái)將簡(jiǎn)要介紹一種對(duì)不等式約束采用一階近似的RO 方法,該方法可以將RO 問(wèn)題2 近似轉(zhuǎn)化為一個(gè)易求解的最優(yōu)控制問(wèn)題,并且轉(zhuǎn)化后的問(wèn)題維度低于PCE 方法.該方法更詳細(xì)的推導(dǎo)過(guò)程見文獻(xiàn)[23].考慮如下一般形式的非線性優(yōu)化問(wèn)題
在上述優(yōu)化問(wèn)題中,通過(guò)等式約束確定一個(gè)X與u,s之間的隱式關(guān)系式,不妨設(shè)為
將F(X,u,s)=0兩端同時(shí)關(guān)于s求微分,可得Xs(u,s)滿足
其中,變量下標(biāo)代表關(guān)于該變量求偏導(dǎo)數(shù).
另一方面,在隱式關(guān)系X=X(u,s)下,令
其中,Sτ由式(8)定義.將(u,s)一階泰勒展開可得
其中,?s代表關(guān)于s的梯度算子,〈·,·〉 代表內(nèi)積.因此
其中,等號(hào)成立是因?yàn)楹諣柕虏坏仁?Holder inequality).為了求取?s(u,),將式(12)關(guān)于s求偏導(dǎo)可得
其中,ei∈Rm是m階單位矩陣的第i列,m是G的維數(shù).將式(14)代入式(13)中,再結(jié)合式(11)和原始問(wèn)題式(10),可以得到原始問(wèn)題在一階近似意義下的魯棒優(yōu)化問(wèn)題
其中,gi代表G的第i個(gè)元素,FX,Fs,GX和Gs均使用參數(shù)標(biāo)稱值計(jì)算,Xs是新的待優(yōu)化變量,代表通過(guò)等式約束F(X,u,s)=0和式(11)所導(dǎo)出的X關(guān)于s的偏導(dǎo)數(shù).
注4.式(15)中不等式約束的第二項(xiàng)可以看作是式(10)的不等式約束對(duì)應(yīng)于參數(shù)變化的“安全裕量”,它代表gi對(duì)于參數(shù)改變的敏感度.雖然在計(jì)算中使用的是參數(shù)標(biāo)稱值,但是通過(guò)在式(15)中添加安全裕量來(lái)盡可能保證實(shí)際的參數(shù)不會(huì)導(dǎo)致優(yōu)化結(jié)果違背不等式約束.
文獻(xiàn)[23]指出,雖然上述魯棒優(yōu)化方法不能保證不等式約束G(X(u,s),u,s)≤0對(duì)?s∈Sτ都一定成立,但是在一定的連續(xù)性假設(shè)下,它能保證超出約束范圍的值具有一個(gè)關(guān)于τ2的上界.
為了敘述簡(jiǎn)便,令
?v代表風(fēng)和火箭的相對(duì)速度,分別代表 ?v在x,y,z軸上的分量.由于風(fēng)場(chǎng)在z軸方向風(fēng)速為零,因此 ?vz=vz.另外,令
那么式(2)中的空氣阻力可以表示為
此外,對(duì)于一個(gè)矩陣P,引入符號(hào)P(i:j,k:l)代表由P的第i行到第j行、第k列到第l列組成的子矩陣,P(i:j,:)代表由P的第i行到第j行構(gòu)成的子矩陣.將上述魯棒優(yōu)化方法應(yīng)用到RO 問(wèn)題2 中,可以得到:
RO 問(wèn)題3.
RO 問(wèn)題3 使用標(biāo)稱風(fēng)場(chǎng)信息進(jìn)行計(jì)算,因此可以直接使用傳統(tǒng)優(yōu)化算法求解.基于RO 問(wèn)題3可以得出如下單階段RO 算法.
單階段RO 算法.輸入X0,rf,vf,Lri,Lvi,i=x,y,z,求解RO 問(wèn)題3,輸出最優(yōu)控制、相應(yīng)的狀態(tài)軌跡與終端時(shí)刻.
第1.3 節(jié)得到了參數(shù)的變化導(dǎo)致的終端約束函數(shù)變化,即注4 所說(shuō)的安全裕量,也即式(17)中的‖rsi(tf)A‖1,‖vsi(tf)A‖1.由式(13)的推導(dǎo)可知,這種方法在每個(gè)時(shí)刻總會(huì)選取最差的參數(shù)(盡管實(shí)際場(chǎng)景中很難出現(xiàn)這種情況),并且隨著時(shí)間推移,這種最差的參數(shù)對(duì)約束函數(shù)帶來(lái)的變化也會(huì)逐漸累積,最終會(huì)導(dǎo)致安全裕量比較大.這就導(dǎo)致RO 問(wèn)題3 出現(xiàn)矛盾: 一方面,終端約束范圍Lri和Lvi,i=x,y,z越小越好,因?yàn)檫@代表著能使火箭更加精準(zhǔn)地以給定速度落到給定地點(diǎn);另一方面,如果Lri和Lvi比安全裕量還小,那么式(17)中的不等式約束將無(wú)法滿足,導(dǎo)致優(yōu)化問(wèn)題無(wú)解,因此獲取更多安全裕量的信息十分重要.事實(shí)上,本節(jié)將會(huì)定量地給出安全裕量的一個(gè)上界.為此,首先給出一個(gè)假設(shè)和一個(gè)引理.
在假設(shè)1 下,可以證明有如下引理成立.
引理1.若假設(shè)1 成立,在RO 問(wèn)題3 中,Xs等式約束內(nèi)的fX(4:6,:)和fs(4:6,:)滿足
由式(26)和假設(shè)1,有
其中,第二個(gè)不等式是因?yàn)槭?7).因此
對(duì)于fX(4:6,4:6),可以得到
對(duì)于fs(4:6,:),在式(23)中,如果記
那么可以得到
引理1 說(shuō)明式(17)中的fX,fs都是有界的,因此可以通過(guò)=fXXs+fs求出rsi,vsi,i=x,y,z的上界,進(jìn)而求出式(17)中安全裕量的界.具體而言,給出如下定理.
定理1.若Xs滿足
那么對(duì)?tf ≥0,式(17)中 的‖rsi(tf)A‖1和‖vsi(tf)A‖1有界,滿足
其中,i=x,y,z,C=1+(tf+1)M1+2M2.
證明.定理中的結(jié)果需要證明rs和vs所有行向量的1 范數(shù)有界,因此只需證明 ‖rs‖∞和‖vs‖∞有界.為了證明 ‖rs‖∞和‖vs‖∞有界,首先證明它們的Frobenius 范數(shù)(F 范數(shù))有界,再通過(guò)矩陣范數(shù)之間的關(guān)系(對(duì)任意n階方陣A,有‖A‖∞≤得到‖rs‖∞和‖vs‖∞的界.為了簡(jiǎn)便,在本證明中,令 |·|代表向量或者矩陣的2范數(shù),對(duì)于標(biāo)量則代表取絕對(duì)值,并記rsj=rs(:,j),vsj=vs(:,j),j=1,2,3 代表rs和vs的各個(gè)列向量,注意和rsi,vsi區(qū)分,i=x,y,z.它們代表rs和vs的各個(gè)行向量.將式(19)拆為三個(gè)列向量對(duì)應(yīng)的微分方程可得
兩端同時(shí)取2 范數(shù)并由引理1 可得
或者等價(jià)地
〈·,·〉 代表內(nèi)積.對(duì)上式兩端在 [0,t] 上積分,注意到sj(0)=0,可得
對(duì)于上式中的 |rsj|2,有
將上式代入式(30)中,得到
由Gronwall 不等式可得
回顧vsj代表vs的列向量,由F 范數(shù)的定義可知
又因?yàn)榫仃嚨? 范數(shù)總是小于矩陣的F 范數(shù),所以
由式(33)即可得知定理1 的式(29)成立.對(duì)于rs,由式(31)可得
類似于式(32)和式(33)的推導(dǎo),有
由式(34)即證定理1 的式(28)成立.
從定理1 可以看出,式(28)~(29)中給出的界是隨時(shí)間指數(shù)增長(zhǎng)的,對(duì)于一些需要較長(zhǎng)時(shí)間才能完成的優(yōu)化場(chǎng)景,顯然會(huì)導(dǎo)致必須選取較大的終端約束范圍優(yōu)化問(wèn)題才能有解,但是大的終端約束范圍又會(huì)降低落點(diǎn)的精度.對(duì)于優(yōu)化時(shí)間較長(zhǎng)的場(chǎng)景,為了通過(guò)縮短優(yōu)化時(shí)長(zhǎng)使問(wèn)題具有可行性,本文提出一種多階段RO 算法.其主要思想是通過(guò)將整個(gè)優(yōu)化進(jìn)程拆分為多個(gè)串聯(lián)的子階段,在每個(gè)子階段中分別求解RO 子問(wèn)題,縮短優(yōu)化時(shí)長(zhǎng),進(jìn)而避免優(yōu)化時(shí)間較長(zhǎng)致使RO 問(wèn)題3 無(wú)解.
具體而言,在如圖1 所示的K階段RO 算法流程圖中,選取每個(gè)階段的預(yù)估終端時(shí)刻t1,t2,···,tK,使得每個(gè)階段的預(yù)估時(shí)長(zhǎng)為標(biāo)稱軌跡時(shí)長(zhǎng)的1/K.將標(biāo)稱軌跡在t1時(shí)刻的狀態(tài)rNO(t1)、vNO(t1)作為第一階段的目標(biāo)終值,將初始狀態(tài)X0作為第一階段的初值,使用單階段RO 算法進(jìn)行優(yōu)化,將最優(yōu)解輸入到火箭中,得到火箭在第一階段末的位置在第二階段中,則以為第二階段的初值,以rNO(t2),vNO(t2)作為第二階段的目標(biāo)終值,同樣使用單階段RO 算法求解,如此進(jìn)行下去,直到第K階段結(jié)束.需要注意的是,最后一個(gè)階段的目標(biāo)終值滿足

圖1 多階段RO 算法流程圖Fig.1 Flow chart of multi-stage RO algorithm
算法1.多階段RO 算法
注5.關(guān)于階段數(shù)的選擇,在實(shí)際場(chǎng)景中應(yīng)綜合考慮多方面因素.一方面,階段數(shù)量越多,每個(gè)階段的優(yōu)化時(shí)長(zhǎng)就越短,根據(jù)本節(jié)的分析,就可以求解更小終端約束范圍的優(yōu)化問(wèn)題,進(jìn)而使最終狀態(tài)更加精確.另一方面,階段數(shù)量進(jìn)一步增加也會(huì)帶來(lái)一些問(wèn)題.首先,階段數(shù)增多會(huì)使每個(gè)階段的優(yōu)化時(shí)長(zhǎng)縮短,縮短優(yōu)化時(shí)長(zhǎng)是為了縮小約束范圍Lri,Lvi,i=x,y,z.但是從定理1 可以看出,當(dāng)優(yōu)化時(shí)長(zhǎng)tf減小并趨于零,式(17)中的安全裕量同樣會(huì)減小并趨于零,若Lri,Lvi同樣減小并趨于零,此時(shí)RO 問(wèn)題3 將近似等價(jià)于NO 問(wèn)題,這表明當(dāng)優(yōu)化時(shí)間過(guò)短時(shí),RO 算法相比于NO 算法便不再具有優(yōu)勢(shì).其次,當(dāng)階段數(shù)增加時(shí),整個(gè)優(yōu)化進(jìn)程需要求解的優(yōu)化問(wèn)題也會(huì)增加,這對(duì)計(jì)算性能有著更大的需求.因此,在實(shí)際場(chǎng)景中,設(shè)計(jì)算法時(shí)應(yīng)綜合考慮火箭初始狀態(tài)、能夠容許的終端約束范圍大小、計(jì)算能力等條件來(lái)選擇階段數(shù).本節(jié)算法將每個(gè)階段的預(yù)估優(yōu)化時(shí)長(zhǎng)取為相等,正是考慮到每個(gè)階段的優(yōu)化時(shí)長(zhǎng)不能過(guò)長(zhǎng)、也不能過(guò)短的因素.
首先說(shuō)明仿真所使用的風(fēng)場(chǎng)模型.本文通過(guò)實(shí)際測(cè)量出的6 個(gè)時(shí)刻風(fēng)場(chǎng)數(shù)據(jù)擬合出6 個(gè)風(fēng)速與高度之間的函數(shù)關(guān)系,將這些擬合函數(shù)各個(gè)參數(shù)對(duì)應(yīng)的最大、最小值的平均值作為參數(shù)標(biāo)稱值.x軸方向的風(fēng)場(chǎng)參數(shù)選為擬合出的標(biāo)稱參數(shù).對(duì)于y軸方向的不確定風(fēng)場(chǎng),將各個(gè)擬合參數(shù)對(duì)應(yīng)的最大、最小值分別作為該參數(shù)取值范圍的上、下界,擬合出的y軸風(fēng)場(chǎng)與高度之間的函數(shù)關(guān)系如圖2 所示.如果令風(fēng)場(chǎng)序號(hào)0 代表標(biāo)稱風(fēng)場(chǎng),風(fēng)場(chǎng)序號(hào)1 到6 代表實(shí)際風(fēng)場(chǎng),那么各個(gè)實(shí)際風(fēng)場(chǎng)的參數(shù)與標(biāo)稱參數(shù)在無(wú)窮范數(shù)意義下的距離如圖3 所示.可以看出在6 個(gè)實(shí)際風(fēng)場(chǎng)中,風(fēng)場(chǎng)1、2、6 與標(biāo)稱參數(shù)差距較大,風(fēng)場(chǎng)4、5 與標(biāo)稱參數(shù)差距較小.因此可以認(rèn)為,相比于標(biāo)稱風(fēng)場(chǎng)而言,風(fēng)場(chǎng)1、2、6 是較差的風(fēng)場(chǎng),風(fēng)場(chǎng)4、5 是較好的風(fēng)場(chǎng).

圖2 y 軸的六個(gè)實(shí)際風(fēng)場(chǎng)和一個(gè)標(biāo)稱風(fēng)場(chǎng)Fig.2 Six actual wind fields and one nominal wind field on the y-axis

圖3 實(shí)際風(fēng)場(chǎng)參數(shù)與標(biāo)稱風(fēng)場(chǎng)參數(shù)的距離Fig.3 Distances between parameters of actual wind fields and parameters of nominal wind field
本節(jié)的仿真參數(shù)如表1 所示,圖4 則給出在表1的仿真參數(shù)條件下隨機(jī)生成的一組風(fēng)場(chǎng),可以看出本節(jié)的參數(shù)條件較好地覆蓋了6 個(gè)實(shí)際風(fēng)場(chǎng).

表1 仿真參數(shù)值Table 1 Simulation parameter values

圖4 在表1 的仿真參數(shù)條件下隨機(jī)生成的一組風(fēng)場(chǎng)Fig.4 A set of wind fields randomly generated under the simulation parameter conditions of Table 1
此外,本節(jié)使用Matlab 軟件包GPOPS 來(lái)求解最優(yōu)控制問(wèn)題.GPOPS 是一個(gè)基于偽譜法和序列二次規(guī)劃方法來(lái)求解最優(yōu)控制問(wèn)題的算法,具有求解速度快、精度高等特點(diǎn).
注6.在表1 中,單階段RO 算法的終端約束范圍大于多階段RO 算法,這是因?yàn)楫?dāng)單階段RO 算法和多階段RO 算法取相同的終端約束范圍參數(shù)時(shí),單階段RO 算法將由于約束范圍太小而無(wú)法求解.
在本節(jié),將通過(guò)對(duì)比單階段RO 算法和多階段RO 算法來(lái)顯示多階段RO 算法在最終落點(diǎn)狀態(tài)精度上的優(yōu)越性.單階段RO 算法和多階段RO 算法在不同風(fēng)場(chǎng)下最終位置和期望位置的誤差如圖5 所示,最終速度和期望速度的誤差如圖6 所示,它們?cè)诟鱾€(gè)實(shí)際風(fēng)場(chǎng)下的軌跡如圖7 所示,控制–時(shí)間關(guān)系則如圖8、圖9 所示.從圖5 和圖6 可以看出,多階段RO 算法無(wú)論是位置還是速度,最終狀態(tài)的精度均高于單階段RO 算法.并且隨著風(fēng)場(chǎng)的變化,多階段RO 算法的最終位置誤差波動(dòng)也明顯小于單價(jià)段RO 算法.這顯示多階段RO 算法相比于單階段RO 算法在最終落點(diǎn)狀態(tài)精度上的優(yōu)越性.

圖5 單階段RO 算法和多階段RO 算法在不同風(fēng)場(chǎng)下的最終位置誤差Fig.5 The terminal position errors of single-stage RO algorithm and multi-stage RO algorithm under different wind fields

圖6 單階段RO 算法和多階段RO 算法在不同風(fēng)場(chǎng)下的最終速度誤差Fig.6 The terminal velocity errors of single-stage RO algorithm and multi-stage RO algorithm under different wind fields

圖7 單階段RO 算法和多階段RO 算法在各個(gè)實(shí)際風(fēng)場(chǎng)下的軌跡Fig.7 The trajectories of single-stage RO algorithm and multi-stage RO algorithm under each actual wind field

圖8 單階段RO 算法的控制–時(shí)間圖Fig.8 Control-time diagram of single-stage RO algorithm

圖9 多階段RO 算法的控制–時(shí)間圖Fig.9 Control-time diagram of multi-stage RO algorithm
此外,單階段RO 算法和多階段RO 算法的運(yùn)行時(shí)間如表2 所示.從表2 可以看出,由于多階段RO 算法需要求解更多優(yōu)化問(wèn)題,所以具有更高的計(jì)算代價(jià),這也表明在設(shè)計(jì)算法時(shí)階段數(shù)不宜過(guò)多.

表2 單階段RO 算法和多階段RO 算法的運(yùn)行時(shí)間Table 2 Running time of single-stage RO algorithm and multi-stage RO algorithm
在本節(jié),將通過(guò)單階段RO 算法與多階段RO算法在不同風(fēng)場(chǎng)下的終端約束滿足情況以及單階段RO 算法與NO 算法的對(duì)比,來(lái)顯示RO 算法的有效性.單階段RO 算法在各個(gè)風(fēng)場(chǎng)下最終位置、最終速度和相應(yīng)的約束范圍如圖10、圖11 所示,其中,虛線代表該方向?qū)?yīng)的約束范圍,圓圈代表對(duì)應(yīng)的位置或速度.多階段RO 算法在各個(gè)風(fēng)場(chǎng)下最終位置、最終速度和相應(yīng)的約束范圍如圖12、圖13所示,其中,虛線代表該方向?qū)?yīng)的約束范圍,叉號(hào)代表對(duì)應(yīng)的位置或速度.從圖10~圖13 可以看出,無(wú)論是單階段RO 算法還是多階段RO 算法,它們?cè)诓煌娘L(fēng)場(chǎng)下均未違背對(duì)應(yīng)的終端狀態(tài)約束.

圖10 單階段RO 算法在不同風(fēng)場(chǎng)下的最終位置與約束范圍Fig.10 The terminal position and its constraint range of single-stage RO algorithm under different wind fields

圖11 單階段RO 算法在不同風(fēng)場(chǎng)下的最終速度與約束范圍Fig.11 The terminal velocity and its constraint range of single-stage RO algorithm under different wind fields

圖12 多階段RO 算法在不同風(fēng)場(chǎng)下的最終位置與約束范圍Fig.12 The terminal position and its constraint range of multi-stage RO algorithm under different wind fields

圖13 多階段RO 算法在不同風(fēng)場(chǎng)下的最終速度與約束范圍Fig.13 The terminal velocity and its constraint range of multi-stage RO algorithm under different wind fields
圖14 則對(duì)比在不同風(fēng)場(chǎng)下單階段RO 算法與NO 算法的最終位置誤差.從圖14 可以看出,除了標(biāo)稱風(fēng)場(chǎng)(風(fēng)場(chǎng)0),單階段RO 算法在其他風(fēng)場(chǎng)下的最終位置精度都高于NO 算法,并且在較差的風(fēng)場(chǎng)條件(風(fēng)場(chǎng)1、2、6)下,兩算法精度差距較大.這表明NO 算法在面對(duì)未知的風(fēng)場(chǎng)時(shí),受風(fēng)場(chǎng)變化的影響較大,而RO 算法則對(duì)變化的風(fēng)場(chǎng)具有更好的魯棒性.上述討論說(shuō)明了RO 算法的有效性.

圖14 單階段RO 算法與NO 算法在不同風(fēng)場(chǎng)下的最終位置誤差Fig.14 The terminal position errors of single-stage RO algorithm and NO algorithm under different wind fields
本文考慮處于不確定風(fēng)場(chǎng)中可回收火箭動(dòng)力下降的魯棒制導(dǎo)問(wèn)題,主要工作包括:
1)建立一個(gè)包含不確定性的風(fēng)場(chǎng)模型用來(lái)刻畫未知風(fēng)場(chǎng),在該模型的基礎(chǔ)上給出火箭動(dòng)力下降的魯棒優(yōu)化問(wèn)題.進(jìn)一步,使用一種對(duì)不等式約束采取一階近似并將一階項(xiàng)作為安全裕量加入約束中的魯棒優(yōu)化方法,得到一個(gè)易于求解的單階段魯棒優(yōu)化算法.
2)從理論上定量地給出安全裕量的上界,借助該上界提出一種多階段魯棒優(yōu)化方法,該方法相比于單階段魯棒優(yōu)化算法,可以求解具有更小的終端約束范圍的優(yōu)化問(wèn)題.
3)通過(guò)仿真對(duì)比分析單階段、多階段魯棒優(yōu)化算法和標(biāo)稱風(fēng)場(chǎng)下的優(yōu)化算法,結(jié)果表明提出的多階段魯棒優(yōu)化方法在落點(diǎn)精度較高的同時(shí),對(duì)不同的風(fēng)場(chǎng)具有較強(qiáng)的魯棒性.