復(fù)雜約束條件下的再入軌跡迭代求解方法

張夢櫻，唐乾剛，韓小軍，張青斌，葛健全

(1. 國防科學(xué)技術(shù)大學(xué) 航天科學(xué)與工程學(xué)院，湖南 長沙410003;2. 第二炮兵裝備部 高新技術(shù)辦公室，北京100085)

0 引言

隨著對(duì)高超聲速滑翔飛行器研究的不斷深入，滿足禁飛區(qū)約束的滑翔軌跡優(yōu)化問題逐漸成為了一個(gè)研究熱點(diǎn)。滿足飛行安全的同時(shí)，需兼顧復(fù)雜航路禁飛區(qū)約束，其軌跡優(yōu)化問題將變得更為復(fù)雜。

在禁飛區(qū)建模方面，文獻(xiàn)［1］建立了基于地理位置的威脅模型，將暴露危險(xiǎn)作為優(yōu)化目標(biāo)設(shè)計(jì)最優(yōu)軌跡。文獻(xiàn)［2］考慮了威脅和目標(biāo)可動(dòng)，建立數(shù)學(xué)模型并采用遺傳算法在線優(yōu)化軌跡。文獻(xiàn)［3］在平面地球假設(shè)下，建立無限長圓柱的禁飛區(qū)模型，設(shè)計(jì)了滿足航路點(diǎn)和禁飛區(qū)約束的最優(yōu)軌跡。文獻(xiàn)［4］的禁飛區(qū)模型建立在球形地球表面，采用球面三角形原理建立無限長圓柱的空間禁飛區(qū)模型。文獻(xiàn)［5］建立了半橢球的禁飛區(qū)模型，但未改變橢球赤道半長軸方向以適用更一般的情況?，F(xiàn)有航路約束或禁飛區(qū)約束建模并不完全適用于實(shí)際情況，為使最優(yōu)軌跡能夠滿足任務(wù)需求，有必要考慮更為復(fù)雜和真實(shí)的禁飛區(qū)模型。

復(fù)雜約束條件下的高超聲速軌跡優(yōu)化問題可歸結(jié)為強(qiáng)非線性、多階段、多約束的最優(yōu)控制問題［6］，對(duì)其求解十分困難。因具有收斂速度快和求解精度高的特點(diǎn)，偽譜法近年來發(fā)展迅速，成為求解復(fù)雜最優(yōu)控制問題的常見工具。文獻(xiàn)［5，7 -9］用Gauss 偽譜法求解包括端點(diǎn)、航路點(diǎn)和禁飛區(qū)約束的通用航空飛行器(CAV)軌跡優(yōu)化問題。國內(nèi)學(xué)者也研究了包括過載、動(dòng)壓、駐點(diǎn)熱流等多約束軌跡優(yōu)化問題［10-14］。在處理復(fù)雜路徑約束問題時(shí)國內(nèi)外研究者多采用分段優(yōu)化的方法［5，8-9，11］，有些會(huì)采用HP自適應(yīng)偽譜方法求解［14-17］此類問題。然而當(dāng)Gauss 節(jié)點(diǎn)較多時(shí)，處理復(fù)雜約束問題會(huì)存在設(shè)計(jì)變量數(shù)目龐大和約束眾多，從而影響求解魯棒性的問題［9，11］。同時(shí)，采用分段方法，也將增加設(shè)定設(shè)計(jì)變量初值的難度。

針對(duì)帶有復(fù)雜禁飛區(qū)約束的高超聲速滑翔飛行器軌跡優(yōu)化問題直接采用偽譜法難以求解的問題，本文提出依序迭代求解策略。

1 數(shù)學(xué)模型

1.1 動(dòng)力學(xué)模型

為簡化計(jì)算，設(shè)地球?yàn)闊o旋圓球，在地心坐標(biāo)系下建立飛行器的三自由度再入動(dòng)力學(xué)模型［18］:

氣動(dòng)力系數(shù)通常是馬赫數(shù)Ma 和攻角α 的函數(shù)，但因高超聲速飛行器的氣動(dòng)遵循馬赫數(shù)無關(guān)原則，可采用(3)式擬合，相應(yīng)氣動(dòng)參數(shù)詳見文獻(xiàn)［15］.

1.2 彈道優(yōu)化模型

1.2.1 目標(biāo)函數(shù)

對(duì)于指定目標(biāo)點(diǎn)的再入軌跡優(yōu)化問題，飛行時(shí)間是影響突防性能一個(gè)重要的優(yōu)化目標(biāo)。

式中:t0和tf分別為起始時(shí)間和終止時(shí)間。

飛行器在大氣中的運(yùn)動(dòng)微分方程采用(1)式。

1.2.2 約束模型

考慮飛行器任務(wù)需求，需要考慮飛行器的指定端點(diǎn)約束，同時(shí)確保熱防護(hù)安全、結(jié)構(gòu)可靠性和姿態(tài)穩(wěn)定性。

1)端點(diǎn)約束:

2)氣動(dòng)熱約束［21］:

3)過載約束

再入段總過載通常較小，主要考慮法向過載ny的約束。

4)動(dòng)壓約束:

1.2.3 控制量選取

攻角α 和側(cè)傾角σ 是影響飛行軌跡的重要變量，選取它們作為優(yōu)化問題的控制變量，主要限制二者的取值范圍，其表達(dá)式為

2 禁飛區(qū)參數(shù)化建模

禁飛區(qū)是指飛行過程中不允許經(jīng)過的空域，如領(lǐng)空、防空識(shí)別區(qū)、雷達(dá)探測區(qū)和反導(dǎo)防御空域等。本文按邊界地面投影類型分為規(guī)則和一般兩類。規(guī)則禁飛區(qū)包括半球形、半橢球、橢圓柱和圓柱形，可建立參數(shù)化數(shù)學(xué)模型描述其覆蓋范圍。對(duì)一般禁飛區(qū)難于數(shù)學(xué)描述，可轉(zhuǎn)化為若干規(guī)則型禁飛區(qū)的疊加組合。因此，典型規(guī)則形狀禁飛區(qū)的數(shù)學(xué)建模是基礎(chǔ)。

如圖1所示，將地球球形表面某點(diǎn)(rnfz，λnfz，φnfz)作為禁飛區(qū)中心，以其為坐標(biāo)原點(diǎn)建立禁飛區(qū)坐標(biāo)系Ox1y1z1，其中Ox1軸指向禁飛區(qū)在過O 點(diǎn)切面上投影的短半軸方向，Oy1指向天，Oz1指向長半軸方向。αnfz為投影短半軸和正北方向的夾角，以順時(shí)針方向?yàn)檎?。OeXeYeZe為地心直角坐標(biāo)系，(r，λ，φ)為飛行器軌跡上某點(diǎn)坐標(biāo)。

將以上定義的位置的球坐標(biāo)轉(zhuǎn)換到地心坐標(biāo)系下:

圖1 禁飛區(qū)模型示意圖Fig.1 Schematic diagram of no-fly zone model

式中:(x，y，z)和(xnfz，ynfz，znfz)分別為飛行器和禁飛區(qū)中心在地心坐標(biāo)系下的坐標(biāo)。

從地心系轉(zhuǎn)化到禁飛區(qū)坐標(biāo)系Ox1y1z1的轉(zhuǎn)化矩陣為

式中:M2、M1、M3表達(dá)式分別為

則飛行器在禁飛區(qū)坐標(biāo)系下的位置坐標(biāo)為

2.1 橢球禁飛區(qū)約束

橢球禁飛區(qū)約束可描述為

式中:a 為橢球赤道長半軸;b 為短半軸;c 為橢球的極半徑。

特別地，當(dāng)a=b=c 時(shí)，橢球禁飛區(qū)退化成圓球禁飛區(qū)。

2.2 無限高橢圓柱禁飛區(qū)約束

當(dāng)c=inf 時(shí)，(18)式退化為

表示一個(gè)無限高橢圓柱禁飛區(qū)約束。

特別地，當(dāng){rnfz=re}∩{αnfz=0}∩{a =b}時(shí)，禁飛區(qū)退化成圓柱禁飛區(qū)，亦可用球面三角表示為

式中:rnfz為禁飛區(qū)在切平面上投影的半徑，rnfz=a =b.

3 軌跡優(yōu)化的迭代求解方法

3.1 偽譜法

近年來，偽譜法因具有收斂速度快和求解精度高的特點(diǎn)，成為求解復(fù)雜最優(yōu)控制問題的常見工具［8］。影響偽譜法收斂速度和計(jì)算精度的主要因素有參數(shù)初值和離散節(jié)點(diǎn)的數(shù)量。由于偽譜法具有以較少的節(jié)點(diǎn)數(shù)獲得較精確的解的優(yōu)點(diǎn)，故參數(shù)初值成為主要因素［9］。偽譜法現(xiàn)廣泛應(yīng)用于軌跡優(yōu)化問題中，但其算法性能仍受初值的影響［9］。

3.2 優(yōu)化流程

考慮多約束的高超聲速飛行器軌跡優(yōu)化問題較為復(fù)雜，直接運(yùn)用偽譜法求解將很難得到優(yōu)化結(jié)果。同時(shí)，第2 節(jié)所建立的禁飛區(qū)模型的復(fù)雜程度也給直接優(yōu)化增添難度。為解決上述問題，彌補(bǔ)算法對(duì)初值的敏感性并提高計(jì)算效率，考慮分段優(yōu)化之外的方法，即采用如圖2所示的基于Gauss 偽譜法的最優(yōu)軌跡依序迭代串行求解流程。

依據(jù)各約束的復(fù)雜程度，逐步施加約束，并將優(yōu)化成功得到的彈道作為下一步的初值彈道進(jìn)行迭代。具體描述如下:

圖2 優(yōu)化流程Fig.2 Flow chart of optimization

步驟1 生成初值彈道。根據(jù)任務(wù)設(shè)定起點(diǎn)和終點(diǎn)的位置速度參數(shù)，擬定初始方向角和速度傾角，在GPOPS 軟件中設(shè)定初值猜測范圍、迭代精度和節(jié)點(diǎn)數(shù)，開始進(jìn)行求解。若滿足約束，則輸出彈道。若不滿足，則返回調(diào)整初值，直到得到優(yōu)化結(jié)果。

步驟2 基于Gauss 偽譜法求滿足端點(diǎn)約束和過程約束的最優(yōu)解。優(yōu)化模型中添加包括動(dòng)壓、駐點(diǎn)熱流以及過載等過程約束，設(shè)置迭代精度和節(jié)點(diǎn)數(shù)，并將上一步所得優(yōu)化結(jié)果作為下一步增加了約束的優(yōu)化問題的初值彈道進(jìn)行迭代求解。若結(jié)果滿足約束，則輸出彈道，若不滿足，首先調(diào)整節(jié)點(diǎn)數(shù)，不能滿足則返回上一步重新調(diào)整計(jì)算。

步驟3 求解滿足端點(diǎn)約束、過程約束和禁飛區(qū)約束的最優(yōu)彈道。如圖3所示在優(yōu)化模型中添加禁飛區(qū)約束。由2.3 節(jié)可得，所提出的禁飛區(qū)模型均可由橢球禁飛區(qū)模型退化變形得到，故首先添加半橢球禁飛區(qū)模型。為提高求解效率，加快收斂，考慮如下逐步擴(kuò)大禁飛區(qū)大小的策略:

圖3 禁飛區(qū)約束處理方法Fig.3 Process of dealing with no-fly zones constraint

設(shè)禁飛區(qū)加權(quán)的長半軸、短半軸及極半徑分別為

將每一步優(yōu)化得到的可行解迭代到下一步，得到滿足禁飛區(qū)約束的軌跡。如在這一步始終不能得到最優(yōu)解，則返回步驟2 調(diào)整初值，重新優(yōu)化彈道。

步驟4 檢查結(jié)果是否滿足所有約束，并判斷收斂性，輸出結(jié)果彈道。

4 仿真計(jì)算

4.1 仿真算例

根據(jù)前述的飛行器運(yùn)動(dòng)模型、優(yōu)化模型和優(yōu)化算法等，設(shè)計(jì)如下3 個(gè)算例:算例1 通過對(duì)經(jīng)典算例的仿真計(jì)算驗(yàn)證了算法的正確性;算例2 對(duì)比不同禁飛區(qū)模型對(duì)優(yōu)化結(jié)果的影響;算例3 對(duì)多組合禁飛區(qū)約束優(yōu)化問題進(jìn)行仿真，驗(yàn)證禁飛區(qū)模型的正確性和所提出計(jì)算方法的有效性。

4.1.1 算例1:算法驗(yàn)證

文獻(xiàn)［19］所提出的最大橫程算例，是典型最優(yōu)控制問題的算例。設(shè)置初始參數(shù)如表1所示，采用所提出優(yōu)化策略進(jìn)行計(jì)算。

表1 仿真參數(shù)設(shè)置Tab.1 The values of simulation parameters

4.1.2 算例2:4 種禁飛區(qū)模型優(yōu)化結(jié)構(gòu)比較

在相同初始運(yùn)動(dòng)參數(shù)和端點(diǎn)約束條件下，以相同地理坐標(biāo)為中心，規(guī)避某一特定區(qū)域，設(shè)置如下4 組禁飛區(qū)，以最小時(shí)間為優(yōu)化目標(biāo)，研究不同禁飛區(qū)模型對(duì)優(yōu)化彈道結(jié)果的影響。為了嚴(yán)格比較4 種禁飛區(qū)模型優(yōu)化結(jié)果，應(yīng)用GPOPS 求解時(shí)，設(shè)置彈道分段數(shù)為10 段，數(shù)值計(jì)算精度為10-3，每段節(jié)點(diǎn)數(shù)最小為4，最大為8.

端點(diǎn)和禁飛區(qū)參數(shù)設(shè)置如表2、表3所示。

表2 端點(diǎn)參數(shù)設(shè)置Tab.2 The values of start point and end point

表3 不同禁飛區(qū)約束分組設(shè)置Tab.3 The values of different group of no-fly zone constraints

4.1.3 算例3:多組合禁飛區(qū)約束優(yōu)化問題

考慮禁飛區(qū)約束的滑翔彈道優(yōu)化問題，低速飛行器將其考慮成一個(gè)低空突防的航路規(guī)劃問題［21］。文獻(xiàn)［22］采用了遺傳算法、航路極坐標(biāo)編碼方式、航跡適應(yīng)度函數(shù)和遺傳操作算子法進(jìn)行仿真，卻不適合高超聲速條件下的軌跡優(yōu)化問題。為滿足禁飛區(qū)約束，采用Gauss 偽譜法進(jìn)行彈道優(yōu)化設(shè)計(jì)。

設(shè)計(jì)一條固定端點(diǎn)的時(shí)間最短飛行軌跡，在飛行過程中需繞過兩個(gè)禁飛區(qū)。第1 個(gè)禁飛區(qū)模擬行政區(qū)域的領(lǐng)空，第2 個(gè)禁飛區(qū)模擬航母戰(zhàn)斗群的防御半徑。具體仿真數(shù)值設(shè)置如表4所示。

表4 多約束條件設(shè)置Tab.4 Multi-constraints

4.2 結(jié)果

4.2.1 算例1

本例仿真運(yùn)算中，設(shè)置彈道分段數(shù)為10 段，精度為10-7，每段節(jié)點(diǎn)數(shù)最小為4，最大為8. 將所得=∞優(yōu)化結(jié)果作為初值，迭代到添加熱流密度約束的問題中，所得優(yōu)化結(jié)果與文獻(xiàn)［19］比較如表5所示。

表5 優(yōu)化結(jié)果比較Tab.5 Comparison of optimal results

本文優(yōu)化所得結(jié)果中，無約束結(jié)果較好復(fù)現(xiàn)了文獻(xiàn)［19］的結(jié)果。比較添加熱流約束的結(jié)果，仿真得到全程駐點(diǎn)最大熱流為795.14 kW/m2，考察熱流極值約束滿足情況，其誤差為0.020 8%，滿足精度要求。與文獻(xiàn)［19］結(jié)果相比，在滿足相同約束的情況下，優(yōu)化目標(biāo)最大橫程更大，因而結(jié)果更優(yōu)。同時(shí)，若采用直接添加約束進(jìn)行優(yōu)化，其仿真時(shí)間為38.70 s，相比迭代結(jié)果要慢。綜上，可驗(yàn)證本文所提出方法的正確性。

4.2.2 算例2

在用偽譜法難以直接求解原問題的情況下，采用本文策略得到4 組優(yōu)化時(shí)間分別為2 021.3 s、1 944.8 s、1 833.8 s 和1 828.5 s，對(duì)應(yīng)圖4第2 ～5 條軌跡結(jié)果。比較可得，球形禁飛區(qū)模型所得優(yōu)化時(shí)間短于柱形禁飛區(qū)模型所得優(yōu)化時(shí)間，且在4 種類型禁飛區(qū)均能滿足規(guī)避如圖陰影區(qū)域的情況下，橢球禁飛區(qū)對(duì)應(yīng)時(shí)間最短，即結(jié)果最優(yōu)。

圖4 優(yōu)化軌跡局部投影Fig.4 Part ground projection view of optimal trajectories

4.2.3 算例3

類似地，由于本例約束復(fù)雜，偽譜法直接求解無法得到結(jié)果，采用3.2 節(jié)提出策略。運(yùn)算中，設(shè)置彈道分段數(shù)為10 段，精度為10-3，每段節(jié)點(diǎn)數(shù)最小為4，最大為8. 取攻角α 和傾側(cè)角σ 為優(yōu)化變量，得到優(yōu)化軌跡及其地面投影和結(jié)果如圖5～圖7所示。

圖5 三維優(yōu)化軌跡Fig.5 3-D optimal trajectory

圖6 優(yōu)化軌跡地面投影Fig.6 Ground projection of optimal trajectories

圖6為無約束和僅添加過程約束的兩條軌跡穿過禁飛區(qū)，而添加禁飛區(qū)約束的優(yōu)化結(jié)果滿足所設(shè)置的路徑約束，即從起點(diǎn)出發(fā)，先后繞過兩個(gè)禁飛區(qū)的側(cè)面實(shí)現(xiàn)規(guī)避。結(jié)果反映了4 條軌跡逐步逼近原問題最優(yōu)解的過程。由圖7(e)可見，控制量的變化滿足了過程約束要求，且在滑翔段，攻角保持在17°附近，即對(duì)應(yīng)最大升阻比攻角。由于所需規(guī)避區(qū)域距離較遠(yuǎn)，軌跡緯度變化平緩(見圖7(c))，其地面投影呈不甚明顯的“S”型，因而傾側(cè)角在飛行過程中并未出現(xiàn)較大翻轉(zhuǎn)，且對(duì)應(yīng)規(guī)避動(dòng)作時(shí)，角度變化率較為穩(wěn)定，在工程上易于實(shí)現(xiàn)。由圖7(h)～圖7(j)可知，飛行器的駐點(diǎn)熱流密度、法向過載及動(dòng)壓均滿足了約束條件。同時(shí)，由圖7可知飛行器終端參數(shù)λf=62.42°，φf=42.28°，hf=20 km，vf=712.62 m/s，γf= -10°，均滿足設(shè)計(jì)要求。

在Windows XP 系統(tǒng)下，Matlab2009 計(jì)算平臺(tái)下，采用第3 節(jié)描述的優(yōu)化策略，從初值計(jì)算到優(yōu)化結(jié)果的得到，仿真計(jì)算總耗時(shí)為4 036.1 s，約1 h，易于仿真實(shí)現(xiàn)。

5 結(jié)論

本文研究了復(fù)雜約束條件下的再入軌跡快速優(yōu)化問題，從復(fù)雜約束的解構(gòu)入手，建立了常用過程約束模型和更符合實(shí)際情況的典型禁飛區(qū)約束模型;針對(duì)采用Gauss 偽譜法所面臨的初值敏感問題，提出依據(jù)復(fù)雜程度逐步添加約束的迭代求解策略;通過與文獻(xiàn)［19］算例進(jìn)行對(duì)比，驗(yàn)證了方法的正確性;分析了不同禁飛區(qū)模型適用的實(shí)際情況，并以仿真算例驗(yàn)證所建立模型的準(zhǔn)確性;以時(shí)間最小為例設(shè)計(jì)了固定端點(diǎn)并考慮復(fù)雜約束條件下的最優(yōu)軌跡。結(jié)果表明，依序逐步添加約束的迭代求解策略，不但能方便解決考慮多種禁飛區(qū)約束的問題，并且避免了Gauss 偽譜法在初值彈道生成上的難題，可作為分段求解軌跡優(yōu)化問題方法的補(bǔ)充和參照。

圖7 優(yōu)化結(jié)果Fig.7 Optimal results

References)

［1］Vian J L，Moore J R. Trajectory optimization with risk minimization for military aircraft［J］. Journal of Guidance，Control，and Dynamics，1989，12(3):311 -317.

［2］Twigg S，Calise A，Johnson E. On-line trajectory optimization including moving threats and targets［C］∥AIAA Guidance，Navigation，and Control Conference and Exhibit. Providence，Rhode Island:American Institute of Aeronautics and Astronautics，2004:2004 -5139.

［3］Jorris T R，Cobb R G. Three-dimensional trajectory optimization satisfying waypoint and no-fly zone constraints［J］. Journal of Guidance，Control，and Dynamics，2009，32(2):551 -572.

［4］陳小慶，侯中喜，劉建霞. 高超聲速滑翔式飛行器再入軌跡多目標(biāo)多約束優(yōu)化［J］. 國防科技大學(xué)學(xué)報(bào)，2009，31(6):77-83.CHEN Xiao-qing，HOU Zhong-xi，LIU Jian-xia. Multi-objective optimization of reentry trajectory for hypersonic glide vehicle with multi-constraints［J］. Journal of National University of Defense Technology，2009，31 (6):77 -83. (in Chinese)

［5］謝愈，劉魯華，湯國建. 多約束條件下高超聲速滑翔飛行器軌跡優(yōu)化［J］. 宇航學(xué)報(bào)，2011，32 (12):2499 -2504.XIE Yu，LIU Lu-hua，TANG Guo-jian. Trajectory optimization for hypersonic glide vehicle with multi-constraints［J］. Journal of Astronautics，2011，32(12):2499 -2504. (in Chinese)

［6］李惠峰. 高超聲速飛行器制導(dǎo)與控制技術(shù)［M］. 北京:中國宇航出版社，2012:242 -243.LI Hui-feng. Guidance and control technology of hypersonic vehicles［M］. Beijing:China Astronautic Publishing House，2012:242 -243.(in Chinese)

［7］Yaple D E. Simulation and application of GPOPS for a trajectory optimization and mission planning tool［D］. Ohio，US:Air Force Institute of Technology，2010.

［8］Jorris T. Common aero vehicle autonomous reentry trajectory optimization satisfying waypoint and no-fly zone constraints［D］.Ohio，US:Graduate School of Engineering and Management，Air Force Institute of Technology，2007.

［9］雍恩米. 高超聲速滑翔式再入飛行器軌跡優(yōu)化與制導(dǎo)方法研究［D］. 長沙:國防科學(xué)技術(shù)大學(xué)，2008.YONG En-mi. Study on trajectory optimization and guidance approach for hypersonic glide-reentry vehicle［D］. Changsha:National University of Defense Technology，2008.(in Chinese)

［10］姚寅偉，李華濱. 基于Gauss 偽譜法的高超聲速飛行器多約束三維再入軌跡優(yōu)化［J］.航天控制，2012，30(2):33 -45.YAO Yin-wei，LI Hua-bin. The generation of three-dimensional constrained entry trajectories for hypersonic vehicle based on the gauss pseudospectral method ［J］. Aerospace Control，2012，30(2):33 -45.(in Chinese)

［11］Zhao J，Zhou R. Reentry trajectory optimization for hypersonic vehicle satisfying complex constraints［J］. Chinese Journal of Aeronautics，2013，26(6):1544 -1553.

［12］Sun Y，Zhang M. Optimal reentry range trajectory of hypersonic vehicle by Gauss pseudospectral method［C］∥International Conference on Intelligent Control and Information Processing. Harbin:IEEE，2011:545 -549.

［13］Liu P，Zhao J S，Gu L X. Rapid approach for three-dimensional trajectory optimization of hypersonic reentry vehicles［J］. Flight Dynamics，2012，30(3):263 -268.

［14］Darby C L，Hager W W，Rao A V. An HP-adaptive pseudospectral method for solving optimal control problems［J］. Optimal Control Applications and Methods，2011，32(4):476 -502.

［15］國海峰，黃長強(qiáng)，丁達(dá)理，等. 多約束條件下高超聲速導(dǎo)彈再入軌跡優(yōu)化［J］.彈道學(xué)報(bào)，2013，25(1):10 -15.GUO Hai-feng，HUANG Chang-qiang，DING Da-li，et al. Reentry trajectory optimization for supersonic missile considering multiple constraints ［J］. Journal of Ballistics，2013，25(1):10 -15.(in Chinese)

［16］王麗英，張友安，趙國榮. 改進(jìn)的HP 自適應(yīng)網(wǎng)格細(xì)化算法及應(yīng)用［J］.彈道學(xué)報(bào)，2013，25(1):16 -21.WANG Li-ying，ZHANG You-an，ZHAO Guo-rong. Improved HP-adaptive mesh refinement algorithm and its application［J］.Journal of Ballistics，2013，25(1):16 -21.(in Chinese)

［17］洪蓓，辛萬青. HP 自適應(yīng)偽譜法在滑翔彈道快速優(yōu)化中的應(yīng)用［J］.計(jì)算機(jī)測量與控制，2012，20(5):1283 -1286.HONG Bei，XIN Wan-qing. Application of HP-adaptive pseudospectral method in rapid gliding trajectory optimization［J］. Computer Measurement ＆ Control，2012，20(5):1283 -1286.(in Chinese)

［18］趙漢元. 飛行器再入動(dòng)力學(xué)和制導(dǎo)［M］. 長沙:國防科學(xué)技術(shù)大學(xué)出版社，1997:84 -85.ZHAO Han-yuan. Reentry dynamics and guidance of flight vehicle［M］. Changsha:National University of Defense Technology Press，1997:84 -85.(in Chinese)

［19］Betts J T. Practical methods for optimal control and estimation using nonlinear programming［M］. Philadelphia:Society for Industrial and Applied Mathematics，2009:247 -251.

［20］賈沛然，陳克俊，何力. 遠(yuǎn)程火箭彈道學(xué)［M］. 長沙:國防科學(xué)技術(shù)大學(xué)出版社，2009:155 -156.JIA Pei-ran，CHEN Ke-jun，HE Li. Long-distance missile ballistics［M］. Changsha:National University of Defense Technology Press，2009:155 -156.(in Chinese)

［21］文葉，姜文志，馬登武. 飛機(jī)低空突防中的航路規(guī)劃技術(shù)研究［J］. 飛行力學(xué)，2006，24 (2):22 -26.WEN Ye，JIANG Wen-zhi，MA Deng-wu. Study on the plane low-altitude penetration route planning technique［J］. Flight Dynamics，2006，24 (2):22 -26.(in Chinese)

［22］熊丹君，蔡滿意，劉宇坤. 多約束條件下飛行器航路規(guī)劃［J］. 彈箭與制導(dǎo)學(xué)報(bào)，2009，29 (2):289 -292.XIONG Dan-jun，CAI Man-yi，LIU Yu-kun. Air vehicle route planning for multi-constraints［J］. Journal of Projectiles，Rockets，Missiles and Guidance，2009，29 (2):289-292.(in Chinese)