梁明
(中國(guó)人民解放軍91851 部隊(duì),遼寧 葫蘆島 125001)
隨著計(jì)算機(jī)技術(shù)的不斷發(fā)展,模型辨識(shí)方法與計(jì)算機(jī)技術(shù)相結(jié)合,使得模型辨識(shí)精度越來(lái)越高,在航天飛行器領(lǐng)域的應(yīng)用越來(lái)越廣泛[1]。本文建立了某型飛行器動(dòng)力學(xué)模型結(jié)構(gòu),對(duì)氣動(dòng)力辨識(shí)輸入?yún)?shù)進(jìn)行了分析;采用迭代算法得出辨識(shí)參數(shù),并對(duì)辨識(shí)精度進(jìn)行了分析,認(rèn)為主要是觀測(cè)量測(cè)量誤差、物理幾何參數(shù)誤差影響辨識(shí)精度。模型辨識(shí)的難點(diǎn)在于參數(shù)精度的分析和確認(rèn),本文在初步分析影響精度主要因素的基礎(chǔ)上,將辨識(shí)得到的氣動(dòng)參數(shù)帶入彈道仿真程序進(jìn)行仿真驗(yàn)證。
在建立模型過(guò)程中,將某型飛行器視為運(yùn)動(dòng)剛體,其在空間的運(yùn)動(dòng)可以分解為3 個(gè)線(xiàn)位移和3 個(gè)角位移,用6 自由度來(lái)進(jìn)行描述,采取輸入激勵(lì)信號(hào)激發(fā)某型飛行器運(yùn)動(dòng)模態(tài),從而辨識(shí)出該型飛行器的相關(guān)氣動(dòng)參數(shù)[2]。
某型飛行器助推器工作時(shí)間比總飛行時(shí)間要短的多,一般只有3~4 s,但助推器的推力很大,通常助推器推力產(chǎn)生的縱向過(guò)載可達(dá)15~18g,在助推器推力持續(xù)作用時(shí),空氣動(dòng)力的影響幾乎可以忽略,且助推段很難獲得滿(mǎn)足參數(shù)辨識(shí)要求的響應(yīng)參數(shù)[3],因此本研究中不考慮助推器的推力。在上述條件下,某型飛行器6 自由度動(dòng)力學(xué)數(shù)學(xué)模型為[4]
式中:Wxd、Wzd為風(fēng)速在地面座標(biāo)系的分量,Wxd=WcosψW,Wzd=-WsinψW,ψW為風(fēng)向角。
氣動(dòng)系數(shù)的模型應(yīng)根據(jù)某型飛行器氣動(dòng)外形的特點(diǎn)來(lái)確定[5],一般情況下可表示為
以上方程再加上控制方程就成為一個(gè)完整的方程組。某型飛行器的調(diào)節(jié)規(guī)律方程通??杀硎緸椋?]
從調(diào)節(jié)規(guī)律的方程還可以看出,某型飛行器3 個(gè)通道的控制是互相獨(dú)立的,正常飛行情況下,可以將6 自由度動(dòng)力學(xué)數(shù)學(xué)模型,簡(jiǎn)化為縱向和側(cè)向2 個(gè)3 自由度的動(dòng)力學(xué)模型[8]。
1)導(dǎo)彈的結(jié)構(gòu)參數(shù):導(dǎo)彈質(zhì)量m,轉(zhuǎn)動(dòng)慣量Jx、Jy、Jz,導(dǎo)彈質(zhì)心XT、YT,參考面積S,平均氣動(dòng)弦長(zhǎng)bA,參考長(zhǎng)度l,傳感器的安裝位置等;
2)外彈道測(cè)量參數(shù):導(dǎo)彈速度V、飛行高度h等;
3)氣象測(cè)量參數(shù):密度ρ(若密度無(wú)實(shí)測(cè)值可根據(jù)飛行高度計(jì)算)、風(fēng)速W、風(fēng)向角ψW等;
4)內(nèi)彈道測(cè)量參數(shù):舵偏角δx、δx、δx,導(dǎo)彈姿態(tài)角ψ、?、γ,導(dǎo)彈轉(zhuǎn)動(dòng)角速度ωx、ωy、ωz,導(dǎo)彈軸向、法向、側(cè)向過(guò)載Nx、Ny、Nz,發(fā)動(dòng)機(jī)推力R,導(dǎo)彈攻角α,側(cè)滑角β等。
采用牛頓-拉夫遜迭代公式[9]
式中的d修正量Δθk可由下列線(xiàn)性代數(shù)方程組計(jì)算[10]:
具體迭代過(guò)程為根據(jù)某型飛行器理論計(jì)算結(jié)果,給出待估氣動(dòng)參數(shù)的初值θ0,由狀態(tài)方程、觀測(cè)方程和靈敏度方程積分出某型飛行器的狀態(tài)值x、觀測(cè)值y和靈敏度陣;然后,解線(xiàn)性代數(shù)方程組求出Δθ0,再以θ1=θ0+Δθ0代替原來(lái)的θ0,重復(fù)以上計(jì)算過(guò)程,每迭代一次都需要計(jì)算似然準(zhǔn)則函數(shù)Jk和Jk/Jk-1[11]。當(dāng)時(shí),則認(rèn)為迭代收斂,此時(shí)的待估參數(shù)θk,即為所求的氣動(dòng)參數(shù),一般情況下取ε=0.01[12]。
影響氣動(dòng)參數(shù)辨識(shí)精度的因素較多,歸納起來(lái)有以下幾個(gè)方面:
1)觀測(cè)量的測(cè)量誤差。觀測(cè)量的測(cè)量誤差是影響辨識(shí)精度的重要因素[13]。利用飛行試驗(yàn)數(shù)據(jù)進(jìn)行參數(shù)辨識(shí)時(shí),觀測(cè)量一般取導(dǎo)彈轉(zhuǎn)動(dòng)角速度ωx、ωy、ωz,姿態(tài)角ψ、?、γ,導(dǎo)彈軸向、法向、側(cè)向過(guò)載Nx、Ny、Nz,導(dǎo)彈攻角α,側(cè)滑角β。分別由安裝在彈體內(nèi)部的角速率傳感器、角位移傳感器、過(guò)載傳感器、攻角傳感器及側(cè)滑角傳感器獲得。通過(guò)遙測(cè)傳感器獲得時(shí),零位漂移、死區(qū)、安裝誤差不可避免,但選擇高精度的傳感器可以使誤差減小,為滿(mǎn)足氣動(dòng)參數(shù)辨識(shí)的需要,遙測(cè)傳感器應(yīng)滿(mǎn)足如下要求:
(1)噪聲強(qiáng)度系數(shù)不高于1%;
(2)角位移傳感器的零位漂移速度應(yīng)小于每分鐘0.1°,角速率傳感器的零位漂移小于峰值的0.2%,過(guò)載傳感器零位漂移應(yīng)小于峰值的1.0%;
(3)傳感器的標(biāo)定誤差應(yīng)小于0.3%;
(4)閉環(huán)飛行試驗(yàn)中,角速率、角位移、過(guò)載傳感器的死區(qū)分別應(yīng)在0.03(°/s,(°),g)的范圍內(nèi);
(5)角速率傳感器安裝角的偏移度要小于0.1°,過(guò)載傳感器安裝角的偏移度要小于0.2°,過(guò)載傳感器質(zhì)心位置及安裝位置誤差要小于3.0 mm。
2)幾何和物理參數(shù)的誤差。幾何參數(shù)包括特征長(zhǎng)度、面積,物理參數(shù)包括慣性力矩、動(dòng)壓、質(zhì)量,這些量中任何一個(gè)量值存在的誤差,都會(huì)對(duì)氣動(dòng)參數(shù)估值產(chǎn)生相同的誤差[14]。
若右端某一項(xiàng)遠(yuǎn)小于左端,則該項(xiàng)是不可辨識(shí)的;若左端接近于零,則該項(xiàng)也是不可辨識(shí)的,其余類(lèi)同。
數(shù)據(jù)取自該型飛行器6 自由度85~86 s 的內(nèi)、外彈道數(shù)據(jù),數(shù)據(jù)包括:光測(cè)數(shù)據(jù)t、V、h、θ、ψc,遙測(cè)采樣數(shù)據(jù)ωx、ωy、ωz、Nx、Ny、Nz、δx、δy、δz、α、β,彈體結(jié)構(gòu)參數(shù)l、bA、S、m、Jx、Jy、Jz,氣象數(shù)據(jù)ρ、W、ψW等。
為了分析這段數(shù)據(jù)的可辨識(shí)性,將動(dòng)力學(xué)方程的各項(xiàng)分別逐項(xiàng)計(jì)算,結(jié)果見(jiàn)表1—表3。可以看出與各項(xiàng)的數(shù)量級(jí)相當(dāng),這表明縱、側(cè)向氣動(dòng)系數(shù)均可進(jìn)行辨識(shí),可以采用6 自由度參數(shù)辨識(shí)的基本方程組進(jìn)行辨識(shí),表中數(shù)據(jù)還表明阻尼參數(shù)的數(shù)量級(jí)較小,會(huì)有較大的辨識(shí)誤差[16]。
表1 辨識(shí)數(shù)據(jù)x表Tab.1 x table of identification data
表1 辨識(shí)數(shù)據(jù)x表Tab.1 x table of identification data
表2 辨識(shí)數(shù)據(jù)y表Tab.2 y table of identification data
表2 辨識(shí)數(shù)據(jù)y表Tab.2 y table of identification data
表3 辨識(shí)數(shù)據(jù)z表Tab.3 z table of identification data
表3 辨識(shí)數(shù)據(jù)z表Tab.3 z table of identification data
采用極大似然法,對(duì)仿真數(shù)據(jù)采取參數(shù)辨識(shí)的基本方法進(jìn)行辨識(shí)[17]。
狀態(tài)變量為
觀測(cè)量為
待估參數(shù)為
共辨識(shí)出21個(gè)氣動(dòng)系數(shù),辨識(shí)結(jié)果見(jiàn)表4~表7。
表4 辨識(shí)結(jié)果(1~5)Tab.4 Identification results(1~5)
表6 辨識(shí)結(jié)果(11~16)Tab.6 Identification results(11~16)
表7 辨識(shí)結(jié)果(17~21)Tab.7 Identification results(17~21)
辨識(shí)得到的氣動(dòng)參數(shù)為計(jì)算飛行器彈道導(dǎo)彈飛行力學(xué)所需要的升力、阻力、側(cè)向力,俯仰、滾轉(zhuǎn)、偏航力矩系數(shù)及其對(duì)迎角、側(cè)滑角、舵偏角的偏導(dǎo)。辨識(shí)得到這些氣動(dòng)參數(shù)之后,便可以對(duì)導(dǎo)彈的彈道進(jìn)行仿真分析[18]。
通過(guò)對(duì)該型飛行器控制系統(tǒng)和動(dòng)力系統(tǒng)進(jìn)行辨識(shí),得出了部分模型辨識(shí)參數(shù),下面驗(yàn)證這些實(shí)辨參數(shù)模型。
根據(jù)理論模型編寫(xiě)彈道仿真程序,可有效地計(jì)算該型號(hào)飛行器的控制彈道,并經(jīng)過(guò)多次實(shí)彈飛行試驗(yàn)驗(yàn)證。驗(yàn)證結(jié)果表明,該理論模型精度較高,滿(mǎn)足要求。
將實(shí)辨模型代入該型飛行器彈道仿真程序,替換原來(lái)的理論模型,進(jìn)行彈道仿真計(jì)算。從仿真結(jié)果看,彈道變化平穩(wěn),與試驗(yàn)數(shù)據(jù)吻合較好,誤差符合規(guī)定的范圍,通過(guò)本方法得到的氣動(dòng)力參數(shù)精度已滿(mǎn)足工程實(shí)際需要,實(shí)際辨識(shí)模型參數(shù)具有可信性[19]。Y方向仿真彈道如圖1 所示,該彈道反映了某型飛行器初始正常爬高,按照預(yù)定程序下滑、穩(wěn)定平飛、91 s 開(kāi)始規(guī)避爬升、躍過(guò)目標(biāo)艦后下滑入水自毀。
圖1 Y 方向仿真彈道Fig.1 Simulation ballistics in the Y-direction
等效飛行高度曲線(xiàn)如圖2 所示,該曲線(xiàn)反映了某型飛行器高度通道輸出變化規(guī)律,從0~90 s 高度通道正常輸出高度表測(cè)量高度,90 s 以后輸出固定高度電平,從而完成安控規(guī)避。X方向仿真彈道曲線(xiàn)如圖3 所示,該曲線(xiàn)反映了該型號(hào)飛行器縱向位移平穩(wěn)。
圖2 等效飛行高度仿真Fig.2 Equivalent flight altitude simulation
圖3 X 方向仿真彈道Fig.3 Simulation ballistics in the X-direction
1~4 舵機(jī)仿真曲線(xiàn)分別如圖4~圖7 所示,4 條仿真曲線(xiàn)反映了某型飛行器舵機(jī)系統(tǒng)的變化規(guī)律,這4 個(gè)舵機(jī)系統(tǒng)動(dòng)作對(duì)稱(chēng)性好,響應(yīng)速度快,滿(mǎn)足該型飛行器控制要求[20]。
圖4 1 舵機(jī)仿真曲線(xiàn)Fig.4 Simulation curve of Servo 1
圖5 2 舵機(jī)仿真曲線(xiàn)Fig.5 Simulation curve of Servo 2
圖6 3 舵機(jī)仿真曲線(xiàn)Fig.6 Simulation curve of Servo 3
圖7 4 舵機(jī)仿真曲線(xiàn)Fig.7 Simulation curve of Servo 4
本文介紹了一種飛行器氣動(dòng)模型參數(shù)辨識(shí)方法,利用該方法進(jìn)行了飛行器的典型氣動(dòng)參數(shù)辨識(shí),并將辨識(shí)得到的氣動(dòng)數(shù)據(jù)代入彈道仿真程序中進(jìn)行驗(yàn)證。驗(yàn)證結(jié)果表明,該氣動(dòng)參數(shù)得到的彈道與實(shí)際情況吻合程度高,辨識(shí)精度較好。同時(shí)作為一種工程方法,該氣動(dòng)參數(shù)辨識(shí)方法使用較為方便,計(jì)算耗時(shí)較少,適合于多種工程場(chǎng)景,在工程領(lǐng)域有廣泛的應(yīng)用前景。