劉 洋,劉 丹,常思江

(1.南京理工大學(xué) 能源與動力工程學(xué)院,江蘇 南京 210094; 2.西北工業(yè)集團(tuán)有限公司 設(shè)計(jì)二所,陜西 西安 710043)

彈丸氣動力系數(shù)是外彈道仿真、彈道設(shè)計(jì)及射表編制等的基礎(chǔ)數(shù)據(jù),如何獲取較高精度的彈丸氣動力系數(shù)是外彈道學(xué)、彈箭空氣動力學(xué)研究人員長期以來一直關(guān)注的研究主題。獲取彈丸氣動力系數(shù)的方式主要包括3種[1]:理論計(jì)算、風(fēng)洞實(shí)驗(yàn)數(shù)據(jù)擬合、根據(jù)飛行試驗(yàn)數(shù)據(jù)辨識。一般情況下,與理論計(jì)算和風(fēng)洞實(shí)驗(yàn)結(jié)果相比,根據(jù)飛行試驗(yàn)數(shù)據(jù)辨識的氣動力系數(shù)精度最高,與實(shí)際狀況最為接近[2]。

國內(nèi)外目前在氣動力系數(shù)辨識領(lǐng)域已開展了較多研究,使用的方法也是多種多樣。文獻(xiàn)[3-4]使用了遺傳算法及自適應(yīng)混沌變異粒子群算法對旋轉(zhuǎn)穩(wěn)定彈進(jìn)行了氣動系數(shù)辨識;文獻(xiàn)[5]提出一種采用雙加速度計(jì)測量的辨識方法,辨識出了自振角頻率、攻角和升力系數(shù)等。文獻(xiàn)[6]中提出了一種基于修正質(zhì)點(diǎn)彈道模型的簡化模型,相比于以往的辨識方法,減小了彈道方程的復(fù)雜程度,并使用了最小二乘法進(jìn)行辨識,但是辨識結(jié)果在精度上仍有提升的空間。

本文借鑒文獻(xiàn)[6]的思路,針對最小二乘法辨識精度較低的問題,提出采用一種改進(jìn)粒子群算法構(gòu)建氣動力系數(shù)辨識的模型,利用彈丸飛行的速度和轉(zhuǎn)速數(shù)據(jù),對旋轉(zhuǎn)穩(wěn)定彈的零升阻力系數(shù)、升力系數(shù)導(dǎo)數(shù)、馬格努斯力系數(shù)導(dǎo)數(shù)及極阻尼力矩系數(shù)導(dǎo)數(shù)進(jìn)行辨識,以期為旋轉(zhuǎn)彈丸氣動力系數(shù)辨識提供一種新的思路和模型。

1 系統(tǒng)模型

根據(jù)外彈道理論,六自由度剛體彈道模型可以全面地描述彈丸的飛行[7],但是該六自由度模型較為復(fù)雜,需要考慮的因素過多[8],從精度和計(jì)算量角度考慮,本文擬使用四自由度修正質(zhì)點(diǎn)彈道模型,并針對該模型做出以下假設(shè):

①彈丸關(guān)于彈軸對稱且飛行過程中無形變;

②彈丸飛行過程中保持穩(wěn)定;

③各個氣動力系數(shù)對攻角的導(dǎo)數(shù)可以認(rèn)為是線性關(guān)系,即:

C=Cαα

(1)

式中:C為各氣動力系數(shù),α為攻角,Cα為各氣動力系數(shù)相對于攻角的導(dǎo)數(shù)。

④忽略對飛行影響較小的力和力矩;

(2)

(3)

(4)

(5)

(6)

(7)

式中:ρ為大氣密度;D為彈丸直徑;l為彈丸特征長度;m為彈丸質(zhì)量;g為重力加速度;Fc為科氏力;u為彈丸的絕對速度矢量;ux、uy和uz為彈丸的絕對速度u在x,y,z方向上的分量;v為彈丸相對于風(fēng)的速度矢量;|v|為相對速度矢量v的模長;I1為極轉(zhuǎn)動慣量;ωc為彈丸的轉(zhuǎn)速;αe為彈丸的攻角;Cx0為彈丸的零升阻力系數(shù);CL,α為彈丸的升力系數(shù)導(dǎo)數(shù),Cm,α為彈丸的靜力矩系數(shù)導(dǎo)數(shù),Cz,α為彈丸的馬格努斯力系數(shù)導(dǎo)數(shù),CMxz,ω為彈丸的極阻尼力矩系數(shù)導(dǎo)數(shù)。

由粒子群算法的特性可知[9],維度數(shù)量直接影響著適應(yīng)度函數(shù)的收斂速度?，F(xiàn)有研究表明[10],利用彈丸紙靶章動試驗(yàn)可單獨(dú)獲得高精度的靜力矩系數(shù)導(dǎo)數(shù),因此,為了提高氣動系數(shù)辨識的精度和效率,本文建立的氣動系數(shù)辨識模型中將靜力矩系數(shù)作為已知量,故可對上述四自由度修正質(zhì)點(diǎn)彈道方程做出一定的簡化。

將攻角方程(7)代入速度方程(2)中,得到簡化后的速度方程(8)為

(8)

從上式可以看出,適應(yīng)度函數(shù)就減少了攻角所帶來的3個維度,在精度不變的基礎(chǔ)上有可能提高收斂的速度。

2 粒子群算法原理與改進(jìn)

2.1 粒子群算法原理

粒子群算法(簡稱PSO算法),與蟻群算法、遺傳算法和神經(jīng)網(wǎng)絡(luò)算法等稱為智能算法。與傳統(tǒng)算法不同,這種智能算法在全局搜索上具有較強(qiáng)的能力,而且這類算法不依賴于模型本身的好壞,其目的在于在解空間中找出使得適應(yīng)度函數(shù)最大或最小的解。這與Chapman-Kirk法[10]類似,但是不會出現(xiàn)Chapman-Kirk法中偶爾發(fā)生的正規(guī)方程系數(shù)矩陣成為準(zhǔn)奇異矩陣的現(xiàn)象。

PSO算法主要分為設(shè)計(jì)變量、適應(yīng)度函數(shù)和數(shù)學(xué)模型等幾個部分,其計(jì)算原理如圖1所示。

①設(shè)計(jì)變量θ。 本文中的設(shè)計(jì)變量即為待辨識的各個系數(shù)。即:

②適應(yīng)度函數(shù)J。 常用的適應(yīng)度函數(shù)有最小二乘準(zhǔn)則、極大似然準(zhǔn)則等。這里選用準(zhǔn)則如下式:

式中:Ui為當(dāng)前粒子的計(jì)算值,Ut0為時間區(qū)間[t0,t1]上t0時刻的測量值,Ut1為t1時刻的測量值。(Ui-Ut1)和(Ut1-Ut0)為2個向量之間的標(biāo)準(zhǔn)歐幾里得距離。由于ux,uy、uz和ωc具有不同的數(shù)量級和量綱,標(biāo)準(zhǔn)歐幾里德距離可以修正數(shù)量級和量綱帶來的差異,使結(jié)果更加精確。

③數(shù)學(xué)模型。 數(shù)學(xué)模型即為簡化后的彈道方程,即式(3)～式(6)和式(8)。

④終止條件設(shè)置。

利用轉(zhuǎn)速辨識:J<0.000 01或者k>200。

利用速度辨識:J<0.000 4或者k>1 000。

⑤粒子群算法主要參數(shù)。 速度vp與位置n主要由下式?jīng)Q定。

其中初始速度與位置在解空間內(nèi)隨機(jī)產(chǎn)生。式中:kn為迭代次數(shù),w為慣性權(quán)重,c1和c2為學(xué)習(xí)因子,gi為當(dāng)前粒子個體最優(yōu)點(diǎn),p為全局最優(yōu)點(diǎn),ξ1和ξ2為[-1,1]之間的一個隨機(jī)數(shù)。

2.2 算法的改進(jìn)

標(biāo)準(zhǔn)PSO算法多用于規(guī)劃、金融、設(shè)計(jì)等多目標(biāo)領(lǐng)域,算法本身具有運(yùn)算簡單、容易實(shí)現(xiàn)等優(yōu)點(diǎn)。但是,在參數(shù)辨識領(lǐng)域,單純的標(biāo)準(zhǔn)PSO算法會有收斂速度慢、迭代方向性差、收斂步數(shù)過長等缺點(diǎn)。因此,針對彈丸氣動力系數(shù)辨識的特點(diǎn),對PSO算法進(jìn)行一定的改進(jìn)與優(yōu)化。

①出界反射。 一般對待辨識的系數(shù)進(jìn)行范圍上的估計(jì),如阻力系數(shù)的取值一般在[0.1,0.5]內(nèi),當(dāng)某個粒子超出該范圍時,需要對該粒子的位置進(jìn)行處理。本文使用出界反射原理,即當(dāng)某一維度超出范圍時,將該系數(shù)反射回來,這個反射僅在該維度上進(jìn)行,不涉及其他的維度。反射公式為

θk=2Mk-θk-1

式中:M為待辨識系數(shù)的范圍,即PSO算法的搜索區(qū)間。

②分步辨識。 由于對轉(zhuǎn)速方程進(jìn)行了簡化,忽略了赤道阻尼力矩系數(shù),使得轉(zhuǎn)速ωc僅是極阻尼力矩系數(shù)導(dǎo)數(shù)CMxz,ω的函數(shù)。故可先用轉(zhuǎn)速方程(6)對極阻尼力矩系數(shù)導(dǎo)數(shù)進(jìn)行辨識,再將辨識后的極阻尼力矩系數(shù)導(dǎo)數(shù)代回方程(6)中,與其他方程一起構(gòu)成完整的彈道方程,對其余的3個系數(shù)進(jìn)行辨識。這樣的辨識順序會提升辨識的速度,提高辨識結(jié)果的精確度。

③變速度與區(qū)間收斂。 由PSO算法的特性可知,過大的速度有可能導(dǎo)致一次次地“越過”全局最優(yōu)解,所以隨著迭代次數(shù)的增多和適應(yīng)度函數(shù)的收斂,速度應(yīng)逐漸減小。一般情況下,可以使用自適應(yīng)學(xué)習(xí)因子和慣性權(quán)重來放慢速度。而在本文中,因總迭代步數(shù)較少,故可依據(jù)適應(yīng)度函數(shù)的值選擇學(xué)習(xí)因子和慣性權(quán)重。同時由于PSO算法是全局的優(yōu)化算法,其搜索的區(qū)間較大,若在迭代過程中隨著全局最優(yōu)解p的變化將區(qū)間半徑R逐漸變小,則能夠有效地提高收斂速度,區(qū)間中心為當(dāng)前全局最優(yōu)解。本文中學(xué)習(xí)因子、慣性權(quán)重的變化和區(qū)間設(shè)置如表1。

表1 PSO參數(shù)設(shè)置

④初始區(qū)間的選擇。 由于彈丸的飛行過程是一個連續(xù)的過程,其氣動系數(shù)的變化也是連續(xù)的,在計(jì)算下一個區(qū)間時,可以以上一個區(qū)間的解θk-1為中心,區(qū)間半徑長度R為

將初始區(qū)間Mk取為

Mk=[θk-1-Rθk-1+R]

β的取值與整體數(shù)據(jù)密度有關(guān),密度越大,β越小;密度越小,β越大。本文中選取β=10。

⑤粒子變異。 為避免陷入局部收斂,根據(jù)迭代的步數(shù)與當(dāng)前J的值對粒子進(jìn)行變異。由于氣動系數(shù)辨識的解空間中存在唯一解,或在某些條件下(如測量不準(zhǔn)確)僅存在少數(shù)解,所以本文的變異選擇在區(qū)間內(nèi)為重置當(dāng)前粒子的位置,即隨機(jī)產(chǎn)生位置n與速度vp。

2.3 算法速度比較

使用改進(jìn)PSO算法與標(biāo)準(zhǔn)PSO算法對同一段數(shù)據(jù)進(jìn)行辨識,其收斂速度如圖2所示,終止條件為適應(yīng)度函數(shù)J<0.000 4。為了讓結(jié)果更加明顯,將J取對數(shù),即lgJ<-3.4。

可以看出,改進(jìn)PSO算法在迭代次數(shù)達(dá)到36次時就已經(jīng)達(dá)到了收斂要求,而標(biāo)準(zhǔn)PSO算法在迭代次數(shù)達(dá)到100次時還與迭代終止條件有著較大的差距?？紤]到所采用的坐標(biāo)是對數(shù)坐標(biāo),這個差距將使迭代次數(shù)成數(shù)十倍地提升。這說明改進(jìn)PSO算法的收斂速度要遠(yuǎn)遠(yuǎn)高于標(biāo)準(zhǔn)PSO算法。

3 數(shù)值驗(yàn)證

3.1 155 mm榴彈仿真數(shù)據(jù)驗(yàn)證

利用某155 mm榴彈的仿真彈道數(shù)據(jù)進(jìn)行氣動力系數(shù)辨識,并將辨識的結(jié)果與仿真使用的氣動力系數(shù)理論值進(jìn)行比較,以驗(yàn)證改進(jìn)PSO算法用于旋轉(zhuǎn)彈氣動系數(shù)辨識的可行性和有效性。在彈道仿真和氣動力系數(shù)辨識的過程中,氣象條件均選擇標(biāo)準(zhǔn)炮兵氣象條件,初速選擇930 m/s,射角ψ為35°,45°和55°,不考慮科氏力,重力加速度為9.81 m/s2。辨識模型的輸入量為:[t0,t1]時間區(qū)間上的3個速度ux,uy,uz和轉(zhuǎn)速ωc;辨識模型的輸入量為氣動力系數(shù)的辨識值。在實(shí)際應(yīng)用中,彈丸速度可由坐標(biāo)雷達(dá)測得,彈丸轉(zhuǎn)速可由彈載地磁傳感器測得。

將辨識后的結(jié)果與理論值對比,并繪制出各個系數(shù)隨馬赫數(shù)的變化曲線。

零升阻力系數(shù)的辨識結(jié)果與理論值對比如圖3所示?？梢钥闯?零升阻力系數(shù)的辨識結(jié)果與理論值基本一致,最大相對誤差約為0.1%,平均相對誤差為0.036%。

CL,α/Cm,α的辨識結(jié)果如圖4所示。其中,圖4中的辨識值乘以對應(yīng)馬赫數(shù)下的靜力矩系數(shù)為升力系數(shù)導(dǎo)數(shù)。平均誤差約為0.05%,最大誤差約為0.15%。多次辨識的誤差并不完全一致,這是因?yàn)镻SO算法位置區(qū)間的初始化與前一個點(diǎn)的辨識值有著直接的關(guān)系,辨識的結(jié)果在辨識過程中產(chǎn)生了一步步的偏移,不同的辨識過程將在總體上產(chǎn)生不同的偏移,但是都在一個較小的范圍內(nèi)。

極阻尼力矩系數(shù)導(dǎo)數(shù)的辨識結(jié)果如圖5所示。辨識結(jié)果的平均誤差為0.041%,最大誤差為0.06%。由于極阻尼力矩系數(shù)是由轉(zhuǎn)速方程(6)進(jìn)行辨識,沒有其他未知量的干擾,所以辨識的精度要比零升阻力系數(shù)和CL,α/Cm,α更高。

Cz,α/Cm,α的辨識結(jié)果如圖6所示。與其他3個系數(shù)的辨識結(jié)果不同,馬格努斯力系數(shù)導(dǎo)數(shù)的辨識結(jié)果誤差較大,但是在不同馬赫段上誤差是不同的?？傮w上,馬赫數(shù)在0.8～1.25范圍內(nèi),辨識值與理論值的相對誤差在30%～50%之間;馬赫數(shù)在1.25～2.7之間,辨識值的相對誤差很大,結(jié)果已不可用。此外,從圖6可知,不同射角下的誤差也有所不同,55°射角下的辨識值要好于其他射角下的辨識值,這可能是由于在55°射角下全彈道上的動力平衡角相對較大,有利于馬格努斯力系數(shù)導(dǎo)數(shù)的辨識。

以上辨識結(jié)果表明,本文所采用的辨識方法對零升阻力系數(shù)Cx0、升力系數(shù)導(dǎo)數(shù)CL,α和極阻尼力矩系數(shù)導(dǎo)數(shù)CMxz,ω具有較高的辨識精度,而馬格努斯力系數(shù)導(dǎo)數(shù)的可辨識性較低[11],辨識值與理論值相差較大。3個射角下Cx0,CL,α和CMxz,ω沒有明顯的差別,這說明發(fā)射條件的改變對改進(jìn)PSO算法在氣動力系數(shù)上的辨識沒有影響。

將辨識結(jié)果代回方程重新計(jì)算彈道，并與辨識所用的仿真值進(jìn)行對比，對比結(jié)果如圖7所示。

從圖7可以看出,由氣動系數(shù)辨識值計(jì)算的彈道軌跡和速度變化均與測量結(jié)果較為吻合,其中彈道落點(diǎn)的測量值為(26 386,0,918),而根據(jù)氣動系數(shù)辨識值計(jì)算的彈道落點(diǎn)是(26 378,0,918),落點(diǎn)差為8 m。這個結(jié)果說明了辨識的結(jié)果具有較高的準(zhǔn)確性。

為了更加貼近工程實(shí)際,將速度ux,uy同時加入相同的高斯白噪聲,共選擇2個工況:①期望E=0,標(biāo)準(zhǔn)差σ=0.01 m/s;②期望E=0,標(biāo)準(zhǔn)差σ=0.03 m/s。由于數(shù)據(jù)密度較密,數(shù)據(jù)初段2個采集點(diǎn)之間的差值在0.8 m/s左右,數(shù)據(jù)末端2個采集點(diǎn)之間的差值在0.04 m/s左右,而在實(shí)際工程運(yùn)用中,需要對數(shù)據(jù)進(jìn)行降噪處理,處理后的數(shù)據(jù)噪聲不會很大,所以這樣的高斯白噪聲是相對符合實(shí)際的。辨識出的零升阻力系數(shù)如圖8所示。

由圖8可以看出,該算法對帶噪聲的飛行數(shù)據(jù)具有較好的辨識能力,其中對噪聲較小數(shù)據(jù)的辨識精度要比噪聲較大數(shù)據(jù)的辨識精度要高。在高馬赫數(shù)下的辨識精度要比低馬赫數(shù)下的辨識精度要高,這是因?yàn)楦唏R赫數(shù)為彈道初段,速度變化較大,所以噪聲影響較小;而低馬赫數(shù)為彈道頂端與彈道末段,速度變化較小,所以噪聲影響較大。而由于算法中依區(qū)間收斂,在辨識過程中減小了噪聲對飛行數(shù)據(jù)的干擾。

4 結(jié)束語

本文以四自由度修正質(zhì)點(diǎn)彈道模型為基礎(chǔ),利用改進(jìn)粒子群算法對零升阻力系數(shù)、升力系數(shù)導(dǎo)數(shù)、極阻尼力矩系數(shù)導(dǎo)數(shù)和馬格努斯力系數(shù)導(dǎo)數(shù)進(jìn)行了辨識,取得了較好的結(jié)果,并得出以下結(jié)論:

①改進(jìn)粒子群算法比標(biāo)準(zhǔn)粒子群算法具有更快的收斂速度。

②通過該方法得到的零升阻力系數(shù)、升力系數(shù)導(dǎo)數(shù)和極阻尼力矩系數(shù)導(dǎo)數(shù)具有較高的精度,滿足實(shí)際工程的需要。

③對于帶有噪聲的數(shù)據(jù)能夠進(jìn)行有效地辨識,且能夠在一定程度上減小噪聲對辨識結(jié)果的影響。