陳 琦，王中原，常思江

（南京理工大學(xué) 能源與動力工程學(xué)院，南京210094）

隨著戰(zhàn)爭理念的進步和發(fā)展，特別是現(xiàn)代非對稱作戰(zhàn)中“遠程壓制、精確打擊、高效毀傷”等概念的提出，傳統(tǒng)火炮射程不遠、精度不高等缺點逐漸暴露，已不能適應(yīng)現(xiàn)代戰(zhàn)爭的需求。在這種背景下，滑翔增程制導(dǎo)炮彈以其反應(yīng)速度快、效費比高、使用靈活并且可打擊靜止和運動目標(biāo)等優(yōu)勢逐漸受到各國的重視?；鑹埑虖椨苫鹋诎l(fā)射后，尾翼張開保持穩(wěn)定飛行，之后火箭發(fā)動機啟動并開始助推，發(fā)動機工作結(jié)束后，彈丸像普通尾翼彈一樣繼續(xù)在升弧段上飛行，到達彈道頂點附近，鴨舵展開，開始對彈丸實施控制，最終命中目標(biāo)。文獻［1－4］闡述了滑翔増程彈的飛行原理；文獻［5］研究了制導(dǎo)炮彈最優(yōu)滑翔彈道與控制問題，分別采用龐特里亞金極大值原理和共軛梯度法求解了控制變量有約束的最優(yōu)滑翔彈道；文獻［6］對滑翔彈道分3段進行了分析，研究了滑翔増程炮彈的彈道模型，并對滑翔彈道進行了設(shè)計。文獻［7］以最大升阻比的方法對滑翔彈道進行了優(yōu)化設(shè)計；文獻［8］分別采用了序列二次規(guī)劃法和最大升阻比法優(yōu)化了俯仰舵偏角規(guī)律。以上文獻主要以射程最遠為性能指標(biāo)優(yōu)化方案彈道，這對于考察彈丸的滑翔能力、打擊距離均具有十分重要的意義。但是這種方法也往往使得優(yōu)化出的彈道飛行時間過長，易被敵方防御系統(tǒng)發(fā)現(xiàn)和攔截；存速較低，削弱了末端攻擊的機動性。這在很大程度上限制了制導(dǎo)炮彈快速精確打擊、及時火力支援等能力的發(fā)揮。此外，上述研究均假設(shè)攻角能瞬時響應(yīng)指令信號，導(dǎo)致設(shè)計出的方案彈道過于理想，給實際應(yīng)用帶來了較大的難度。針對以上問題，本文以飛行時間為性能指標(biāo)對滑翔彈道進行了優(yōu)化設(shè)計，考慮了動力學(xué)滯后，將其簡化為一階慣性環(huán)節(jié)，引入了虛擬控制量，并將其作為優(yōu)化變量，同時為了保證攻擊效果，對末端存速、落角進行一定的限制，在縱向平面內(nèi)建立了彈道優(yōu)化模型；引入了偽譜理論對模型進行離散，這在很大程度上降低了彈道優(yōu)化耗時，縮短了滑翔彈的作戰(zhàn)反應(yīng)時間，同時也為今后實現(xiàn)算法的實時性提供了理論依據(jù)。

1 彈道優(yōu)化問題描述

1.1 滑翔増程彈運動方程

為研究問題的本質(zhì)，本文只考慮炮彈在縱向平面內(nèi)的運動，并將動力學(xué)滯后簡化為一階慣性環(huán)節(jié)，大氣為標(biāo)準(zhǔn)氣象條件。在這種情況下，滑翔増程制導(dǎo)炮彈的滑翔段運動方程可簡化為

Fx，F(xiàn)y分別為阻力和升力，有：

式中：θ為彈道傾角，q為動壓，Sref為參考面積，k和k′分別為彈翼組合體和舵面的誘導(dǎo)阻力系數(shù)，Cx0BW和Cx0δ分別為彈翼組合體和舵面的零升阻力系數(shù)，和Cδz

y分別為彈翼組合體和舵面的升力系數(shù)導(dǎo)數(shù)，以上各參數(shù)的值由吹風(fēng)實驗獲得。αB和δzB分別為平衡攻角和平衡舵偏角，Td為一階慣性環(huán)節(jié)的時間常數(shù)，uα為虛擬控制量，ε（uα）＝0為虛擬量uα的控制方程。方程中其它變量的意義參見文獻［9］。

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

為了實現(xiàn)快速打擊的要求，滑翔増程制導(dǎo)炮彈在彈道最高點處張開鴨舵，對彈體實施控制，經(jīng)過滑翔飛行，在最短的時間內(nèi)命中目標(biāo)，同時，為了保證對目標(biāo)的毀傷效果，還要對彈道末端速度和落角進行一定的限制。因此，滑翔彈道優(yōu)化問題可以具體描述為：在［0，tf］時間內(nèi)（tf為未知的彈道末端時刻），確定出控制量uα，使得由式（1）確定的動態(tài)系統(tǒng)從初始狀態(tài)轉(zhuǎn)移到終端狀態(tài)，在滿足規(guī)定約束的條件下，使得參數(shù)tf最小。由以上分析可知，彈道優(yōu)化問題實際上是一個動態(tài)優(yōu)化問題，其優(yōu)化模型可表示如下：

①性能指標(biāo)。為了實現(xiàn)快速打擊的要求，性能指標(biāo)取為攻擊時間最短，即：

②初始狀態(tài)為無控段末端狀態(tài)，本文取為彈道最高點處狀態(tài)。

③終端約束條件。為了保證最終的攻擊效果，滑翔段終端狀態(tài)應(yīng)滿足一定的約束，包括終端位置約束、速度約束和落角約束，即：

④攻角約束。考慮到彈體的控制性能，同時避免彈體失速，需要對攻角的幅值進行約束，即：

⑤系統(tǒng)狀態(tài)方程為式（1）描述的彈體縱向平面運動模型。

可見，滑翔彈道優(yōu)化問題是一個起始狀態(tài)固定，終端狀態(tài)受約束，終端時刻自由的非線性動態(tài)優(yōu)化問題。

2 Gauss偽譜法

作為一種求解最優(yōu)控制問題的直接方法，Gauss偽譜法的主要思想是將連續(xù)的無限維最優(yōu)控制問題轉(zhuǎn)換為非線性規(guī)劃問題進行求解，在一系列的Legendre－Gauss（LG）節(jié)點上將狀態(tài)變量和控制變量進行離散，并且利用這些離散點構(gòu)造全局Lagrange插值多項式來近似系統(tǒng)的動力學(xué)方程。該方法具有收斂速度快、對初值不敏感且無需猜測協(xié)態(tài)變量等優(yōu)點，2005年經(jīng)Benson D在其博士論文中介紹并完善后［10］，便受到了廣泛的關(guān)注和應(yīng)用。

2.1 系統(tǒng)動力學(xué)方程離散

考慮積分形式的系統(tǒng)動力學(xué)方程：

式中：X（t）∈Rn為狀態(tài)向量，U（t）∈Rm為控制向量，f：Rn×Rm×R→Rn為連續(xù)向量函數(shù)。采用Gauss偽譜法需要將時間區(qū)域［t0，tf］轉(zhuǎn)換到［－1，1］上，為此，引入變量τ對時間t進行變換［6］，即：

因此，系統(tǒng)動力學(xué)方程變?yōu)?/p>

由數(shù)值分析的知識可知，選取N個Lagrange插值基函數(shù)Lk（ζ）近似式（4）中的積分項，其代數(shù)精度可達2N－1次［11］，有：

這樣，結(jié)合式（4）和式（5），第i個節(jié)點處的狀態(tài)向量便可表示為

若記X0＝X（－1），Ui＝U（τi），Xi＝X（τi），1≤i≤N，根據(jù)式（6），便可將積分形式的系統(tǒng)動力學(xué)方程在N個LG節(jié)點上轉(zhuǎn)換為代數(shù)方程，即：

式中：Wi為第i個節(jié)點處的積分權(quán)重，Wi＝

2.2 邊界條件和性能指標(biāo)離散

式（8）僅在區(qū)間的內(nèi)點計算狀態(tài)量，并未包含終端時刻節(jié)點，終端狀態(tài)應(yīng)滿足動力學(xué)方程約束，即：

將終端條件離散并利用Gauss積分來近似，可得：

由于本文所選取的性能指標(biāo)為Mayer型指標(biāo)，不包含積分項。因此無需進行特殊處理，即：

根據(jù)以上的數(shù)學(xué)變換，本文的問題可以描述為：在［t0，tf］時間內(nèi)（tf為未知的彈道末端時刻），確定離散點上的狀態(tài)向量Xi（i＝1，2，…，N）、控制向量Ui（i＝1，2，…，N）和終端時刻tf，使得性能指標(biāo)（11）最小，并滿足系統(tǒng)動力學(xué)方程約束（8）、終端狀態(tài)約束（10），以及原問題的邊界條件和過程約束：

從而將連續(xù)的無限維最優(yōu)控制問題轉(zhuǎn)換為一般的非線性規(guī)劃問題，寫成標(biāo)準(zhǔn)形式為

3 數(shù)值仿真算例

以某滑翔增程制導(dǎo)炮彈為對象進行數(shù)值仿真。經(jīng)火箭發(fā)動機助推后彈丸質(zhì)量m＝37.27kg，炮彈達到彈道最高點時的高度為Hp0＝14 746.6m，速度為vp0＝299.3m／s，彈道傾角為θp0＝0；滑翔段舵偏角約束為｜αB（t）｜≤15°，滑翔末端高度取為yf＝0，彈道傾角取為｜θ（tf）｜≥θf＝89°，同時，為了保證攻擊效果，也要對末端存速進行限制，v（tf）≥vf＝200m／s，一階慣性環(huán)節(jié)的時間常數(shù)Td取為1.5s。取彈道最高點處x0＝0，對目標(biāo)xT＝40km，50km，60km分別進行仿真計算，LG節(jié)點個數(shù)取為35，優(yōu)化算法工具包為SNOPT，仿真結(jié)果如圖1～圖5所示。

從圖1可見，制導(dǎo)炮彈在經(jīng)過彈道最高點后，近乎以直線滑翔，當(dāng)快達到末端時，彈道開始俯沖，這是由于為了保證最后的攻擊效果，對末端條件進行約束所致。圖2描述了滑翔過程中速度的變化規(guī)律，可以發(fā)現(xiàn)，速度在中間階段變化非常緩慢，但當(dāng)接近彈道末端時開始衰減，這是因為彈丸在末段以較大的負(fù)攻角進行俯沖，增大了飛行阻力，因此速度衰減較快。

圖1 不同射程情況下的彈道曲線

圖2 不同射程情況下的速度曲線

圖3 也很明顯地反映出彈道傾角在經(jīng)過短暫的減小后，在較長的時間內(nèi)近乎保持不變，驗證了彈丸近似直線飛行，當(dāng)接近彈道末端時，彈道傾角先增加后急劇減小，說明彈丸開始俯沖，且最終落角均大于89°，滿足設(shè)計指標(biāo)。

圖3 不同射程情況下的彈道傾角曲線

圖4 為不同射程條件下虛擬控制量uα和實際控制量αB的對比曲線，可以看出，因為考慮了動力學(xué)滯后，為了使攻角迅速增加，uα在初始時刻的值較大，之后慢慢減小，逐漸與αB重合。由于在構(gòu)建優(yōu)化模型時對末端彈道傾角進行了約束，｜θ（tf）｜≥θf，為了滿足這種約束，升力需在彈道末端變?yōu)樨?fù)值，根據(jù)升力計算公式，此時攻角也應(yīng)該為負(fù)值，且θf越大，攻角的絕對值也越大，從而攻角曲線在末端出現(xiàn)了圖4所示的陡峭現(xiàn)象。此外，由于動力學(xué)滯后的存在，uα在末端的下降過程中要超前于αB。

圖4 不同射程情況下的實際控制量和虛擬控制量曲線

圖5 展示了考慮動力學(xué)滯后的攻角改善效果?？梢悦黠@地看出，考慮動力學(xué)滯后的攻角變化較為平緩，以xT＝60km為例，不考慮動力學(xué)滯后時，攻角在初始時刻從0瞬間增加到接近6°，這種突變往往會給今后控制系統(tǒng)的準(zhǔn)確跟蹤帶來較大的困難。而在考慮動力學(xué)滯后時，這種現(xiàn)象便不再出現(xiàn)，優(yōu)化得到的攻角曲線變化非常平緩。

圖5 考慮動力學(xué)滯后攻角改善效果

為了進一步驗證Gauss偽譜法的高效性，利用直接打靶法［14］在同樣的計算平臺上分別對xT＝40km，50km，60km的彈道進行了優(yōu)化計算，離散控制量所需的節(jié)點設(shè)置為均勻分布，且個數(shù)取為35，積分步長取為0.5s，對于優(yōu)化模型中的各種約束，采用罰函數(shù)法進行處理，優(yōu)化結(jié)果對比見表1。從中可看出，直接打靶法優(yōu)化出的彈道飛行時間大于Gauss偽譜法，這說明其尋優(yōu)效果要低于Gauss偽譜法。此外，直接打靶法的優(yōu)化結(jié)果并沒有完全滿足落角約束，可見，在終端位置、速度和落角等諸多約束的條件下，相比于Gauss偽譜法，直接打靶法的精度有所欠缺。最后，直接打靶法的優(yōu)化耗時tcost要遠大于Gauss偽譜法，這是因為直接打靶法在每次迭代時都會利用控制量積分彈道，雖然選取了較大的積分步長（0.5s），但是大量的迭代還是產(chǎn)生了巨大的計算耗時，明顯地降低了算法的效率，而Gauss偽譜法無需在每次迭代中積分彈道，因此在計算效率上具有很大的優(yōu)勢。圖6為兩種方法優(yōu)化得到的攻角的對比圖，從中可看出二者基本吻合，但Gauss偽譜法的結(jié)果要更為光滑。從以上的分析可以看出，Gauss偽譜法不論在優(yōu)化效率還是精度上均強于直接打靶法。

表1 直接打靶法和Gauss偽譜法優(yōu)化結(jié)果對比

圖6 直接打靶法和Gauss偽譜法攻角對比曲線

4 結(jié)束語

本文基于Gauss偽譜法求解了遠程制導(dǎo)炮彈滑翔彈道快速優(yōu)化問題。以快速精確打擊為指標(biāo)對滑翔彈道進行了優(yōu)化設(shè)計，考慮動力學(xué)滯后，引入虛擬控制量，并將其作為優(yōu)化變量，在縱向平面內(nèi)建立了彈道優(yōu)化模型；利用Gauss偽譜法將優(yōu)化模型在一系列LG節(jié)點上離散，將連續(xù)的無限維最優(yōu)控制問題轉(zhuǎn)換為一般的非線性規(guī)劃問題。仿真結(jié)果表明，彈丸能以一定的速度和大著角攻擊目標(biāo)，所有約束條件均能較好地滿足，引入虛擬控制量后，優(yōu)化得到的攻角曲線更為平滑。與傳統(tǒng)的直接打靶法相比，Gauss偽譜法計算時間短，優(yōu)化效率高，具有在線優(yōu)化的潛力，促進了制導(dǎo)炮彈快速精確打擊、及時火力支援等能力的有效發(fā)揮，可為遠程制導(dǎo)炮彈的方案彈道優(yōu)化設(shè)計提供一定的參考。

［1］COSTELLO M.Extended range of a gun launched smart projectile using controllable canards［J］.Journal of Guidance，Control，and Dynamics，2006，29（6）：1 404－1 414.

［2］COSTELLO M.Extended range of a gun launched smart projectile using controllable canards［J］.Journal of Shock and Vibration，2001，8：203－213.

［3］史金光，王中原，易文俊.滑翔増程彈方案彈道特性研究［J］.彈道學(xué)報，2003，15（1）：51－54.SHI Jin－guang，WANG Zhong－yuan，YI Wen－jun.A study on the projectile trajectory characteristics of gliding extended range projectile［J］.Journal of Ballistics，2003，15（1）：51－54.（in Chinese）

［4］史金光，王中原，易文俊.滑翔増程彈飛行彈道［J］.火力與指揮控制，2007，32（11）：88－90.SHI Jin－guang，WANG Zhong－yuan，YI Wen－jun.Study of fight trajectory for gliding extended range projectile［J］.Fire Control and Command Control，2007，32（11）：88－90.（in Chinese）

［5］卜奎晨，劉莉.復(fù)合制導(dǎo)炮彈最優(yōu)滑翔彈道與控制［J］.彈道學(xué)報，2007，19（4）：23－25.BU Kui－chen，LIU Li.Optimal glide trajectory and control for combined guidance munition［J］.Journal of Ballistics，2007，19（4）：23－25.（in Chinese）

［6］汪小娜，王樹宗，朱華兵.火箭助推滑翔増程炮彈彈道模型研究［J］.艦船科學(xué)技術(shù)，2005，27（1）：77－79.WANG Xiao－na，WANG Shu－zong，ZHU Hua－bing.Study on ballistic models of rocket－aid gliding extended projectiles［J］.Ship Science and Technology，2005，27 （1）：77 － 79.（in Chinese）

［7］易文俊，王中原，史金光，等.帶鴨舵滑翔増程炮彈方案彈道研究［J］.南京理工大學(xué)學(xué)報，2008，32（3）：322－326.YI Wen－jun，WANG Zhong－yuan，SHI Jin－guang.Trajectory of glide extended－range projectile with canards configuration［J］.Journal of Nanjing University of Science and Technology，2008，32（3）：322－326.（in Chinese）

［8］史金光，王中原，涂泗華.求解滑翔増程彈較優(yōu)舵偏角方法［J］.彈箭與制導(dǎo)學(xué)報，2005，25（4）：717－719.SHI Jin－guang，WANG Zhong－yuan，TU Si－h(huán)ua.Calculating method on a suitable angle of lift－canards of gliding extended range projectile［J］.Journal of Projectiles，Rockets，Missiles and Guidance，2005，25（4）：717－719.（in Chinese）

［9］錢杏芳，林瑞雄，趙亞男.導(dǎo)彈飛行力學(xué)［M］.北京：北京理工大學(xué)出版社，2008.QIAN Xing－fang，LIN Rui－xiong，ZHAO Ya－nan.Missile flight dynamics［M］.Beijing：Beijing Institute of Technology Press，2008.（in Chinese）

［10］BENSON D.A Gauss pseudospectral transcription for optimal control［D］.Cambridge：Massachusetts Institute of Technology，2005.

［11］關(guān)冶，陸金甫.數(shù)值分析基礎(chǔ)［M］.北京：高等教育出版社，2005.GUAN Ye，LU Jin－fu.Fundamentals of numerical analysis［M］.Beijing：High Education Press，2005.（in Chinese）

［12］GILL P E，MURRAY W，SAUNDERS M A.SNOPT：an SQP algorithm for large－scale constrained optimization［J］.SIAM Review，2005，47（1）：99－131.

［13］ACHTER A W，BIEGLER L T.On the implementation of a primal－dual interior point filter line search algorithm for largescale nonlinear programming［J］.Mathematical Programming，2006，106（1）：25－57.

［14］GERDTS M.Direct shooting method for the numerical solution of higher－index DAE optimal control problems［J］.Journal of Optimization Theory and Applications，2003，17（2）：267－294.