用于戰(zhàn)術(shù)導(dǎo)彈氣動(dòng)參數(shù)辨識(shí)的輸入設(shè)計(jì)研究*

楊闖，張弫

(北京電子工程總體研究所，北京 100854)

用于戰(zhàn)術(shù)導(dǎo)彈氣動(dòng)參數(shù)辨識(shí)的輸入設(shè)計(jì)研究*

楊闖，張弫

(北京電子工程總體研究所，北京 100854)

為了研究h=8 km，Ma=2條件附近導(dǎo)彈飛行試驗(yàn)中輸入設(shè)計(jì)對(duì)氣動(dòng)參數(shù)辨識(shí)精度的影響，驗(yàn)證氣動(dòng)力模型的有效性，以方波信號(hào)為基礎(chǔ)設(shè)計(jì)了3種開環(huán)舵偏指令。通過開環(huán)飛行彈道仿真試驗(yàn)，采用遞推最小二乘法，對(duì)多項(xiàng)式非線性氣動(dòng)力模型中的氣動(dòng)參數(shù)進(jìn)行了辨識(shí)。氣動(dòng)參數(shù)辨識(shí)結(jié)果表明，單級(jí)方波，偶極方波，“211”多級(jí)方波3種輸入用于氣動(dòng)參數(shù)辨識(shí)是可行的，辨識(shí)精度是可接受的，辨識(shí)所用的縱向三自由度氣動(dòng)力模型是有效的。

戰(zhàn)術(shù)導(dǎo)彈；氣動(dòng)參數(shù)辨識(shí)；輸入設(shè)計(jì)；遞推最小二乘；氣動(dòng)模型；飛行試驗(yàn)

0 引言

采用飛行試驗(yàn)數(shù)據(jù)進(jìn)行氣動(dòng)參數(shù)辨識(shí)是驗(yàn)證導(dǎo)彈氣動(dòng)特性的一種有效方法。與固定翼飛機(jī)不同，對(duì)于戰(zhàn)術(shù)導(dǎo)彈而言，其飛行試驗(yàn)成本較高，無法進(jìn)行大批次的飛行試驗(yàn)，且每次飛行試驗(yàn)時(shí)間短暫，用于氣動(dòng)參數(shù)辨識(shí)的時(shí)間更是十分有限[1]。為了能在有限時(shí)間內(nèi)充分激發(fā)導(dǎo)彈運(yùn)動(dòng)模態(tài)，使得飛行試驗(yàn)為導(dǎo)彈氣動(dòng)參數(shù)辨識(shí)提供所需信息量，提高氣動(dòng)參數(shù)辨識(shí)的精度，對(duì)辨識(shí)時(shí)間段內(nèi)的輸入信號(hào)進(jìn)行精心設(shè)計(jì)是十分必要的[2-4]。

文獻(xiàn)[5-7]的研究表明以方波信號(hào)為基礎(chǔ)的輸入設(shè)計(jì)對(duì)于有人駕駛飛機(jī)，閉環(huán)控制無人機(jī)，開環(huán)控制飛行器的氣動(dòng)參數(shù)辨識(shí)是有效的。戰(zhàn)術(shù)導(dǎo)彈飛行試驗(yàn)中氣動(dòng)參數(shù)隨著試驗(yàn)馬赫數(shù)，高度，攻角等眾多因素劇烈變化，為了提高氣動(dòng)參數(shù)辨識(shí)精度，在參考文獻(xiàn)[8]中導(dǎo)彈縱向三自由度線性氣動(dòng)模型基礎(chǔ)上，結(jié)合文獻(xiàn)[9-13]中非線性氣動(dòng)力模型的構(gòu)型，本文提出一種改進(jìn)的縱向三自由度多項(xiàng)式形式非線性氣動(dòng)模型。

針對(duì)該氣動(dòng)力模型，在Ma=2，h=8 km的條件附近設(shè)計(jì)3種方波形式的輸入指令，通過開環(huán)飛行彈道仿真試驗(yàn)獲取試驗(yàn)數(shù)據(jù)，采取遞推最小二乘法從仿真試驗(yàn)數(shù)據(jù)中辨識(shí)出相應(yīng)氣動(dòng)參數(shù)，對(duì)3種輸入的氣動(dòng)參數(shù)辨識(shí)效果進(jìn)行分析研究，為飛行試驗(yàn)中用于氣動(dòng)參數(shù)辨識(shí)的輸入設(shè)計(jì)提供一定的參考。

1 參數(shù)辨識(shí)最小二乘遞推算法

戰(zhàn)術(shù)導(dǎo)彈動(dòng)力學(xué)系統(tǒng)參數(shù)辨識(shí)問題的一般性描述為[9]

(1)

式中：x(t)為n維狀態(tài)向量；y(t)為m維輸出向量；z(t)為m維觀測(cè)向量；u(t)為l維輸入向量；θ為p維參數(shù)向量；η(t)為q維隨機(jī)噪聲向量；Γ為n×q維系統(tǒng)噪聲分布矩陣；F為已知實(shí)值函數(shù)；H為觀測(cè)矩陣。

參數(shù)的遞推估計(jì)，每取得一次新觀測(cè)數(shù)據(jù)后，在前一次估計(jì)結(jié)果基礎(chǔ)上，采用新引入觀測(cè)數(shù)據(jù)對(duì)前次估計(jì)結(jié)果進(jìn)行修正，從而遞推得出新的參數(shù)估計(jì)值[14-15]。這樣，隨著新觀測(cè)數(shù)據(jù)的不斷引入，參數(shù)估計(jì)值不斷更新，直到估計(jì)值達(dá)到滿意的精度為止。其基本思想[15]可概括為

(2)

最小二乘參數(shù)估計(jì)的遞推算法具體算式[15]如下：

(3)

(4)

式中：α為充分大的實(shí)數(shù)；ε為充分小的實(shí)向量。按上述算法進(jìn)行遞推計(jì)算，直到辨識(shí)參數(shù)向量滿足以下收斂條件[15]：

(5)

2 縱向氣動(dòng)參數(shù)辨識(shí)模型

早期飛行器外作用力模型一般采用線性模型，隨著飛行器氣動(dòng)性能的不斷提升，為了更精確描述外作用力，現(xiàn)代飛行器尤其是戰(zhàn)術(shù)導(dǎo)彈一般多采用非線性氣動(dòng)力模型。

常用非線性模型包括多項(xiàng)式模型，樣條函數(shù)模型，階躍過渡函數(shù)模型等[9]。其中多項(xiàng)式模型是最簡(jiǎn)單方便的一種模型，由于其各項(xiàng)物理意義明確，模型有較高的光滑性，便于進(jìn)行參數(shù)辨識(shí)，因此在工程實(shí)踐中被廣泛采用[9-10]。

2.1氣動(dòng)力及氣動(dòng)力矩模型

在對(duì)飛行器常用非線性氣動(dòng)力模型進(jìn)行了分析研究后，本文針對(duì)戰(zhàn)術(shù)導(dǎo)彈，提出了一種改進(jìn)的縱向三自由度多項(xiàng)式形式的非線性氣動(dòng)力模型，具體表達(dá)式為[9-13]

(6)

2.2狀態(tài)方程組

在彈體坐標(biāo)系中，導(dǎo)彈縱向三自由度狀態(tài)方程組為

(7)

式中：vx，vy分別為導(dǎo)彈軸向速度和法向速度；ωz為導(dǎo)彈俯仰角速度；?為導(dǎo)彈俯仰角；q為動(dòng)壓；S為參考面積；Jz為導(dǎo)彈俯仰轉(zhuǎn)動(dòng)慣量。

補(bǔ)充方程組：

(8)

式中：Nx，Ny分別為軸向過載和法向過載；?為俯仰角；α為攻角；G為重力；R為發(fā)動(dòng)機(jī)推力；φR，ψR(shí)為發(fā)動(dòng)機(jī)推力偏心角。

2.3觀測(cè)方程組

取導(dǎo)彈被動(dòng)段俯仰角速度，軸向過載，法向過載，俯仰角，攻角作為氣動(dòng)參數(shù)辨識(shí)的觀測(cè)量。觀測(cè)方程組為

(9)

式中：vi(i=1,2,…,5)為互不相關(guān)的零均值高斯分布隨機(jī)測(cè)量噪聲。

2.4待辨識(shí)參數(shù)

待估計(jì)參數(shù)為式(6)氣動(dòng)模型中的未知參數(shù)，向量形式表達(dá)式為

(10)

3 仿真試驗(yàn)設(shè)計(jì)

方波式輸入是導(dǎo)彈操縱機(jī)構(gòu)一種典型的偏轉(zhuǎn)方式，由于波形簡(jiǎn)單，工程上容易實(shí)現(xiàn)，且輸入后導(dǎo)彈響應(yīng)強(qiáng)烈，引起過渡過程中超調(diào)量較大，因此在用于氣動(dòng)參數(shù)辨識(shí)的飛行試驗(yàn)中被廣泛采用[16-18]。

仿真試驗(yàn)共設(shè)計(jì)3條彈道，第1條彈道在t=30 s時(shí)斷開穩(wěn)定控制系統(tǒng)反饋回路，在辨識(shí)窗口內(nèi)加載如圖1a)所示的單極方波開環(huán)舵偏指令，圖1b)為對(duì)應(yīng)的攻角響應(yīng)曲線。所加載的具體指令為

(11)

圖1 單極方波輸入及其攻角響應(yīng)Fig.1 Single square wave input and the response of attack angle

第2條彈道在t=31 s時(shí)斷開穩(wěn)定控制系統(tǒng)反饋回路，在辨識(shí)窗口內(nèi)加載如圖2a)所示的偶極方波開環(huán)舵偏指令，攻角響應(yīng)曲線如圖2b)所示。加載具體指令為

(12)

圖2 偶極方波輸入及其攻角響應(yīng)Fig.2 Double square wave input and the response of attack angle

第3條彈道在t=31 s時(shí)斷開穩(wěn)定控制系統(tǒng)反饋回路，在辨識(shí)窗口內(nèi)加載如圖3a)所示的“211”多級(jí)方波開環(huán)舵偏指令，對(duì)應(yīng)的攻角響應(yīng)曲線如圖3b)所示。具體指令為

(13)

圖3 “211”多級(jí)方波輸入及其攻角響應(yīng)Fig.3 “211” multistep square wave input and the response of attack angle

在仿真試驗(yàn)中，3條彈道理論上的辨識(shí)窗口為(8 km,2Ma)附近。由于加載的3種開環(huán)舵偏指令不同，在理論辨識(shí)窗口附近彈道會(huì)產(chǎn)生偏離，因此3條彈道實(shí)際的辨識(shí)窗口均取在(8±0.3 km,2±0.3Ma)的區(qū)間內(nèi)，辨識(shí)窗口持續(xù)時(shí)間均為6 s，所加載的開環(huán)舵偏指令寬度均為2 s，幅值均為7°。圖4給出了辨識(shí)窗口附近馬赫數(shù)及高度的變化情況。表1為3條仿真彈道辨識(shí)窗口中馬赫數(shù)與海拔高度的變化范圍。

表1 辨識(shí)窗口對(duì)應(yīng)馬赫數(shù)與高度區(qū)間Table 1 Mach range and altitude range in identification window

圖4 辨識(shí)窗口附近馬赫數(shù)及高度變化Fig.4 Variation of mach and altitude near identification window

4 辨識(shí)結(jié)果分析

根據(jù)彈道仿真試驗(yàn)數(shù)據(jù)，采用遞推最小二乘法對(duì)導(dǎo)彈縱向各氣動(dòng)參數(shù)進(jìn)行辨識(shí)，將各參數(shù)辨識(shí)結(jié)果匯總后分類進(jìn)行比較分析。

4.1軸向力系數(shù)辨識(shí)結(jié)果

3種輸入對(duì)軸向力系數(shù)辨識(shí)結(jié)果分別如圖5～7所示。

圖5 單極方波輸入軸向力系數(shù)辨識(shí)結(jié)果Fig.5 Identification results of axial force coefficient using single square wave input

圖6 偶極方波輸入軸向力系數(shù)辨識(shí)結(jié)果Fig.6 Identification results of axial force coefficient using double square wave input

圖7 “211”多級(jí)方波輸入軸向力系數(shù)辨識(shí)結(jié)果Fig.7 Identification results of axial force coefficient using “211”multistep square wave input

3種輸入對(duì)軸向力系數(shù)辨識(shí)精度的定量分析由表2給出。在辨識(shí)窗口內(nèi)取120個(gè)辨識(shí)點(diǎn)計(jì)算其相對(duì)誤差，經(jīng)過統(tǒng)計(jì)后發(fā)現(xiàn)，對(duì)于軸向力系數(shù)，單極方波輸入的辨識(shí)精度較高，平均相對(duì)誤差為0.664%，偶極方波與“211”多級(jí)方波辨識(shí)相對(duì)誤差次之。

表2 軸向力系數(shù)辨識(shí)相對(duì)誤差Table 2 Relative error of axial force coefficient identification

4.2法向力系數(shù)辨識(shí)結(jié)果

圖8～10分別為3種輸入在辨識(shí)窗口內(nèi)對(duì)法向力系數(shù)的辨識(shí)結(jié)果。

圖8 單極方波輸入法向力系數(shù)辨識(shí)結(jié)果Fig.8 Identification results of normal force coefficient using single square wave input

圖9 偶極方波輸入法向力系數(shù)辨識(shí)結(jié)果Fig.9 Identification results of normal force coefficient using double square wave input

圖10 ”211”方波輸入法向力系數(shù)辨識(shí)結(jié)果Fig.10 Identification results of normal force coefficient using “211” multistep square wave input

表3為法向力系數(shù)辨識(shí)的相對(duì)誤差。在辨識(shí)窗口內(nèi)，對(duì)法向力系數(shù)辨識(shí)效果較好的是偶極方波輸入，平均相對(duì)誤差為5.774%，單極方波和“211”多級(jí)方波對(duì)法向力系數(shù)辨識(shí)的平均相對(duì)誤差分別為6.110%和6.255%。

表3 法向力系數(shù)辨識(shí)相對(duì)誤差Table 3 Relative error of normal force coefficient identification

4.3俯仰力矩系數(shù)辨識(shí)結(jié)果

圖11～13分別為單極方波，偶極方波，“211”多級(jí)方波在辨識(shí)窗口內(nèi)對(duì)俯仰力矩系數(shù)的辨識(shí)結(jié)果。

圖11 單極方波俯仰力矩系數(shù)辨識(shí)結(jié)果Fig.11 Identification results of pitching moment coefficient using single square wave input

圖12 偶極方波俯仰力矩系數(shù)辨識(shí)結(jié)果Fig.12 Identification results of pitching moment coefficient using double square wave input

圖13 “211”多級(jí)方波俯仰力矩系數(shù)辨識(shí)結(jié)果Fig.13 Identification results of pitching moment coefficient using “211”multistep square wave input

表4給出了各輸入信號(hào)對(duì)俯仰力矩系數(shù)辨識(shí)的相對(duì)誤差。在持續(xù)時(shí)間6 s的辨識(shí)窗口中，3種輸入辨識(shí)結(jié)果的平均相對(duì)誤差均超過了8%。俯仰力矩系數(shù)總體辨識(shí)精度略遜于軸向力系數(shù)和法向力系數(shù)，這說明俯仰力矩系數(shù)的氣動(dòng)力模型有一定的改進(jìn)潛力。

表4 俯仰力矩系數(shù)辨識(shí)相對(duì)誤差Table 4 Relative error of pitching moment coefficient identification

5 結(jié)束語

本文設(shè)計(jì)了3種輸入指令，在(8 km,2Ma)條件附近，通過開環(huán)飛行彈道仿真試驗(yàn)獲取了相關(guān)試驗(yàn)數(shù)據(jù)，針對(duì)一種改進(jìn)的縱向三自由度非線性多項(xiàng)式氣動(dòng)力模型，采用遞推最小二乘法對(duì)氣動(dòng)參數(shù)進(jìn)行了辨識(shí)，比較了3種不同輸入的氣動(dòng)參數(shù)辨識(shí)精度。

辨識(shí)結(jié)果表明，對(duì)于軸向力系數(shù)，單極方波指令辨識(shí)效果較好。對(duì)于法向力系數(shù)和俯仰力矩系數(shù)，3種輸入的辨識(shí)效果接近。在開環(huán)飛行試驗(yàn)中，單級(jí)方波，偶極方波，“211”多級(jí)方波3種輸入用于氣動(dòng)參數(shù)辨識(shí)是可行的。

辨識(shí)結(jié)果驗(yàn)證了縱向三自由度氣動(dòng)力模型的有效性，參數(shù)辨識(shí)的相對(duì)誤差表明法向力系數(shù)模型和俯仰力矩系數(shù)模型仍有改進(jìn)的空間。

[2] 蔡金獅．飛行器氣動(dòng)參數(shù)辨識(shí)進(jìn)展[J]．力學(xué)進(jìn)展，1987，17(4)：467-478．

CAI Jin-shi.Advances in Identification of Aircraft Aerodynamic Parameters[J].Advances in Mechanics，1987，17(4)：467-478.

[3] 于君彩．飛行數(shù)據(jù)處理與氣動(dòng)參數(shù)辨識(shí)[D]．西安：西北工業(yè)大學(xué)，2007．

YU Jun-cai.Flight Data Processing and Aerodynamic Parameter Identification[D].Xi’an :Northwestern Polytechnical University,2007.

[4] 張?zhí)戽]環(huán)控制戰(zhàn)術(shù)導(dǎo)彈氣動(dòng)參數(shù)辨識(shí)技術(shù)研究[D]．綿陽：中國空氣動(dòng)力研究與發(fā)展中心，2010．

ZHANG Tian-jiao.Aerodynamic Parameter Identification Technology of Closed-Loop Controlled Tactical Missiles[D].Mianyang：China Aerodynamics Research and Development Center，2010.

[5] SHAFER M F.Flight Investigation of Various Control Inputs intended for Parameter Estimation[C]∥The 11th Atmospheric Flight Mechanics Conference，Seattle，USA，August，1984：21-23．

[6] KRINGS M，HENNING K，THIELECKE F.Flight Test Oriented Autopilot Design for Improved Aerodynamic Parameter Identification[C]∥The 2nd CEAS Specialist Conference on Guidance，Navigation and Control，Delft，Netherland，April，2013：11-12.

[7] ALBISSER M，DOBRE S，BERNER C，et al.Identifiability Investigation of the Aerodynamic Coefficients from Free Flight Tests[C]∥AIAA Atmospheric Flight Mechanics Conference，Boston，USA，August，2013：19-22.

[8] 錢杏芳，林瑞雄，趙亞男．導(dǎo)彈飛行力學(xué)[M]．北京：北京理工大學(xué)出版社，2011．

QIAN Xing-fang，LIN Rui-xiong，ZHAO Ya-nan.Missile Flight Dynamics[M].Beijing：Beijing Institute of Technology Press,2011．

[9] 蔡金獅．動(dòng)力學(xué)系統(tǒng)辨識(shí)與建模[M]．北京：國防工業(yè)出版社，1991．

CAI Jin-shi.Dynamic System Identification and Modeling[M].Beijing：National Defence Industry Press，1991.

[10] James E Kain，Charles M Brown，Jang G Lee．Missile Aerodynamic Parameter and Structure Identification from Flight Test Data[R].Florida：Air Force Armament Laboratory，1977：ADA056343．

[11] HALL W E，VINCENT J R.Flight Test Design for CH-47 Parameter Identification[R].NASA Contractor Report，1978:141-148．

[12] TLIFF K W.Parameter Estimation for Flight Vehicle[J].Guidance and Control，1989，12(5)：608-622．

[13] 王曉鵬，焦天峰，李響．戰(zhàn)術(shù)導(dǎo)彈氣動(dòng)參數(shù)辨識(shí)仿真[C]∥第14屆中國系統(tǒng)仿真技術(shù)及其應(yīng)用學(xué)術(shù)年會(huì)，三亞，2012:802-806.

WANG Xiao-peng，JIAO Tian-feng，LI Xiang.Simulation of Aerodynamic Parameter Identification for Tactical Missile[C]∥The 14th Chinese Conference on System Simulation Technology & Application，Sanya，2012:802-806.

[14] 強(qiáng)明輝，張京娥．基于MATLAB的遞推最小二乘法辨識(shí)與仿真[J]．自動(dòng)化與儀器儀表，2008(6)：4-5．

ZHANG Ming-hui，ZHANG Jing-e.Recursive Least Squares Identification and Simulation Based on MATLAB[J].Automation and Instrumentation，2008(6):4-5.

[15] 候媛彬，汪梅．系統(tǒng)辨識(shí)及其MATLAB仿真[M]．北京：科學(xué)出版社，2004．

HOU Yuan-bin，WANG Mei.System Identification and MATLAB Simulation[M].Beijing：Science Press，2004.

[16] 汪清，錢煒祺，何開鋒．導(dǎo)彈氣動(dòng)參數(shù)辨識(shí)與優(yōu)化輸入設(shè)計(jì)[J]．宇航學(xué)報(bào)，2008，29(3)：789-793．

WANG Qing，QIAN Wei-qi，HE Kai-feng.AerodynamicParameter Identification and Optimal Input Design for Missile[J].Journal of Astronautics，2008，29(3):789-793.

[17] 徐秀峰．高超聲速飛行器閉環(huán)氣動(dòng)參數(shù)辨識(shí)[D]．哈爾濱：哈爾濱工業(yè)大學(xué)，2014．

XU Xiu-feng.Closed-Loop Identification of Aerodynamic Parameters for Hypersonic Vehicle[D].Harbin：Harbin Institute of Technology，2014.

[18] 李超．高超聲速飛行器氣動(dòng)參數(shù)辨識(shí)和輸入設(shè)計(jì)技術(shù)研究[D]．廈門：廈門大學(xué)，2014．

LI Chao.Research on Aerodynamic Parameter Identification and Input Design for Hypersonic Vehicle[D].Xiamen：Xiamen University，2014.

InputDesignfortheAerodynamicParametersIdentificationofTacticalMissile

YANG Chuang，ZHANG Zhen

(Beijing Institute of Electronic System Engineering，Beijing 100854，China)

In order to study the influence of input design on accuracy of aerodynamic parameters identification whenh=8 km,Ma=2 in flight test, three open-loop commands of control surface deflection based on square wave signal are designed. In the open-loop flight trajectory simulation test, recursive least squares method is used to estimate the parameters of nonlinear aerodynamic model. The aerodynamic parameters identification results show that it is feasible to use three inputs including single square wave, double square wave, “211” multistep square wave. The identification accuracy is acceptable. The longitudinal aerodynamic model of three degrees of freedom used for identification is effective.

tactical missile; aerodynamic parameter identification; input design; recursive least squares; aerodynamic model; flight test

2017-01-03；

2017-02-06

有

楊闖(1992-)，男，陜西咸陽人。 碩士生，主要研究方向?yàn)轱w行器總體設(shè)計(jì)。

通信地址：100854 北京142信箱30分箱E-mail：505511790@qq.com

10.3969/j.issn.1009-086x.2017.05.007

TJ761.1;V417+.1

A

1009-086X(2017)-05-0035-07