李金晟,常思江,陳升富

(1.北京航天長征飛行器研究所,北京 100076;2.南京理工大學(xué) 能源與動力工程學(xué)院,江蘇 南京 210094)

氣動參數(shù)對于常規(guī)彈箭的飛行穩(wěn)定性及彈道軌跡具有重要影響。目前,彈箭的氣動參數(shù)的估計方法已有許多,其中經(jīng)典常用的有基于最小二乘準(zhǔn)則的C-K法[1]、基于極大似然準(zhǔn)則的輸出誤差法[2]以及結(jié)合了濾波思想的濾波誤差法[3]。但是,這些方法在應(yīng)用于實測飛行數(shù)據(jù)時,無一例外都會出現(xiàn)一些復(fù)雜的問題。其中較為普遍的一個難題便是對于真實系統(tǒng)缺少精確的數(shù)學(xué)模型。大多數(shù)傳統(tǒng)的方法利用動力學(xué)模型來近似描述真實系統(tǒng),這種近似會帶來建模誤差的問題。但對于建模誤差的處理目前還沒有成熟的理論可用,通常是將其簡單的歸為狀態(tài)或測量噪聲來處理,這種處理方式?jīng)]有考慮到建模誤差可能是確定性的而非隨機性的[4]。為了解決上述傳統(tǒng)方法存在的建模誤差問題,相關(guān)人員開始采用人工神經(jīng)網(wǎng)絡(luò)(artificial neural networks,ANN)進行氣動參數(shù)辨識,其完全基于對輸入輸出的測量而不考慮系統(tǒng)內(nèi)部的物理機制。

基于神經(jīng)網(wǎng)絡(luò)的氣動參數(shù)辨識方法主要有三類:一是基于動力學(xué)模型的辨識方法[5],其優(yōu)點是無需預(yù)先計算基準(zhǔn)氣動數(shù)據(jù),缺點是所辨識的氣動模型仍是一個“黑箱”,無法得到具體表達式;二是基于穩(wěn)定或控制導(dǎo)數(shù)的Delta法和Zeros法[6],其優(yōu)點是可以辨識出氣動模型,缺點是需預(yù)先得知較為準(zhǔn)確的氣動參數(shù)數(shù)據(jù),從而利用其訓(xùn)練神經(jīng)網(wǎng)絡(luò);三是將前饋神經(jīng)網(wǎng)絡(luò)(feed forward neural networks,FFNNs)與優(yōu)化算法高斯-牛頓法(Gauss-Newton)相結(jié)合的神經(jīng)網(wǎng)絡(luò)-牛頓(neural Gauss-Newton,NGN)辨識方法[7],其獨特之處在于無需預(yù)先對系統(tǒng)進行建模,不需求解運動方程。這三類方法已經(jīng)成功應(yīng)用于飛機、導(dǎo)彈、再入體等飛行器的氣動參數(shù)辨識,但尚沒有應(yīng)用于彈丸。

本文針對某無控旋轉(zhuǎn)彈丸的阻力系數(shù)辨識問題,采用基于誤差反向傳播(back propagation,BP)神經(jīng)網(wǎng)絡(luò)的NGN法進行研究,并對結(jié)果進行了分析。

1 BP網(wǎng)絡(luò)辨識動力學(xué)模型

1.1 BP網(wǎng)絡(luò)原理簡介

BP神經(jīng)網(wǎng)絡(luò)是一種單向傳播的多層前向網(wǎng)絡(luò),其分為輸人層、隱含層和輸出層。隱含層可以有若干層,每一層的神經(jīng)元只接收前一層神經(jīng)元的輸出,同一層神經(jīng)元之間不連接。該網(wǎng)絡(luò)的拓?fù)浣Y(jié)構(gòu)如圖1所示[8]。

1.2 辨識動力學(xué)模型

NGN算法區(qū)別于傳統(tǒng)辨識算法的獨特之處在于,利用神經(jīng)網(wǎng)絡(luò)強大的非線性逼近能力來減小建模誤差。由于本文主要是辨識彈丸的零升阻力系數(shù)Cx0,且文獻[2]已表明通過量測速度數(shù)據(jù)v便可實現(xiàn)對Cx0的辨識。故利用網(wǎng)絡(luò)辨識出速度狀態(tài)滿足的關(guān)系式即可。

BP神經(jīng)網(wǎng)絡(luò)在對狀態(tài)進行預(yù)測之前首先要對網(wǎng)絡(luò)進行多次訓(xùn)練,使其具備聯(lián)想記憶和逼近能力,具體的訓(xùn)練過程包括以下幾個步驟[9]。

第1步:網(wǎng)絡(luò)初始化。由質(zhì)點彈道模型可得,加速度主要與速度、密度及零升阻力系數(shù)有關(guān),而在標(biāo)準(zhǔn)條件下密度又是高度的函數(shù)。因此本文的輸入向量為

X(t)=(vx(t)vy(t)vz(t)y(t)Cx0(t))

(1)

式中:t為離散的采樣時間,Cx0(t)是通過對風(fēng)洞吹風(fēng)或理論計算所得的基準(zhǔn)氣動參數(shù)進行線性插值所得?；赬(t),輸出向量為

Y(t+1)=(vx(t+1)vy(t+1)vz(t+1)y(t+1))

(2)

由輸入、輸出向量可得輸入層節(jié)點n=5,輸出層節(jié)點數(shù)m=4。

第2步:求隱含層輸出H。由X(t),輸入層和隱含層間連接權(quán)值wij及隱含層閾值a,求H。

(3)

式中:l為隱含層節(jié)點數(shù);f(*)為隱含層激勵函數(shù),本文設(shè)置其為Sigmoid函數(shù),表達式為

(4)

式中:β>0,本文取β=1。

第3步:求輸出層輸出Y(t+1)。由H,隱含層與輸出層連接權(quán)值wjk及輸出層閾值b,求預(yù)測輸出Y(t+1)。

(5)

第4步:誤差計算。由Y(t+1)和期望輸出O(t+1),求出預(yù)測誤差e及誤差性能指標(biāo)函數(shù)E。

(6)

第5步:權(quán)值與閾值更新?；谧钏傧陆捣?由e更新連接權(quán)值wij,wjk以及閾值a和b。

(7)

(8)

η(g)=ηmax-[g(ηmax-ηmin)/gmax]

(9)

式中:ηmax,ηmin分別為最大、最小學(xué)習(xí)速率;gmax為最大迭代次數(shù);g為當(dāng)前迭代次數(shù)。

第6步:判斷算法迭代是否結(jié)束,若沒有結(jié)束,返回第2步。為了提高迭代效率,本文設(shè)置了2個迭代終止準(zhǔn)則,準(zhǔn)則一為達到最大迭代次數(shù)gmax,準(zhǔn)則二為誤差性能函數(shù)E小于給定精度。

神經(jīng)網(wǎng)絡(luò)訓(xùn)練結(jié)束后,所得到的BP網(wǎng)絡(luò)就可以看成一個非線性動力學(xué)模型,其表達了從n個輸入自變量到m個輸出因變量的函數(shù)映射關(guān)系。

2 NGN算法應(yīng)用過程

當(dāng)神經(jīng)網(wǎng)絡(luò)的相關(guān)參數(shù)調(diào)試恰當(dāng)后,就可以對t時刻的任意測量輸入變量X(t)進行預(yù)測,得到t+1時刻的輸出變量Y(t+1)。接下來便可根據(jù)神經(jīng)網(wǎng)絡(luò)模型辨識所需的氣動參數(shù)。

2.1 分段辨識思想

本文基于分段辨識思想提取阻力系數(shù),即將雷達所測的速度數(shù)據(jù)劃分為若干個小區(qū)間,在每一個小區(qū)間內(nèi)由于速度變化很小,故可把該區(qū)間內(nèi)的Cx0視為常數(shù)。

2.2 辨識準(zhǔn)則

本文采用極大似然準(zhǔn)則做為目標(biāo)準(zhǔn)則函數(shù)。在實際應(yīng)用中由于彈丸自由飛行試驗要求較為嚴(yán)格,這時系統(tǒng)的過程噪聲很小,故目標(biāo)函數(shù)可表示為

(10)

V(q)=ym(q)-y(q)

(11)

式中:ym(q)=(vm,x(q)vm,y(q)vm,z(q))T為已經(jīng)過坐標(biāo)轉(zhuǎn)換的雷達量測速度數(shù)據(jù);y(q)=(vx(q)vy(q)vz(q))T是基于當(dāng)前輸入,通過BP神經(jīng)網(wǎng)絡(luò)所計算的理論速度值。

2.3 靈敏度參數(shù)計算

θ(r+1)=θ(r)+Δθ(r)

(12)

傳統(tǒng)極大似然法是采用數(shù)值積分法聯(lián)立求解彈道微分方程組與靈敏度方程組,從而獲得靈敏度系數(shù)值。由于需要建立靈敏度方程組,故其較為繁瑣,尤其是當(dāng)模型比較復(fù)雜時,計算量會急劇增大。而NGN憑借其獨特優(yōu)勢,極大地簡化了計算靈敏度系數(shù)的過程,其方式如下:

(13)

式中:d為辨識系數(shù)個數(shù);對于本文而言d=1;Ypk(t+1)為擾動輸出,即給當(dāng)前輸入θd加一任意小擾動Δθd,利用已訓(xùn)練好的神經(jīng)網(wǎng)絡(luò)計算得到的輸出值。

至此,利用式(13)計算的靈敏度值可求出修正量ΔCx0。重復(fù)以上的計算過程,每迭代一次都需要由似然準(zhǔn)則函數(shù)計算Qr,直到當(dāng)|1-Qr/Qr-1|≤0.01時,迭代收斂,計算結(jié)束。此時輸出的Cx0就是該彈道段內(nèi)辨識計算所得到的阻力參數(shù)。

3 數(shù)值結(jié)果分析

3.1 仿真驗證

為了檢驗NGN算法的可行性,本文針對無控旋轉(zhuǎn)穩(wěn)定彈對其進行仿真分析,相關(guān)參數(shù)設(shè)計如下。①訓(xùn)練神經(jīng)網(wǎng)絡(luò)所需數(shù)據(jù)是由質(zhì)點彈道模型計算所得,初始彈道諸元設(shè)置為:初速v0=1 100 m/s,射角θa0=45°,射向角ψ20=0,積分步長Δt=0.01 s;②BP神經(jīng)網(wǎng)絡(luò)調(diào)整參數(shù):隱層層數(shù)為1,隱層節(jié)點數(shù)為l=10,學(xué)習(xí)速率ηmax=0.7,ηmin=0.4,動量因子α=0.1,權(quán)值wij、wjk和閾值a、b在[-1,1]間隨機初始化,迭代終止準(zhǔn)則一為最大迭代次數(shù)gmax=1 000,準(zhǔn)則二為誤差性能函數(shù)E≤0.001。

基于上述參數(shù)的設(shè)計,首先分析BP神經(jīng)網(wǎng)絡(luò)對動力學(xué)模型的逼近精度。本文基于質(zhì)點彈道模型共生成7 713組彈道數(shù)據(jù),隨機選擇其中5 000組數(shù)據(jù)來訓(xùn)練神經(jīng)網(wǎng)絡(luò),利用其余的數(shù)據(jù)來測試神經(jīng)網(wǎng)絡(luò)的逼近能力。質(zhì)點彈道的速度規(guī)律如圖2所示。

本文采用均方誤差(mean square error,MSE)及均方根誤差(root mean square error,RMSE)來反映神經(jīng)網(wǎng)絡(luò)的預(yù)測精度,表達式為

(14)

BP神經(jīng)網(wǎng)絡(luò)的MSE及RMSE值如表1所示。

表1 神經(jīng)網(wǎng)絡(luò)模型的MSE及RMSE值

由表1可得,本文所建立的神經(jīng)網(wǎng)絡(luò)模型預(yù)測的3個速度分量值、高度值與真實值基本一致,且RMSE值相比于狀態(tài)值而言非常小,這表明辨識效果較好。此外,神經(jīng)網(wǎng)絡(luò)的性能,即預(yù)測能力及收斂速度,則取決于網(wǎng)絡(luò)的調(diào)整參數(shù)。經(jīng)過調(diào)試分析可知,隱層層數(shù)、隱層節(jié)點數(shù)、最大迭代次數(shù)主要決定了BP神經(jīng)網(wǎng)絡(luò)的預(yù)測精度,而初始權(quán)值及閥值、學(xué)習(xí)速率和動量因子主要影響B(tài)P神經(jīng)網(wǎng)絡(luò)的收斂速度。

網(wǎng)絡(luò)訓(xùn)練成功后,便可采用NGN法提取零升阻力系數(shù)。本文選擇7 500組仿真雷達測速數(shù)據(jù)(其是由質(zhì)點彈道模型與高斯白噪聲疊加而成)。為了滿足每一彈道段內(nèi)阻力系數(shù)可以近似為常數(shù)這一條件,本文共設(shè)置500個彈道段,每個彈道段內(nèi)最大速度降約為20 m/s。選定雷達量測噪聲為σvx=σvy=σvz=0,辨識結(jié)果如圖3所示。圖中,Cd為阻力系數(shù),Ma為馬赫數(shù)。

由圖2可得,基于NGN法辨識的阻力系數(shù)與基準(zhǔn)阻力系數(shù)基本一致,辨識精度較高,可滿足要求。相對而言,在跨聲速段的辨識誤差較大,這是分段小區(qū)間常數(shù)法的固有缺陷。經(jīng)過調(diào)試可知,神經(jīng)網(wǎng)絡(luò)的預(yù)測精度決定了NGN算法的辨識精度,而靈敏度參數(shù)的設(shè)置決定了NGN算法的收斂速度。

3.2 實測數(shù)據(jù)處理

接下來分析NGN法對雷達實測飛行數(shù)據(jù)的辨識效果。本文所用雷達的測量誤差:σγ=5 m,σβ=0.001 rad,σε=0.001 rad,σv=2 m/s。首先直接采用雷達量測數(shù)據(jù)訓(xùn)練BP神經(jīng)網(wǎng)絡(luò),訓(xùn)練并調(diào)試成功后,其MSE及RMSE值如表2所示。

表2 神經(jīng)網(wǎng)絡(luò)模型的MSE及RMSE值

將表1與表2進行對比可得,BP網(wǎng)絡(luò)基于實測數(shù)據(jù)的預(yù)測精度較低,其原因可能如下:一是神經(jīng)網(wǎng)絡(luò)參數(shù)設(shè)置仍不是最優(yōu);二是由于訓(xùn)練數(shù)據(jù)是雷達實測數(shù)據(jù),結(jié)合六自由度剛體彈道模型可得,此時的狀態(tài)量還受到其他參量的影響,如速度狀態(tài)還受到攻角、彈道傾角、升力系數(shù)導(dǎo)數(shù)等的影響,而這些參量由于傳感器的原因,并沒有被作為輸入輸出量進行考慮;三是訓(xùn)練網(wǎng)絡(luò)使用的基準(zhǔn)阻力系數(shù)并不是該發(fā)試驗彈丸最精確的系數(shù)。

NGN法辨識的零升阻力系數(shù)Cx0與馬赫數(shù)Ma的關(guān)系如圖4所示。

由圖4可得,NGN法辨識的阻力系數(shù)與基準(zhǔn)阻力系數(shù)在馬赫數(shù)為[2,2.8]之間時,曲線基本一致,偏差較小;在馬赫數(shù)為[1.2,2]之間時,辨識曲線已難以跟隨基準(zhǔn)曲線,誤差較大;在馬赫數(shù)為[1,1.2]之間時,辨識曲線圍繞基準(zhǔn)曲線振蕩不斷。出現(xiàn)上述現(xiàn)象的原因可能是:一是神經(jīng)網(wǎng)絡(luò)辨識的模型精度不高,這一點已在前文進行了分析;二是分段小區(qū)間常數(shù)法的固有缺陷使然;三是雷達量測數(shù)據(jù)存在野值,可考慮針對濾波數(shù)據(jù)進行訓(xùn)練和辨識。

為了進一步分析NGN法的應(yīng)用特點,本文又采用目前工程應(yīng)用已成熟的擴展卡爾曼濾波(extended Kalman filter,EKF)方法處理該發(fā)彈丸的實測數(shù)據(jù),將分段辨識所得的阻力符合系數(shù)與對應(yīng)彈道段內(nèi)的基準(zhǔn)系數(shù)相乘,便得到全彈道的阻力系數(shù)圖,如圖5所示。

由圖5可得,EKF前期辨識精度較差,其原因是EKF法需要一定的數(shù)據(jù)量才能收斂于穩(wěn)定值。

綜上所述,在實際應(yīng)用中可將NGN和EKF這2種方法進行結(jié)合,取長補短,從而達到工程要求。基于上述思想,本文在馬赫數(shù)為[2.2,2.8]時,采用NGN法辨識的阻力系數(shù);在馬赫數(shù)為[1,2.2]時,采用EKF辨識的結(jié)果。然后對雷達測量彈道、基于理論基準(zhǔn)阻力系數(shù)計算的彈道、基于辨識阻力系數(shù)計算的彈道三者的速度進行比較,如圖6所示。

理論參數(shù)計算彈道、辨識參數(shù)計算彈道與雷達測量彈道的速度及位置的偏差均值如表3所示。

表3 狀態(tài)偏差均值表

由圖6及表3可得,分段辨識參數(shù)計算的彈道相比理論基準(zhǔn)參數(shù)計算的彈道而言,其與雷達測量彈道的接近程度更高。因此,可認(rèn)為分段辨識所得的阻力系數(shù)比基準(zhǔn)阻力系數(shù)更為準(zhǔn)確,從而說明結(jié)合使用NGN法與EKF法具有一定的應(yīng)用意義。

4 結(jié)束語

本文首次將基于神經(jīng)網(wǎng)絡(luò)的參數(shù)辨識方法應(yīng)用于彈丸的氣動參數(shù)辨識。通過對仿真結(jié)果和實際辨識結(jié)果進行分析可得,當(dāng)輸入輸出量測數(shù)據(jù)較為完整且網(wǎng)絡(luò)的調(diào)整參數(shù)設(shè)置恰當(dāng)時,NGN法能夠以較高的精度及可靠性辨識出阻力系數(shù);但當(dāng)所獲量測狀態(tài)數(shù)據(jù)并不完整時,NGN法可用于補償EKF法,提高其前期的辨識精度,從而滿足工程應(yīng)用要求。本文研究結(jié)果可為神經(jīng)網(wǎng)絡(luò)在彈丸氣動參數(shù)辨識中的實際應(yīng)用提供一定的理論參考。