劉昊東，張慶振，郭云鶴，茅佳雯

(1. 北京航空航天大學(xué)自動(dòng)化科學(xué)與電氣工程學(xué)院，北京 100083； 2. 上海機(jī)電工程研究所，上海 201109)

0 引 言

高馬赫數(shù)飛行器是指飛行速度超過5Ma的有翼或無翼飛行器，如各類導(dǎo)彈、飛機(jī)、炮彈等[1]，具有巨大的軍事價(jià)值和經(jīng)濟(jì)價(jià)值，是21世紀(jì)航空航天領(lǐng)域的一個(gè)重要研究方向[2-3]。相比于傳統(tǒng)的飛行器，高馬赫數(shù)飛行器具有以下4個(gè)明顯的優(yōu)勢(shì)：飛行速度快；突防成功率高，攔截困難；射程遠(yuǎn)，殺傷力大；稍縱即逝，探測(cè)難度大[4-5]?；谝陨蟽?yōu)勢(shì)，各國都相繼對(duì)高馬赫數(shù)飛行器展開了深入的研究，并且取得了一定的研究成果[6]。變體飛行器是指一類能夠依據(jù)具體的飛行環(huán)境和任務(wù)要求實(shí)時(shí)調(diào)整外形結(jié)構(gòu)的飛行器[7]。變體飛行器將外形參數(shù)作為可控變量，利用外形參數(shù)對(duì)氣動(dòng)特性的影響來改變飛行器的性能，使其能夠適應(yīng)更寬的飛行空域和速域，從而能夠適應(yīng)更復(fù)雜的飛行任務(wù)和環(huán)境；同時(shí)，它在飛行過程中能夠根據(jù)飛行環(huán)境及任務(wù)實(shí)時(shí)調(diào)整外形，以獲得最優(yōu)的氣動(dòng)和操縱性能，從而達(dá)到減小能耗的目的。高馬赫數(shù)變體飛行器則是指具有變形能力的高馬赫數(shù)飛行器，憑借高馬赫數(shù)以及變形帶來的雙重優(yōu)點(diǎn)，大大提升了其適用范圍，在未來軍事或商業(yè)航天中將具有巨大的應(yīng)用潛力[8-9]。

目前高馬赫數(shù)變體飛行器大多采用機(jī)體/發(fā)動(dòng)機(jī)一體化設(shè)計(jì)技術(shù)，這種設(shè)計(jì)會(huì)導(dǎo)致高馬赫數(shù)變體飛行器的各個(gè)狀態(tài)量之間存在非線性、強(qiáng)耦合等特點(diǎn)。并且，相比于傳統(tǒng)飛行器，高馬赫數(shù)變體飛行器自身的氣動(dòng)特性更為復(fù)雜，這些因素就導(dǎo)致高馬赫數(shù)變體飛行器在制導(dǎo)、控制等方面存在著巨大的困難。同時(shí)，高馬赫數(shù)變體飛行器相比于傳統(tǒng)飛行器在飛行過程中馬赫數(shù)和飛行高度的變化范圍更大，環(huán)境因素變化情況更為復(fù)雜，這就會(huì)導(dǎo)致高馬赫數(shù)變體飛行器中某些模型參數(shù)發(fā)生較大變化[10-11]。因此，有必要采用在線辨識(shí)的方法，實(shí)時(shí)地獲取高馬赫數(shù)變體飛行器模型的各個(gè)參數(shù)，為飛行控制系統(tǒng)的在線更新、飛行能力的在線評(píng)估以及飛行器故障檢測(cè)等提供更加準(zhǔn)確的模型[12]。

系統(tǒng)辨識(shí)是利用系統(tǒng)的輸入輸出數(shù)據(jù)獲取系統(tǒng)模型的技術(shù)，主要包括試驗(yàn)設(shè)計(jì)、模型辨識(shí)、參數(shù)估計(jì)和系統(tǒng)驗(yàn)證4部分。目前系統(tǒng)辨識(shí)中常用的辨識(shí)算法有遞推最小二乘法(recursive least squares，RLS)[13]、遞推極大似然法[14]和卡爾曼濾波方法[15]等。

高馬赫數(shù)變體飛行器模型的非線性較強(qiáng)、耦合性較高且導(dǎo)致模型階次較高且未知參數(shù)較多，因此高馬赫數(shù)變體飛行器模型參數(shù)辨識(shí)比傳統(tǒng)飛行器模型參數(shù)辨識(shí)更加復(fù)雜。對(duì)于非線性模型的辨識(shí)問題，可以通過辨識(shí)非線性模型里的線性子系統(tǒng)來得到整個(gè)非線性模型的辨識(shí)結(jié)果。本文選擇對(duì)高馬赫數(shù)變體飛行器的動(dòng)力學(xué)模型進(jìn)行簡(jiǎn)化，忽略耦合因素，選取縱向小擾動(dòng)線性化模型，分析設(shè)計(jì)辨識(shí)激勵(lì)信號(hào)，采用經(jīng)典的遞推最小二乘法對(duì)多輸入多輸出的線性運(yùn)動(dòng)模型進(jìn)行辨識(shí)，通過仿真驗(yàn)證了遞推最小二乘算法的有效性。

1 高馬赫數(shù)變體飛行器辨識(shí)模型

辨識(shí)模型的建立是系統(tǒng)辨識(shí)的一個(gè)重要環(huán)節(jié)，本章建立了高馬赫數(shù)變體飛行器的縱向小擾動(dòng)線性運(yùn)動(dòng)模型，并以此模型為辨識(shí)模型，對(duì)其進(jìn)行了分解。

1.1 縱向小擾動(dòng)線性運(yùn)動(dòng)模型

在對(duì)高馬赫數(shù)變體飛行器模型進(jìn)行參數(shù)辨識(shí)之前首先需要建立辨識(shí)模型?？紤]到高馬赫數(shù)變體飛行器的六自由度運(yùn)動(dòng)方程過于復(fù)雜，多個(gè)參數(shù)之間存在耦合，辨識(shí)難度較大，本文對(duì)高馬赫數(shù)變體飛行器的運(yùn)動(dòng)方程進(jìn)行簡(jiǎn)化，將飛行器的運(yùn)動(dòng)分解為縱向運(yùn)動(dòng)和橫側(cè)向運(yùn)動(dòng)，忽略耦合項(xiàng)?？紤]到飛行器的縱向運(yùn)動(dòng)模型中包含升力系數(shù)、阻力系數(shù)、俯仰力矩系數(shù)等表征飛行能力的主要參數(shù)，本文采用下反角可以變化的高馬赫數(shù)變體飛行器縱向運(yùn)動(dòng)模型為基礎(chǔ)建立辨識(shí)模型。

通過對(duì)高馬赫數(shù)變體飛行器縱向受力進(jìn)行分析，可得到式(1)所示微分方程組。

(1)

式中：m為飛行器質(zhì)量；g為重力加速度；v為速度；X為阻力；Y為升力；α為攻角；?為俯仰角；θ為航跡角；Jz為繞z軸的轉(zhuǎn)動(dòng)慣量；ωz為俯仰角速度；Mz為俯仰力矩；x,y為沿x軸和y軸的位移。

基于建立縱向運(yùn)動(dòng)方程的假設(shè)，忽略二階及三階小量，同時(shí)采用以下基本假設(shè)，將縱向擾動(dòng)運(yùn)動(dòng)方程組線性化：

1) 忽略橫側(cè)向參數(shù)的基準(zhǔn)值及縱向角速度。即假定在未擾動(dòng)運(yùn)動(dòng)中，側(cè)向運(yùn)動(dòng)參數(shù)β、γ、γc、ψ、ψc、ωx、ωy、z，舵偏角δx、δy，縱向運(yùn)動(dòng)參數(shù)中的ωz、α都是小量，可以在方程中忽略其乘積以及這些參數(shù)與其他微量的乘積。此外還假定側(cè)向氣動(dòng)力及其力矩對(duì)縱向參數(shù)的偏導(dǎo)數(shù)在基準(zhǔn)狀態(tài)下均為零。

2) 干擾只改變縱向運(yùn)動(dòng)參數(shù)，不改變側(cè)向參數(shù)，即側(cè)向運(yùn)動(dòng)參數(shù)的偏量為零。

利用偏導(dǎo)數(shù)的簡(jiǎn)略表示法，可得航向靜不穩(wěn)定飛行器縱向小擾動(dòng)線性化模型為

(2)

1.2 辨識(shí)模型的分析

式(2)為一個(gè)MIMO系統(tǒng)，由4個(gè)多輸入單輸出(multiple input single output，MISO)子系統(tǒng)組成。由于這4個(gè)MISO子系統(tǒng)的輸出都毫無關(guān)聯(lián)，并且每個(gè)子系統(tǒng)的結(jié)構(gòu)相同，通過分析研究其中一個(gè)MISO子系統(tǒng)的辨識(shí)問題就可以得到整個(gè)MIMO系統(tǒng)的辨識(shí)結(jié)果。一個(gè)MISO系統(tǒng)可以分解為兩個(gè)單輸入單輸出(single input single output，SISO)系統(tǒng)，并且這兩個(gè)SISO子系統(tǒng)滿足疊加原理，通過分別研究每一個(gè)SISO子系統(tǒng)，最終可以得到整個(gè)MIMO系統(tǒng)的辨識(shí)結(jié)果。

2 參數(shù)辨識(shí)方法

目前，最小二乘法在參數(shù)辨識(shí)領(lǐng)域中是應(yīng)用最為廣泛的方法之一。最小二乘法的改進(jìn)算法——遞推最小二乘法，相比于傳統(tǒng)的最小二乘法，由于狀態(tài)方程非線性較強(qiáng)，可以得到更加精確的辨識(shí)結(jié)果，有效地解決線性系統(tǒng)的模型參數(shù)辨識(shí)問題。對(duì)于非線性系統(tǒng)，可以通過辨識(shí)非線性系統(tǒng)中的線性子系統(tǒng)來獲得整個(gè)非線性系統(tǒng)的辨識(shí)結(jié)果。本文采用時(shí)域遞推最小二乘法對(duì)高馬赫數(shù)變體飛行器的縱向小擾動(dòng)線性模型進(jìn)行參數(shù)辨識(shí)。

2.1 最小二乘法

(3)

基于殘差建立最小二乘指標(biāo)函數(shù)

(4)

(5)

取殘差展開式的一階項(xiàng)得

(6)

聯(lián)立式(4)和式(6)得

(7)

(8)

聯(lián)立式(7)和式(8)得

(9)

將式(9)簡(jiǎn)記為

A(χ-χ0)+B=0

(10)

式中，

(11)

(12)

2.2 遞推最小二乘法

由于系統(tǒng)的狀態(tài)方程和觀測(cè)方程之間的非線性較強(qiáng)，對(duì)式(12)進(jìn)行一次求解難以得到參數(shù)的精確值。因此，利用迭代最小二乘方法進(jìn)行參數(shù)辨識(shí)，迭代格式為

(13)

式中，

(14)

3 辨識(shí)激勵(lì)信號(hào)設(shè)計(jì)

在對(duì)高馬赫數(shù)變體飛行器進(jìn)行系統(tǒng)辨識(shí)之前，需要設(shè)計(jì)適合該系統(tǒng)的辨識(shí)激勵(lì)信號(hào)，合適的激勵(lì)信號(hào)可以有效地激發(fā)系統(tǒng)的各種運(yùn)動(dòng)模態(tài)，從而獲得足夠多的信息，來滿足系統(tǒng)辨識(shí)的需要。不同的輸入激勵(lì)形式會(huì)激發(fā)高馬赫數(shù)變體飛行器不同的運(yùn)動(dòng)模態(tài)。工程上最常見的激勵(lì)信號(hào)為方波信號(hào)，包括3-2-1-1方波、偶極方波、2-1-1方波以及階躍信號(hào)等。在相同的基頻條件下，相比于正弦信號(hào)，方波激勵(lì)信號(hào)的邊緣可以產(chǎn)生更寬的信號(hào)頻譜，是較為理想的系統(tǒng)辨識(shí)激勵(lì)輸入，因此本文采用3-2-1-1方波作為激勵(lì)信號(hào)。

3-2-1-1方波信號(hào)的特點(diǎn)是幅值呈正負(fù)交替出現(xiàn)，寬度比為3∶2∶1∶1，其激勵(lì)形式及功率譜分布如圖1所示。

圖1 3-2-1-1方波輸入激勵(lì)信號(hào)及其頻譜圖Fig.1 3-2-1-1 square wave input excitation signal and its spectrum diagram

由圖1可以看出，3-2-1-1方波信號(hào)的頻帶較寬，并且在一定的頻帶范圍內(nèi)具有較高的能量；3-2-1-1方波信號(hào)正負(fù)交替的特點(diǎn)使飛行器不會(huì)過于偏離原有航跡，并且可以有效地激發(fā)飛行器的各種運(yùn)動(dòng)模態(tài)，是飛行器系統(tǒng)辨識(shí)較為理想的輸入信號(hào)。

4 仿真驗(yàn)證

本章在已建立的高馬赫數(shù)變體飛行器辨識(shí)模型的基礎(chǔ)上把狀態(tài)空間轉(zhuǎn)化成傳遞函數(shù)，將得到的傳遞函數(shù)離散化，最后轉(zhuǎn)化為差分方程，對(duì)高馬赫數(shù)變體飛行器進(jìn)行參數(shù)辨識(shí)。

4.1 由狀態(tài)空間模型構(gòu)建辨識(shí)差分方程

根據(jù)式(2)所示的已建立的飛行器縱向小擾動(dòng)線性運(yùn)動(dòng)模型，定義

(15)

(16)

則式(2)可以表示為

(17)

式(17)描述的是多輸入多輸出系統(tǒng)的模型，本文將其分解為8個(gè)單輸入單輸出的子系統(tǒng)，對(duì)8個(gè)子系統(tǒng)分別進(jìn)行最小二乘遞推辨識(shí)，這種思路便于理解和算法實(shí)現(xiàn)。

將式(17)所示的縱向小擾動(dòng)線性狀態(tài)模型離散化，可以得到其脈沖傳遞函數(shù)為

(18)

傳遞函數(shù)G(z)可由輸出量與輸入量之比表示，即

(19)

式中：Y(z)為輸出量；U(z)為輸入量。

可得系統(tǒng)的差分方程為

b13δe(k-3)+b14δe(k-4)+c11δr(k-1)+c12δr(k-2)+c13δr(k-3)+c14δr(k-4)

(20)

α(k)=-a1α(k-1)-a2α(k-2)-a3α(k-3)-a4α(k-4)+b21δe(k-1)+b22δe(k-2)+

b23δe(k-3)+b24δe(k-4)+c21δr(k-1)+c22δr(k-2)+c23δr(k-3)+c24δr(k-4)

(21)

?(k)=-a1?(k-1)-a2?(k-2)-a3?(k-3)-a4?(k-4)+b31δe(k-1)+b32δe(k-2)+

b33δe(k-3)+b34δe(k-4)+c31δr(k-1)+c32δr(k-2)+c33δr(k-3)+c34δr(k-4)

(22)

ωz(k)=-a1ωz(k-1)-a2ωz(k-2)-a3ωz(k-3)-a4ωz(k-4)+b41δe(k-1)+b42δe(k-2)+

b43δe(k-3)+b44δe(k-4)+c41δr(k-1)+c42δr(k-2)+c43δr(k-3)+c44δr(k-4)

(23)

選取式(20)～式(23)為參數(shù)辨識(shí)模型，采用遞推最小二乘法對(duì)其差分方程中的未知參數(shù)進(jìn)行辨識(shí)。

4.2 仿真實(shí)例

本文以采樣頻率為100 Hz、飛行高度為25 000 m、飛行速度為2 200 m/s的某變下反角飛行器(下反角可以變化的飛行器，下反角為翼梢小翼與機(jī)翼基準(zhǔn)面的夾角)為例對(duì)其縱向小擾動(dòng)線性運(yùn)動(dòng)模型進(jìn)行了遞推最小二乘參數(shù)辨識(shí)。本文仿真實(shí)驗(yàn)中數(shù)據(jù)來源為文獻(xiàn)[9]。

各參數(shù)取值為A11=-0.064 1，A12=0.004 3，A13=-0.004 4，A21=-0.040 4，A22=-0.015 9，A41=1.039 2×10-14，A42=8.806 6×10-16，A44=8.103 1×10-21，B21=0.000 437，B22=0.000 063，B41=0.000 809，B42=0.000 422。

升降舵和下反角輸入激勵(lì)信號(hào)均選取3-2-1-1方波信號(hào)，如圖2所示。

圖2 輸入激勵(lì)信號(hào)Fig.2 Input excitation signal

高馬赫數(shù)變體飛行器飛行過程中的噪聲干擾大部分來自觀測(cè)噪聲，過程噪聲的影響相對(duì)較小，因此本文在仿真過程中只考慮輸出的觀測(cè)噪聲的影響。基于遞推最小二乘辨識(shí)算法，對(duì)式(20)所示的縱向模型子系統(tǒng)進(jìn)行參數(shù)辨識(shí)。

在辨識(shí)過程中在輸出觀測(cè)量上加入期望為0、方差為10×10-7數(shù)量級(jí)的白噪聲。

圖3 攻角系統(tǒng)辨識(shí)過程Fig.3 Attack angle system identification process

在加入噪聲后，由圖3可以看出，在約4 s內(nèi)，攻角系統(tǒng)參數(shù)辨識(shí)結(jié)果收斂，仿真結(jié)果表明遞推最小二乘算法計(jì)算效率較高、收斂速度較快，能夠滿足在線辨識(shí)算法的時(shí)間要求。同時(shí)，由表1可知，辨識(shí)結(jié)果與參數(shù)真值之間誤差小于6%，辨識(shí)結(jié)果達(dá)到10×10-7數(shù)量級(jí)，辨識(shí)精度較高，基本滿足工程要求。

5 結(jié)束語

本文研究了基于遞推最小二乘法的高馬赫數(shù)變體飛行器動(dòng)力學(xué)模型中一些未知參數(shù)的在線識(shí)別和計(jì)算方法。首先，根據(jù)高馬赫數(shù)變體飛行器的六自由度運(yùn)動(dòng)模型，建立了適用于飛行品質(zhì)評(píng)估、控制律設(shè)計(jì)的高馬赫數(shù)變體飛行器縱向小擾動(dòng)線性運(yùn)動(dòng)模型。將多輸入多輸出系統(tǒng)分解為多個(gè)單輸入單輸出子系統(tǒng)，并以此作為飛行器在線辨識(shí)的基礎(chǔ)模型。之后，對(duì)高馬赫數(shù)變體飛行器參數(shù)辨識(shí)所需要的輸入激勵(lì)信號(hào)進(jìn)行了分析與設(shè)計(jì)，選取了可以有效激發(fā)飛行器各種運(yùn)動(dòng)模態(tài)并且能有效避免飛行器過于偏離原有飛行航跡的3-2-1-1方波信號(hào)。最后，以某高馬赫數(shù)變體飛行器為例，采用遞推最小二乘算法，對(duì)飛行器縱向小擾動(dòng)線性運(yùn)動(dòng)模型里的單輸入單輸出線性子系統(tǒng)的在線辨識(shí)進(jìn)行了仿真驗(yàn)證。仿真的過程和結(jié)果表明，本文基于遞推最小二乘法對(duì)飛行器模型參數(shù)進(jìn)行在線辨識(shí)的方法收斂速度快、計(jì)算效率高，辨識(shí)的平均誤差小于6%，辨識(shí)精度較高，整體辨識(shí)精度滿足工程要求。