朱大林，唐勝景，郭 杰，熊俊輝

（北京理工大學(xué) 宇航學(xué)院，北京100081）

對(duì)于具有無控起飛段的導(dǎo)彈，實(shí)際飛行中由于受到各種隨機(jī)干擾的影響，其啟控點(diǎn)彈道參數(shù)是隨機(jī)分布的，從而形成相對(duì)于理論啟控點(diǎn)的散布。為了保證導(dǎo)彈受到良好的控制，啟控點(diǎn)散布應(yīng)滿足制導(dǎo)系統(tǒng)開始導(dǎo)引的要求。如對(duì)于激光駕束制導(dǎo)的反坦克導(dǎo)彈，一般要求啟控點(diǎn)要落入規(guī)定的光束截面內(nèi)。導(dǎo)彈啟控點(diǎn)的選擇一般是期望啟控點(diǎn)落入規(guī)定的區(qū)域內(nèi)，并期望啟控點(diǎn)波動(dòng)越小越好，這兩點(diǎn)一般可以通過啟控點(diǎn)彈道參數(shù)的期望和方差進(jìn)行描述。為達(dá)到這2個(gè)目標(biāo)，通常需要不斷地進(jìn)行統(tǒng)計(jì)分析和篩選，以確定合適的啟控距離或時(shí)間以及相關(guān)總體參數(shù)（如初速度、發(fā)射角等）。文獻(xiàn)[1－2]采用Monte Carlo法對(duì)反坦克導(dǎo)彈的啟控點(diǎn)散布進(jìn)行了研究，其中文獻(xiàn)[2]通過多次仿真確定了最佳的發(fā)射角；文獻(xiàn)[3]采用協(xié)方差分析描述函數(shù)法建立彈箭外彈道的隨機(jī)優(yōu)化模型，并通過優(yōu)化算法確定了偏角方差最小時(shí)有關(guān)總體參數(shù)的取值。上述方法雖然都可以用于導(dǎo)彈的啟控點(diǎn)選擇，但其不足之處為：統(tǒng)計(jì)分析效率比較低，篩選過程比較費(fèi)機(jī)時(shí)；優(yōu)化目標(biāo)中沒有同時(shí)考慮優(yōu)化條件對(duì)散布均值和方差的影響。

針對(duì)上述問題，本文采用穩(wěn)健優(yōu)化的概念和方法處理導(dǎo)彈啟控點(diǎn)的選擇問題，綜合考慮啟控距離（時(shí)間）和總體參數(shù)的取值對(duì)啟控點(diǎn)散布的影響。穩(wěn)健優(yōu)化（Robust Optimization）是產(chǎn)品質(zhì)量工程中保證產(chǎn)品質(zhì)量穩(wěn)定的一種重要的設(shè)計(jì)方法，其基本思想就是通過優(yōu)化算法在可行域內(nèi)尋找最佳的設(shè)計(jì)方案，力圖使產(chǎn)品的質(zhì)量特性均值盡可能地接近目標(biāo)值，同時(shí)對(duì)干擾因素引起的波動(dòng)盡可能地小[4]。可見，穩(wěn)健優(yōu)化與導(dǎo)彈啟控點(diǎn)選擇的目的是一致的。本文根據(jù)穩(wěn)健優(yōu)化響應(yīng)面建模的思想首先對(duì)導(dǎo)彈擾動(dòng)運(yùn)動(dòng)系統(tǒng)建立高斯型徑向基網(wǎng)絡(luò)（Gaussian Radial Basis Network，GRBN），并以GRBN為基礎(chǔ)給出不確定性分析的解析表達(dá)式；然后以某反坦克導(dǎo)彈啟控點(diǎn)選擇的要求為例建立穩(wěn)健優(yōu)化模型，并采用小生境Pareto遺傳算法獲得該模型的Pareto最優(yōu)解集。對(duì)最優(yōu)解進(jìn)行Monte Carlo模擬打靶驗(yàn)證表明，優(yōu)化結(jié)果比較令人滿意，從而為彈道設(shè)計(jì)、制導(dǎo)系統(tǒng)設(shè)計(jì)等提供參考和依據(jù)。

1 可控因素和噪聲因素的確定

影響導(dǎo)彈啟控點(diǎn)散布的因素有很多種，如初始擾動(dòng)、推力偏心、彈體質(zhì)量分布不對(duì)稱、風(fēng)、制造和安裝工藝誤差等。一般認(rèn)為這些因素是隨機(jī)分布的，事先可以根據(jù)試驗(yàn)數(shù)據(jù)或經(jīng)驗(yàn)確定其分布類型。為了區(qū)別對(duì)待這些因素，根據(jù)穩(wěn)健優(yōu)化的概念可將其分為兩大類：可控因素和噪聲因素。可控因素就是可以選擇其水平即名義值的因素，也是設(shè)計(jì)變量，通常由設(shè)計(jì)名義值和設(shè)計(jì)容差進(jìn)行描述；噪聲因素就是不可控的因素，通常由其分布類型及數(shù)字特征進(jìn)行描述。因此，導(dǎo)彈啟控點(diǎn)選擇的穩(wěn)健優(yōu)化就是在約束范圍內(nèi)確定可控因素的最佳水平組合。

以某反坦克導(dǎo)彈的啟控點(diǎn)選擇為例，假設(shè)初步考慮的因素有：初速度v0，發(fā)射角θ0，初始滾轉(zhuǎn)角速度ωx0，加速發(fā)動(dòng)機(jī)點(diǎn)火時(shí)間tb，啟控距離xc，初始彈道偏角ψv0，慣性中心主軸在彈體坐標(biāo)系下的方位角λ和μ以及橫風(fēng)wz。如果其水平可以選擇的因素有初速度v0，發(fā)射角θ0，初始滾轉(zhuǎn)角速度ωx0，加速發(fā)動(dòng)機(jī)點(diǎn)火時(shí)間tb以及啟控距離xc，則確定的可控因素x和噪聲因素ω分別為x＝（v0θ0ωx0tbxc），ω＝（ψv0λμ wz）。這些因素的詳細(xì)描述分別如表1和表2所示。

表1 可控因素

表2 噪聲因素

2 擾動(dòng)運(yùn)動(dòng)系統(tǒng)建模和不確定性分析

確定了影響因素之后，接下來就是對(duì)導(dǎo)彈擾動(dòng)運(yùn)動(dòng)系統(tǒng)進(jìn)行統(tǒng)計(jì)分析，從而為第3節(jié)的穩(wěn)健優(yōu)化做準(zhǔn)備。導(dǎo)彈擾動(dòng)運(yùn)動(dòng)系統(tǒng)的數(shù)學(xué)模型可參考文獻(xiàn)[5]，常用的分析方法有 Monte Carlo法和協(xié)方差分析描述函數(shù)法。Monte Carlo法是啟控點(diǎn)散布分析方法中一種簡單而有效的方法，但為了獲得精確的統(tǒng)計(jì)值，通常需要幾百次甚至上千次的模擬打靶，比較費(fèi)機(jī)時(shí)。協(xié)方差分析描述函數(shù)法是運(yùn)用描述函數(shù)理論將原來的非線性隨機(jī)系統(tǒng)進(jìn)行統(tǒng)計(jì)線性化，通過求解導(dǎo)出的關(guān)于狀態(tài)變量均值和協(xié)方差的傳播方程，從而獲得相應(yīng)的統(tǒng)計(jì)值[6]。這種方法可以明顯提高統(tǒng)計(jì)分析的效率，但目前統(tǒng)計(jì)線性化的理論主要是建立在狀態(tài)變量服從聯(lián)合正態(tài)分布的假設(shè)上，當(dāng)原系統(tǒng)非線性程度比較高或狀態(tài)變量嚴(yán)重偏離正態(tài)分布時(shí)，協(xié)方差分析描述函數(shù)法的分析精度會(huì)明顯降低；而且，這種方法求解的微分方程比較多，不利于后續(xù)的穩(wěn)健優(yōu)化?；谏鲜鲇懻?，借鑒穩(wěn)健優(yōu)化響應(yīng)面建模的思想，對(duì)導(dǎo)彈擾動(dòng)運(yùn)動(dòng)系統(tǒng)采用GRBN進(jìn)行擬合，并在此基礎(chǔ)上給出不確定性分析的解析表達(dá)式，這種方法可以顯著提高統(tǒng)計(jì)分析的效率。

2.1 GRBN建模

徑向基函數(shù)網(wǎng)絡(luò)是一種簡單而有效的建模工具，其中GRBN由于可以寫成便于概率積分的張量積函數(shù)形式[7]，因而在文中被選用。GRBN的數(shù)學(xué)模型可寫為

式中，‖·‖為歐氏范數(shù)，x＝（x1x2…xn），為n維輸入向量；x1，x2…，xm為徑向基網(wǎng)絡(luò)中心點(diǎn)；α0為常數(shù)項(xiàng)；α＝（α1α2…αm），為權(quán)系數(shù)向量，c為正的常實(shí)數(shù)。

GRBN的構(gòu)造就是要確定網(wǎng)絡(luò)的結(jié)構(gòu)和參數(shù)，其訓(xùn)練過程可參考一般徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)的學(xué)習(xí)算法[8－9]。對(duì)于本文的導(dǎo)彈擾動(dòng)運(yùn)動(dòng)系統(tǒng)，GRBN的構(gòu)造過程為：首先選用合適的試驗(yàn)設(shè)計(jì)獲得以可控因素和噪聲因素為試驗(yàn)點(diǎn)、以啟控點(diǎn)參數(shù)為輸出的樣本；其次，對(duì)樣本進(jìn)行分析，剔除一些離群點(diǎn)；最后，隨機(jī)選擇一部分樣本試驗(yàn)點(diǎn)作為GRBN的中心點(diǎn)并通過Levenberg－Marquardt優(yōu)化算法極小化樣本的誤差平方和，從而確定GRBN的α0，α和c。GRBN的擬合精度可以通過重采樣樣本進(jìn)行檢驗(yàn)，常用的檢驗(yàn)準(zhǔn)則有相對(duì)平均絕對(duì)誤差檢驗(yàn)準(zhǔn)則和相對(duì)最大誤差檢驗(yàn)準(zhǔn)則[10]。

2.2 基于GRBN的不確定分析

GRBN滿足擬合精度要求后，就可以對(duì)其進(jìn)行不確定性分析。將式（1）改寫為

式中：xj為點(diǎn)x的第j個(gè)分量，xij為點(diǎn)xi的第j個(gè)分量。令

則

假設(shè)可控因素、噪聲因素彼此相互獨(dú)立，則根據(jù)文獻(xiàn)[7]的推導(dǎo)，式（3）輸出的均值和方差可表示為

對(duì)于形式為式（1）的GRBN，經(jīng)過推導(dǎo)可得到如下結(jié)果：

若因素xj服從均勻分布，即

則

式中：Φ（·）為標(biāo)準(zhǔn)正態(tài)分布的分布函數(shù)。

若因素xj服從正態(tài)分布，即

則

3 啟控點(diǎn)選擇的穩(wěn)健優(yōu)化模型與求解

與確定性優(yōu)化模型不同，穩(wěn)健優(yōu)化模型一般是采用概率論的術(shù)語進(jìn)行描述的，比較常用的就是隨機(jī)變量的均值和方差。下面以第1節(jié)提到的反坦克導(dǎo)彈為例，建立啟控點(diǎn)選擇的穩(wěn)健優(yōu)化模型。

通常，反坦克導(dǎo)彈的發(fā)射角不太大，初始段彈道比較低。在選擇啟控距離和總體參數(shù)時(shí)，一般期望啟控點(diǎn)的方向偏差和散布半徑比較小，并期望啟控點(diǎn)有一定的高度且彈道角落入一定的范圍內(nèi)，以防止其掉地或啟控后飛離視場范圍。根據(jù)這些要求，建立的穩(wěn)健優(yōu)化模型為

式中：E（·）表示隨機(jī)變量的期望；μxk為可控因素xk的名義值；yl，yu為高度界限；θl，θu為彈道角界限為可控因素名義值的界限；RCEP為圓概率誤差半徑，可采用如下公式進(jìn)行計(jì)算[11]：

式中：σy，σz分別為啟控點(diǎn)的高度和方向均方差。

需要指出的是，如果希望穩(wěn)健優(yōu)化模型中某個(gè)約束具有一定的穩(wěn)健性，可以通過引入方差進(jìn)行調(diào)整。例如：如果希望高度約束具有一定的穩(wěn)健性，可以將其修改為yl＋nlσy≤E（y）≤yu－nuσy，式中：nl，nu稱為σ水平。

求解上面建立的穩(wěn)健優(yōu)化模型，其實(shí)就是求解一個(gè)多目標(biāo)優(yōu)化問題。多目標(biāo)優(yōu)化問題的求解策略大致分為兩類，可將其轉(zhuǎn)化為單目標(biāo)優(yōu)化問題或者直接求解多目標(biāo)優(yōu)化問題，每種優(yōu)化策略又有不同的優(yōu)化方法。這里采用小生境Pareto遺傳算法直接求解建立的穩(wěn)健優(yōu)化模型。該算法采用Pareto支配關(guān)系進(jìn)行錦標(biāo)賽選擇，并采用小生境技術(shù)實(shí)現(xiàn)群體的多樣性，能較好地獲得穩(wěn)健優(yōu)化模型的Pareto最優(yōu)解集，從而方便設(shè)計(jì)者根據(jù)自己的偏好進(jìn)行挑選，算法的具體描述可參考文獻(xiàn)[12]。

4 算例仿真

下面仍以第1節(jié)提到的反坦克導(dǎo)彈啟控點(diǎn)選擇為例，假設(shè)瞄準(zhǔn)線水平，根據(jù)第3節(jié)建立的具體穩(wěn)健優(yōu)化模型為

為了利用式（4）和式（5）進(jìn)行不確定性分析，首先需要對(duì)導(dǎo)彈擾動(dòng)運(yùn)動(dòng)系統(tǒng)構(gòu)造GRBN。這里選用拉丁超方試驗(yàn)設(shè)計(jì)選擇500個(gè)試驗(yàn)點(diǎn)，通過對(duì)擾動(dòng)運(yùn)動(dòng)系統(tǒng)的仿真獲得啟控點(diǎn)坐標(biāo)y，z和彈道角θ的取值。隨后通過分析樣本剔除一些離群點(diǎn)，形成樣本容量為491的樣本。在構(gòu)造以啟控點(diǎn)坐標(biāo)y，z和彈道角θ為輸出的GRBN時(shí)隨機(jī)選擇75%的試驗(yàn)點(diǎn)作為GRBN的中心點(diǎn)。GRBN的擬合精度通過重采樣300個(gè)樣本點(diǎn)進(jìn)行檢驗(yàn)，其相對(duì)平均絕對(duì)誤差均小于7%。

GRBN滿足使用要求后，接著采用式（11），優(yōu)化算法設(shè)置5個(gè)子種群，共500個(gè)個(gè)體，種群進(jìn)化100代，最后共獲得181個(gè)Pareto非支配解，形成的Pareto前沿如圖1所示。

圖1 Pareto前沿

為驗(yàn)證優(yōu)化結(jié)果的可信性，從上述Pareto非支配解解集中隨機(jī)挑選出一個(gè)解對(duì)原擾動(dòng)運(yùn)動(dòng)系統(tǒng)進(jìn)行500次Monte Carlo模擬打靶驗(yàn)證。挑選的解為μv0＝24.15m/s，μθ0＝12.94°，μω0＝9.18r/s，μtb＝0.10s，μxc＝21.93m。

驗(yàn)證結(jié)果如表3所示，通過 Monte Carlo模擬打靶獲得的散布圖如圖2所示。

表3 驗(yàn)證結(jié)果

圖2 散布圖

由表3的驗(yàn)證結(jié)果可以看出，除了E（θ）的數(shù)值相對(duì)有所偏大之外，其余的統(tǒng)計(jì)值基本上比較接近，啟控點(diǎn)高度和彈道傾角的均值均在約束范圍之內(nèi)。從圖2的散布分布來看，這組解達(dá)到了預(yù)期的目的。

5 結(jié)論

本文采用穩(wěn)健優(yōu)化的思想和方法處理導(dǎo)彈啟控點(diǎn)選擇的問題，從仿真驗(yàn)證的結(jié)果來看，這種做法是可行的。文中采用小生境Pareto遺傳算法可以獲得啟控點(diǎn)穩(wěn)健優(yōu)化模型的Pareto最優(yōu)解集，從而為彈道設(shè)計(jì)者提供更多的參考和選擇。本文仍需改進(jìn)的地方就是在優(yōu)化過程中需要提高GRBN的局部擬合精度從而進(jìn)一步提高穩(wěn)健最優(yōu)解的可信度。

[1]LIU Yang，LI Kun，SONG Xiao－dong，et al.Monte Carlo ballistic simulation applied to dispersion analysis of starting control points for antitank missile[C]//Proceedings of 2011 International Conference on Electronic and Mechanical Engineer and Information Technology.Piscataway：IEEE Computer Society，2011：2 929－2 932.

[2]閆智強(qiáng)，周鳳歧，周軍.基于 Monte Carlo的激光駕束反坦克導(dǎo)彈 啟 控 點(diǎn) 散 布 研 究 [J].彈 箭 與 制 導(dǎo) 學(xué) 報(bào)，2005，26（1）：689－691.YAN Zhi－qiang，ZHOU Feng－qi，ZHOU Jun.The study to the spread of starting control points for laser beam－rider guidance anti－tank missile on Monte Carlo method[J].Journal of Projectiles，Rockets，Missiles and Guidance，2005，26（1）：689－691.（in Chinese）

[3]徐明友，丁松濱，徐直軍.外彈道優(yōu)化設(shè)計(jì)的隨機(jī)模型[J].南京理工大學(xué)學(xué)報(bào)，1997，21（5）：403－406.XU Ming－you，DING Song－bin，XU Zhi－jun.Random model for optimum design of exterior ballistics[J].Journal of Nanjing U－niversity of Science and Technology，1997，21（5）：403－406.（in Chinese）

[4]陳立周，于曉紅，翁海珊.基于隨機(jī)優(yōu)化的工程穩(wěn)健設(shè)計(jì)[J].北京科技大學(xué)學(xué)報(bào)，1999，21（1）：57－59.CHEN Li－zhou，YU Xiao－h(huán)ong，WENG Hai－shan.Robust engineering design by stochastic optimization[J].Journal of University of Science and Technology Beijing，1999，21（1）：57－59.（in Chinese）

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

[6]袁子懷，錢杏芳.有控飛行力學(xué)與計(jì)算機(jī)仿真[M].北京：國防工業(yè)出版社，2001：139－174.YUAN Zi－h(huán)uai，QIAN Xing－fang.Control flight mechanics and computer simulation[M].Beijing：National Defence Industry Press，2001：139－174.（in Chinese）

[7]CHEN Wei，JIN Rui－chen，SUDJIANTO A.Analytical variance－based global sensitivity analysis in simulation－based design under uncertainty[J].Journal of Mechanical Design，2005，127（5）：875－886.

[8]UYKAN Z，GUZELIS C，CELEBI M E.Analysis of inputoutput clustering for determining centers of RBFN[J].IEEE Transactions on Neural Networks，2000，11（4）：851－858.

[9]CHEN S，COWAN C F N，GRANT P M.Orthogonal least squares learning algorithm for radial basis function network[J].IEEE Transactions on Neural Networks，1991，2（2）：302－309.

[10]龍騰.飛行器多學(xué)科設(shè)計(jì)方法與集成設(shè)計(jì)平臺(tái)研究[D].北京：北京理工大學(xué)宇航學(xué)院，2009.LONG Teng.Research on methods of multidisciplinary design optimization and integrated design[D].Beijing：School of Aerospace Engineering，Beijing Institute of Technology，2009.（in Chinese）

[11]韓子鵬.彈箭外彈道學(xué)[M].北京：北京理工大學(xué)出版社，2008：113－116.HAN Zi－peng.Projectile and rocket exterior ballistics[M].Beijing：Beijing Institute of Technology Press，2008：113－116.（in Chinese）

[12]HORN J，NAFPLIOTIS N，GOLDBERG D E.A niched Pareto genetic algorithm for multiobjective optimization[C]//Proceedings of the First IEEE Conference on Evolutionary Computation.New York：IEEE，1994：82－87.