基于自適應混沌變異粒子群優化算法的旋轉彈丸氣動參數辨識

管軍，周家勝，易文俊，劉世平，常思江，史繼剛

(1.南京理工大學 瞬態物理國家重點實驗室, 江蘇 南京 210094；2.海軍駐沈陽彈藥專業軍事總代表室，遼寧 沈陽 110045；3.南京理工大學 能源與動力工程學院, 江蘇 南京 210094)

基于自適應混沌變異粒子群優化算法的旋轉彈丸氣動參數辨識

管軍1，周家勝2，易文俊1，劉世平3，常思江3，史繼剛1

(1.南京理工大學 瞬態物理國家重點實驗室, 江蘇 南京 210094；2.海軍駐沈陽彈藥專業軍事總代表室，遼寧 沈陽 110045；3.南京理工大學 能源與動力工程學院, 江蘇 南京 210094)

將最大似然準則應用于高速旋轉彈丸的氣動參數辨識問題中，提出一種新的自適應混沌變異粒子群算法求解該準則下的氣動參數最優解，進而得到彈丸的氣動參數。該算法通過自適應調整慣性權重、利用混沌優化的思想產生初始粒子、設定早熟判別機制來判斷是否陷入局部最優解，并通過粒子變異的策略使其跳出局部最優解等方法進一步優化基本粒子群算法。通過常用的測試函數對該算法進行了測試，測試結果表明：相比于基本粒子群算法，該算法具有收斂速度快、尋優精度高、應用范圍廣等優點。利用系統仿真的方法模擬彈丸的自由飛行數據，并利用該數據結合所提算法對彈丸的主要氣動參數進行辨識，辨識結果表明：該算法可以有效辨識彈丸的氣動參數，且精度高，收斂速度快，可以應用于工程實際問題。

兵器科學與技術；彈丸；氣動參數辨識；粒子群優化算法；最大似然準則

0 引言

圖1 全彈道坐標雷達Fig.1 Full range trajectory tracking radar

高速旋轉穩定彈的氣動參數辨識問題一直是系統辨識領域的難點，此類彈丸目前主要是通過風洞吹風和理論計算的方法來確定彈丸的氣動參數。對于傳統榴彈由于彈丸內沒有彈載測量傳感器，無法對彈丸的飛行姿態、位置、速度等信息進行實時測量和記錄，所以很難通過實際飛行數據辨識此類彈丸的氣動參數。但在實際工程應用中，可以通過彈道跟蹤雷達測量彈丸坐標、速度等數據，通過攻角- 紙靶的方法測量彈丸攻角運動數據。但通過紙靶所獲取的攻角數據通常含有較大的誤差，所以在確定攻角數據之前，通常要采用最優估計等技術對紙靶判讀數據進行平滑、濾波等預處理，從而得到更加精確的攻角數據。由此，可以得到位置、速度、姿態等離線飛行數據，并通過相關方法辨識彈丸氣動參數[1]。圖1和圖2分別為全彈道坐標跟蹤雷達以及攻角- 紙靶試驗彈丸穿孔后的實際情況。目前工程上應用最為廣泛的氣動力辨識準則是最大似然方法(MLE)，該準則廣泛應用于飛機、導彈、航天器等飛行器的氣動辨識問題中[2-5]，將氣動參數辨識問題轉換成最優化問題，通過優化求取氣動力模型參數，使模型的輸出與實測值之間的偏差達到最小。在旋轉穩定彈丸的氣動參數辨識中，Chapman等[6]將基于最小二乘思想的Champ-Kirk方法應用與彈丸的氣動參數辨識中。相比于最小二乘準則，最大似然準則更多地應用于多參數同時辨識問題，以求得多參數的全局最優解。利用最大似然準則辨識旋轉穩定彈丸的氣動參數目前鮮有公開報道，本文將最大似然準則應用于該類型彈丸的氣動參數辨識問題中。在求解該優化問題時，一般采用基于目標函數梯度構造的優化方法，如最速下降法、共軛梯度法等等，但這類方法具有相當的局限性，當系統存在非線性或時滯特性時，根本無法求取梯度值[7]。盡管Champ-Kirk方法可以處理非線性問題，但該方法在非線性問題線性化的過程中，存在一定的誤差，且該方法在旋轉穩定彈丸6自由度方程敏感系數的推導過程中非常復雜。同時，基于梯度構造的優化方法在優化迭代的過程中，始于一組特定的參數，使得優化出的解極易陷入起始點附近的局部最優解，從而給全局最優解的求解問題帶來一定的麻煩。粒子群優化(PSO)算法是由Kennedy等[8]提出的一種進化算法，它是通過追隨當前搜索到的最優值來尋找全局最優解。本文利用PSO算法求解氣動參數辨識問題，并通過對基本PSO算法的改進進一步提高PSO算法的優化效率。

根據彈箭外彈道學[9]的基礎理論可知，彈丸的線運動和氣動力系數相關性較大，而彈丸的角運動和氣動力矩系數相關性較大。本文通過模擬雷達測量的位置、速度數據來辨識主要的氣動力系數：零升阻力系數cx0、誘導阻力系數cxi和升力系數導數c′y；通過模擬紙靶試驗得到攻角數據來辨識主要氣動力矩系數：翻轉力矩(旋轉彈丸中為靜不穩定力矩)系數導數m′z. 其中系統模型中的極阻尼力矩系數導數m′xz已經利用轉速數據對其進行了辨識，赤道阻尼力矩系數導數m′zz、馬氏力系數導數的導數c″z、馬氏力矩系數導數的導數m″y通過敏感性分析以及相關試驗數據的分析認為可辨識性很小，在參數辨識過程中，利用風洞吹風數據代替其真實值，所以在本文中將m′xz、m′zz、c″z、m″y看作是已知量。

1 系統模型

為了更加精確地描述彈丸在空中的實際運動情況，并考慮到高速旋轉彈丸的縱橫向運動相互耦合，本文采用6自由度彈道模型作為參數辨識的理論模型。其具體表達式為

(1)

(2)

(3)

(4)

(5)

(6)

(7)

(8)

(9)

(10)

(11)

(12)

式中：v為速度；vr為相對速度；vrε、vrη、vrξ分別為相對速度在彈軸坐標系下的分量；x、y、z為坐標數據；θv為速度高低角；ψv為速度方向角，ωε、ωη、ωξ為3個方向轉速分量；θa為彈軸高低角；ψa為彈軸方向角；m為彈丸質量；S為參考面積；l為參考彈長；d為參考彈徑；g為重力加速度；ρ為空氣密度；δud、δlr分別為高低攻角和方向攻角；δr為相對攻角；δM、δN分別為附加力矩的氣動偏心角和附加力的氣動偏心角；wxv、wyv、wzv為風速在速度坐標系Ovxvyvzv下的分量；Λ為設計點緯度；ΩE為地球自轉角速率；αd為射向；cx為阻力系數；cy為升力系數；cz為馬氏力系數；mz為翻轉力矩系數；m′y為馬氏力矩系數導數；γ1為氣動偏心角所在平面與彈軸坐標系η軸之間的夾角；γ2為附加力矩偏心角所在平面與彈軸坐標系η軸之間的夾角。

2 PSO算法原理及其改進

2.1 PSO算法

PSO算法是由Kennedy等[8]提出的基于群智能的優化算法，它采用位置- 速度的數學模型，每個粒子代表解空間中的一個候選解，通過適應度函數來判斷候選解的優劣，不同的工程應用對應于不同的適應度函數，通過更新粒子的數值，最終得到最優解。其基本的步驟[10]可描述為：

1) 根據實際情況確定相關的參數。主要有：種群大小N，粒子的個數及維數D，慣性權重w，個體認知系數cI，社會認知系數cS，最大迭代代數kmax，粒子位置和速度的取值范圍等。

2) 初始化粒子群。隨機產生N個D維初始粒子Xi=(xi1xi2…xiD)T，i=1,…,N，并隨機產生每個粒子所對應的粒子速度vi=(vi1vi2…viD)T.

3) 位置和速度更新。速度和位置更新的表達式如(13)式和(14)式所示：

(13)

(14)

(15)

f(·)表示實際問題中的目標函數。

(16)

5)判斷是否滿足停止條件，若滿足，停止計算輸出最優解，若不滿足轉到步驟3，繼續迭代。

2.2 自適應混沌變異PSO算法

標準PSO算法廣泛應用與函數優化、參數辨識、控制系統設計等領域[11-13]。該算法具有進化運算簡單、容易實現等優點。但也存在一些問題，針對這些問題，文中對標準PSO算法進行了改進得到自適應混沌變異PSO算法(ACM-PSO)。

2.2.1 自適應慣性權重

(13)式中的慣性權重w的取值對算法的尋優性能有著重要的影響，算法在迭代初期應該增大慣性權重以快速逼近最優解，而在迭代的后期應該減小慣性權重來減小搜索范圍，本文提出一種新的自適應更新策略對慣性權重進行實時更新，其表達式為

(17)

2.2.2 混沌初始化策略

在PSO算法中，種群的初始化通常是隨機產生，初始位置的散布程度及其在搜索空間中的位置是否均勻，將直接影響整個搜索過程的收斂速度和算法的尋優效率[14]。當解空間較大時，粒子數目不夠多時，有可能造成初始化的粒子不能夠均勻分布。混沌是一種無規則的運動狀態，具有非常強烈的非線性特點，且其運動狀態具有很好的遍歷性和隨機性。本文利用混沌優化策略對粒子群進行初始化，立方映射具有較均勻的遍歷性分布區間，而且只要初始值不為0，就會產生分布均勻的混沌序列，有效地保證了種群在解空間中的均勻分布，從而提高算法的搜索效率。本文選取立方映射的方法產生初始混沌序列。該映射的表達式為

(18)

圖3為隨機法和混沌序列法產生范圍為0～5之間的初始粒子，從圖中可以明顯地看出：混沌序列法產生的粒子分布更均勻，更廣泛，充滿了整個解空間，粒子分布的越均勻，越容易快速找到最優解。

圖3 混沌初始化和隨機初始化粒子比較Fig.3 Comparison of chaos and random initialized particle generating methods

2.2.3 早熟判別機制

(19)

(20)

式中：μ為高斯白噪聲。

2.2.4 ACM-PSO算法流程

1) 確定相關參數。主要有：種群大小N，粒子的個數及維數D，慣性權重最大值wmax和最小值wmin，個體認知系數cI，社會認知系數cS，最大迭代代數kmax，粒子位置和速度的取值范圍等。

2) 初始化粒子。根據2.2.2節的方法產生混沌初始化粒子。

3) 更新慣性權重。根據(17)式來更新慣性權重。

4) 速度位置更新。根據(13)式、(14)式進行速度和位置更新。

5) 判斷是否早熟。根據2.2.3節判斷是否進入早熟，若進入早熟按照式(20)式以概率Pm對群體最優位置進行變異操作。若沒有，轉入步驟6.

6) 判斷是否滿足停止條件，若不滿足轉入步驟3，若滿足，停止計算，輸出最優結果。

3 算法測試

為了驗證文中ACM-PSO算法的有效性，利用標準測試函數對該算法進行了測試，并與標準PSO算法進行了對比。

選取3個基本測試函數分別為Sphere函數fSp，Rastrigrin函數fR和Schaffer函數fSc.

Sphere函數：

(21)

Rastrigrin函數：

(22)

Schaffer函數：

(23)

算法驗證實驗中相關參數設置為：最大迭代代數1 000，最大速度vmax=0.3，最小速度vmin=-0.3，最大慣性權重wmax=0.9，最小慣性權重wmin=0.5，個體認知系數cI=1.496 18，社會認知系數cS=1.496 18，方差判斷閾值ε=0.05. 分別利用PSO算法和ACM-PSO算法對上述函數進行測試，分別取粒子數為50、100和200對函數fSp、fR和fSc進行100次獨立實驗，其計算結果如表1所示。

從表1可以看出：在不同的測試函數下，粒子數越多，求解精度越高。在相同粒子數的情況下，ACM-PSO算法比PSO算法具有明顯的優越性，ACM-PSO的計算結果更加接近真實值，測試的最小值、平均值以及方差均比PSO算法要小。圖4是Sphere函數的目標函數值隨進化代數的變化關系。從圖4中還可以看出，ACM-PSO在第300代進入收斂狀態，而PSO大約在第500代進入收斂狀態，前者的收斂速度要比后者快得多。綜上所述：本文使用的ACM-PSO算法可以有效地解決函數的最優化問題，且有較高的求解精度，可以滿足實際工程需要。

4 彈丸氣動參數辨識

利用ACM-PSO算法處理氣動參數辨識問題，本文辨識的主要氣動參數有：零升阻力系數cx0、誘導阻力系數cxi、升力系數導數c′y和翻轉力矩系數導數m′z. 在本文氣動參數的辨識過程中還需要其他的氣動參數，其中極阻尼力矩系數m′xz已經通過彈丸的實測轉速數據對其進行了辨識和處理，由于篇幅有限，本文中不再詳述。馬氏力系數導數的導數、馬氏力矩系數導數的導數和赤道阻尼力矩系數導數均采用風洞吹風數據，在本文中認為是已知量。參數辨識的本質是尋找一組待辨識參數，使得目標函數最小。

4.1 目標函數

由系統辨識理論可知，參數的最大似然估計是漸進無偏、漸進一致和漸進有效的；辨識實踐也證明，最大似然法是有效、實用的參數辨識方法。本文中待辨識參數向量為

θ=[m′zcx0cxic′y]Τ.

(24)

實際測量值矢量為

ym(i)=[δr(ti)v(ti)x(ti)y(ti)z(ti)]Τ.

(25)

相對應的由動力學方程計算得到的觀測量矢量為y(i)，則計算值和測量值之間的誤差為

e(i)=y(i)-ym(i).

(26)

表1 不同測試函數的目標函數值Tab.1 Objective values of different test functions

圖4 目標函數值隨進化代數的變化關系Fig.4 Relationship between objective value and generation number

最大似然判據，即目標函數可表示為

(27)

式中：M為數據個數；R為觀測噪聲的協方差矩陣，當觀測噪聲的統計特性未知時，采用R的最優估計值：

(28)

4.2 仿真校驗

以某高速旋轉彈丸為背景，利用6自由度彈道方程生成仿真的飛行試驗數據，再用仿真數據進行氣動參數辨識，并與仿真所用的氣動參數值進行比較，從而驗證算在解決旋轉彈丸氣動參數辨識問題方面的有效性。在生成仿真數據的過程中，氣象諸元采用某次實際試驗的氣象條件帶入方程進行計算，能夠更好的模擬實際彈道數據。圖5為實際測量的風速。

圖5 實際風速Fig.5 Actual air speed

生成飛行試驗數據，所需要的主要初始數據為：初速930 m/s，射角65°；PSO算法和ACM-PSO算法的相關參數的設置同第3節。仿真結果如表2所示。

工程上通常將辨識出的氣動參數帶入彈道方程組(1)式～(12)式計算出各個觀測物理量的值，通過和原始觀測值相比較來判斷參數辨識的準確與否以及其誤差大小。本文同樣采用這種方法進行反驗算。 圖6～圖10分別給出了攻角、速度、射程、高程和側偏測量值和用辨識出的氣動參數進行計算的計算值之間的關系，可以看出：由辨識出的參數進行計算的結果和測量值之間差距非常小，測量值落點坐標(21 622.21 m, 7 365.19 m, 2 055.59 m)，計算值落點坐標(21 612.32 m, 7 371.39 m, 2 061.71 m)，兩者之間的差值為13.18 m，對于榴彈來講，該辨識結果已經具有了很高的精度。

表2 主要氣動參數辨識結果Tab.2 Identified result of aerodynamic parameters

注：工況1中，仿真值作為測量值；工況2中，仿真值加上高斯白噪聲作為測量值，信噪比為30 dB；工況3中，仿真值加上高斯白噪聲作為測量值，信噪比為20 dB；工況4中，仿真值加上高斯白噪聲作為測量值，信噪比為10 dB.

圖6 攻角對比圖Fig.6 Comparison of angles of attack

圖7 速度對比圖Fig.7 Comparison of velocities

圖8 射程對比圖Fig.8 Comparison of firing ranges

圖9 高程對比圖Fig.9 Comparison of ballistic heights

圖10 側偏對比圖Fig.10 Comparison of ballistic side-range

5 結論

對于旋轉彈丸的氣動參數辨識問題，利用PSO算法對其進行參數辨識，在公開報道的資料中是很少見的，將最大似然準則應用到該類型的彈丸參數辨識的工作中，亦是很少見的，本文將這兩種方法綜合起來應用到該類彈丸的氣動參數辨識中，并對基本PSO算法進行了改進，取得了很好的辨識結果并得出相關結論：

1) 提出的ACM-PSO算法比基本PSO算法收斂速度更快、精度更高。

2) 通過該方法得到高速旋轉彈丸的氣動參數精度可以滿足工程需要。

3) 該方法可以輔助其他方法進行氣動參數辨識工作，提高辨識的魯棒性。

References)

[1] 閆章更, 祁載康. 射表技術[M]. 北京:國防工業出版社, 2000. YAN Zhang-geng, QI Zai-kang. Firing table technology[M]. Beijing:National Defense Industry Press, 2000.(in Chinese)

[2] Chase A T, Mcdonald R A. Flight testing small UAVs for aerodynamic parameter estimation[C]∥AIAA Atompspheric Flight Mechanics Conference. National Harbor, MD, US:AIAA, 2014.

[3] Peyada N K, Sen A, Ghosh A K. Aerodynamic characterization of HANSA-3 aircraft using equation error, maximum likelihood and filter error methods[C]∥International Multi Conference on Engineers and Computer Scientists. Hong Kong, China: IMECS, 2008.

[4] 趙亞男, 錢杏芳. 反坦克導彈氣動參數辨識的最大似然法[J]. 兵工學報, 1991, 12(4):13-21. ZHAO Ya-nan, QIAN Xing-fang. Aerodynamic identification for antitank missile using maximum likelihood estimator[J]. Acta Armamentarii, 1991, 12(4):13-21.(in Chinese)

[5] 王永驥, 劉莎, 劉磊, 等. 基于粒子群優化算法的氣動參數在線辨識方法[J]. 華中科技大學學報:自然科學版, 2016, 44(3):116-120. WANG Yong-ji, LIU Sha, LIU Lei, et al. Aerodynamic parameters online identification method based on PSO[J]. Journal of Huazhong University of Science & Technology:Natual Science Edition, 2016, 44(3):116-120. (in Chinese)

[6] Chapman G T, Kirk D B. A new method for extracting aerodynamic coefficients from free-flight data[J]. AIAA Journal, 1970, 8(4):753-758.

[7] 張天姣, 汪清, 何開鋒. 粒子群算法在氣動力參數辨識中的應用[J]. 空氣動力學學報, 2010, 28(6): 633-638. ZHANG Tian-jiao, WANG Qing, HE Kai-feng. Application of particle swarm optimization for aerodynamic parameter estimation[J]. Acta Aerodynamic Sinica, 2010, 28(6):633-638.(in Chinese)

[8] Kennedy J, Eberhart R. Particle swarm optimization[C]∥1995 IEEE International Conference on Neural Networks. Perth, WA, US:IEEE, 1995.

[9] 韓子鵬. 彈箭外彈道學[M]. 北京:北京理工大學出版社, 2014. HAN Zi-peng. Exterior ballistics of projectiles and rockets[M]. Beijing:Beijing Institute of Technology Press, 2014.(in Chinese)

[10] 崔志華, 曾建潮. 微粒群優化算法[M]. 北京: 科學出版社, 2011. CUI Zhi-hua, ZENG Jian-chao. Particle swarm optimization[M]. Beijing:Science Press, 2011.(in Chinese)

[11] 李軍偉, 程詠梅, 陳克喆, 等. 基于AIWCPSO算法的三次樣條氣動參數插值方法[J]. 控制與決策. 2014(1):129-134. LI Jun-wei, CHENG Yong-mei, CHEN Ke-zhe, et al. Cubic spline interpolation method of aerodynamic parameters based on AIWCPSO algorithm[J]. Control and Decision, 2014(1):129-134.(in Chinese)

[12] 劉奇治. 高超聲速飛行器氣動參數辨識與LPV混合控制方法研究[D]. 天津:天津大學, 2014. LIU Qi-zhi. Research on aerodynamic parameter identification and LPV mixing control for hypersonic vehicle[D]. Tianjin:Tianjin University, 2014.(in Chinese)

[13] Boundjou O D, Jacobsen L. Adaptive particle swarm optimization applied to aircraft control[C]∥20th AIAA International Space Plances and Hypersonic Systems and Technologies Conference. Glasgow, Scotland:AIAA, 2015.

[14] 程澤, 董夢男, 楊添剴, 等. 基于自適應混沌粒子群算法的光伏電池模型參數辨識[J]. 電工技術學報, 2014, 29(9):245-252. CHENG Ze, DONG Meng-nan, YANG Tian-kai, et al. Extraction of solar cell model parameters based on self-adaptive chaos particle swarm optimization alogorithm[J]. Tramsactions of China Electro-Technical Society, 2014, 29(9):245-252. (in Chinese)

Identification of Spinning Projectile Aerodynamic Parameters Using Adaptive Chaotic Mutation Particle Swarm Optimization

GUAN Jun1，ZHOU Jia-shen2, YI Wen-jun1，LIU Shi-ping3，CHANG Si-jiang3，SHI Ji-gang1

(1.National Key Laboratory of Transient Physics, Nanjing University of Science and Technology, Nanjing 210094, Jiangsu, China;2.Military Representation Office of Navy Ammunition in Shenyang Area, Shenyang 110045, Liaoning, China;3.School of Energy and Power Engineering, Nanjing University of Science and Technology, Nanjing 210094, Jiangsu, China)

The maximum likelihood estimation is applied to the identification of spinning projectile aerodynamic parameters. A new algorithm called adaptive chaotic mutation particle swarm optimization is proposed to solve the optimal solution of aerodynamic parameters, thus obtaining the aerodynamic parameters of a spinning projectile. The proposed algorithm is to use an adaptive weight function, generate the initial particles based on chaos theory, and set a discriminant mechanism which judges whether the algorithm falls into the local optimum. If the algorithm falls into the local convergence, the mutation operator is used to make the algorithm jump out of local. The common test function is used to test this algorithm. The test result shows that the proposed algorithm has the advantages of more quick convergence, higher optimization precision and wide range of application compared to basic PSO. Simulated ballistic data is used to test the algorithm. The result shows that the proposed algorithm can identify the aerodynamic parameters effectively with high precision and quickly converging velocity.

ordnance science and technology; projectile; aerodynamic parameters identification; particle swarm optimization algorithm; maximum likelihood estimation

2016-05-03

國家自然科學基金項目(11472136、11402117)

管軍(1987—)，男，博士研究生。E-mail:guanjun8710@163.com

易文俊(1970—)，男，教授，博士生導師。E-mail:yiwenjun0@163.com

TJ012.3+4

A

1000-1093(2017)01-0073-08

10.3969/j.issn.1000-1093.2017.01.010