基于MPSP算法的炮彈分段氣動參數辨識

柴 勁,張鵬飛,張 意,齊竹昌,王天明,韓旭東

(西安現代控制技術研究所, 西安 710000)

0 引言

升阻特性作為飛行器的基本氣動特性,對總體射程指標以及滑翔方案彈道設計等影響較大,因此在項目前期論證以及飛行試驗中被重點關注[1-4]。在得到飛行試驗數據之后,工程實踐中常用2種方法進行參數辨識:一種方法是彈道復合法,即在彈道程序中對氣動插值結果乘以修正系數,以使速度、彈道傾角、彈道等外彈道參數接近實測值,這種方法工作量較大且過程繁瑣,辨識過程依賴工程師的經驗判斷,辨識結果精度有限;另一種辨識方法是通過外彈道數據求得彈體的受力信息,根據質心動力學方程直接求解升力阻力,從而得到阻力系數與升力系數的辨識值,該方法簡單方便,但需要對觀測數據直接進行微分,辨識噪聲及辨識誤差較大。

目前常見的氣動辨識算法有最小二乘法[5-7]、極大似然估計[8-9]、卡爾曼濾波[10-11]以及神經網絡[12-13]、蜂群算法[14]、基于學生心理的啟發式優化算法[15]等智能算法。這些智能算法模型復雜且對觀測量及其誤差分布等信息需求較多,不利于工程應用。

本文中的辨識算法受模型預測靜態規劃(MPSP)[16-17]啟發,該算法多用于帶有終端約束的制導律設計,通過狀態靈敏度矩陣建立了控制變量與終端誤差的聯系。本文中將待辨識氣動參數類比制導律中的控制量,將辨識誤差類比制導律中的終端約束誤差,建立起辨識氣動參數與外彈道數據誤差之間的聯系,通過算法自動調整氣動參數,從而實現仿真彈道與外彈道數據誤差小于指定精度。本文中提出的氣動辨識方法簡單新穎,可通過編制程序自動執行,無需通過人工手動復合彈道,減少辨識工作量,提高辨識精度,為飛行試驗的參數辨識提供了一種新思路。

1 建立辨識模型

1.1 辨識動力學模型

暫時不考慮有助推發動機的情況,使用簡化的縱向彈體質心動力學方程作為辨識模型,即:

(1)

式中:V、θ、x、y、q、S、cx、cy、m、g分別為速度、彈道傾角、發射系射距、發射系高度、動壓、參考面積、阻力系數、升力系數、質量、重力加速度。

1.2 歸一化處理

由于狀態變量存在量級差異,先進行歸一化處理。使用最大、最小標準化處理,歸一化處理算子為:

(2)

將歸一化后的變量代入式(1)中,即得到歸一化后的縱向質心運動模型。

2 MPSP氣動參數辨識算法

2.1 算法原理

首先確定辨識步長,在每個辨識步長內應用MPSP算法進行迭代計算,使狀態誤差收斂至誤差閾值以下,即得到待辨識參數的估計值。算法流程圖如圖1所示。

圖1 算法流程Fig.1 Algorithm flow chart

2.2 提供辨識輸入

提供外彈道數據、理論氣動參數以及彈體結構參數等辨識輸入。

2.3 設置辨識步長

辨識步長設置要根據飛行狀態確定。若彈體全程飛行狀態變化范圍不大,辨識步長可以適當取長;若飛行狀態變化劇烈,辨識步長需減小。

2.4 MPSP算法用于氣動參數辨識

(3)

采用歐拉法對式(3)進行離散處理(為簡化表示,忽略步長上標p,不失一般性,下同),得:

Xk=Xk-1+hf(Xk-1,λk-1)=Fk-1(Xk-1,λk-1)

(4)

式中:h為數值積分步長;k=1,2,…,N為離散點。

考慮式(4)全微分形式為:

(5)

式中,

(6)

(7)

(8)

(9)

令k=N,N-1,…,2,將式(5)循環代入式(9)中可得終端誤差近似為:

(10)

式(10)中,

(11)

由于初始狀態量為常量,故dX1=0。式(10)可簡化為:

(12)

式(12)通過靈敏度矩陣Bk建立了終端誤差ΔYN與待辨識參數變化量dλk之間的聯系。其物理含義是,當實際外彈道數據與辨識模型的狀態量間存在誤差時,可通過調整辨識模型中的待辨識參數去減小該誤差,直到其滿足精度要求。理論上滿足誤差要求的辨識參數組合并不唯一。但存在一組既滿足誤差要求,同時與辨識初值相比變化量又最小的一組解作為“較優”的辨識結果。因此可以通過構造帶有誤差要求約束、最小化辨識參數變化量的優化問題。其目標函數為:

(13)

式中,Rk為權重矩陣,可對辨識參數的調整量影響進行加權。對式(12)和式(13)組成的優化問題,可通過拉格朗日乘子法轉化為無約束優化問題求解。省去求解過程,直接給出最優解的解析形式為:

(14)

待辨識參數更新為:

(15)

至此,基于MPSP的氣動辨識算法推導完畢。將其整理為:

步驟1提供辨識算法需要的外彈道數據、初始氣動參數及彈體相關參數作為辨識輸入。

步驟2設置辨識步長。

步驟4計算矩陣Bk、Aλ、Bk。

3 仿真結果與分析

3.1 仿真條件設置

以某型155 mm制導炮彈為例,驗證本文中辨識算法的有效性。主要彈體參數及仿真參數設置如表1所示。其他仿真條件設置為:大氣模型使用理想大氣。使用真實阻力系數和升力系數進行彈道仿真,生成的速度、彈道傾角、射程、彈道高作為辨識輸入。阻力系數初值使用零攻角和馬赫數進行插值計算,并將插值結果拉偏+30%作為辨識初值。升力系數的初始迭代值設為零。權重矩陣Rk設為:

表1 參數設置

3.2 仿真結果與分析

3.2.1辨識精度分析

根據以上設置的仿真條件進行參數辨識。辨識步長與每個辨識步長對應的迭代次數如圖2所示?？偙孀R步長為47步。每個步長內迭代次數在4～6次,且迭代5次的情況占多數?？梢娝惴ㄔ趶椀廊尉哂休^好的收斂性與一致性。

圖2 辨識步長及每個步長內的迭代次數

辨識后的彈道狀態量與實際彈道狀態量對比如圖3—圖6所示。可見辨識后的彈道速度、彈道傾角、射距、彈道高與真實實際值基本重合,收斂誤差小于指定精度要求。

圖3 真實速度與參數辨識后的速度對比

圖4 真實彈道傾角與參數辨識后的彈道傾角對比

圖5 真實射距與參數辨識后的射距對比

圖6 真實射高與參數辨識后的射高對比

阻力系數與升力系數辨識結果如圖7、圖8所示。

圖7 阻力系數(真實、初始值以及辨識值)對比

圖8 升力系數(真實、初始值以及辨識值)對比

由圖7、圖8可見,在初始阻力系數、升力系數存在誤差的情況下,經過迭代,辨識結果可以較好地收斂至真值。辨識誤差如圖9所示,無控段阻力系數誤差除,在開舵等個別點達到2%,其余絕大部分飛行段誤差均在1%以內。有控段阻力系數誤差在0.5%以內甚至更小。升力系數誤差基本在0.5%～1%以內。

圖9 阻力系數、升力系數辨識誤差

3.2.2算法在線計算能力分析

從算法復雜度的角度來看,假設每個辨識步長內算法迭代次數為n,每個辨識步長內的節點數記為N,則本文中算法的復雜度為O(n·N)。由于MPSP算法收斂性較快,通過幾次迭代一般即可滿足要求,因此迭代次數n通常較小。本文中的仿真案例對阻力系數初始值進行了30%的拉偏,升力系數初始值直接賦值為零,此時迭代次數為4～6次?？紤]實際情況氣動參數初值誤差一般小于文中的仿真設置,因此迭代次數會更少,算法復雜度可降至O(N),具備在線計算的潛力。

從算法耗時的角度來看,在CPU為AMD Ryzen 7 4800H@2.9Hz、內存為16G的計算機上,采用Matlab 2019a編制程序對本文中算法的耗時進行了統計。結果如圖10所示。橫軸表示全彈道的辨識步長序號,總共47步。其中前10步為無控段辨識,每個步長為4 s,平均耗時約0.34 s;后37步為有控段辨識,每個步長為2 s,平均耗時約0.17 s?？梢娝惴ê臅r較低。

圖10 每個步長的耗時

4 結論

本文中提出的基于MPSP算法的氣動辨識方法,是一種類比帶有終端約束制導律設計思想的一種辨識新思路。通過在某型155mm制導炮彈的升阻力系數的辨識應用中證實了其有效性,并有如下結論:

1) 在僅依賴輸入外彈道數據及彈體質量等必要信息的前提下,本文中算法可對升、阻力系數進行辨識,且辨識精度較高,阻力系數辨識誤差小于2%,升力系數辨識誤差小于1%。

2) 本文算法對待辨識氣動參數初值精度要求不高。在阻力系數拉偏30%、升力系數未知的情況下,通過數次迭代即可收斂至真值附近。

3) 本文辨識算法設計簡單、運算速度快,有利于工程實現。