劉 洋,劉 丹,常思江
(1.南京理工大學 能源與動力工程學院,江蘇 南京 210094; 2.西北工業(yè)集團有限公司 設計二所,陜西 西安 710043)
彈丸氣動力系數(shù)是外彈道仿真、彈道設計及射表編制等的基礎數(shù)據(jù),如何獲取較高精度的彈丸氣動力系數(shù)是外彈道學、彈箭空氣動力學研究人員長期以來一直關注的研究主題。獲取彈丸氣動力系數(shù)的方式主要包括3種[1]:理論計算、風洞實驗數(shù)據(jù)擬合、根據(jù)飛行試驗數(shù)據(jù)辨識。一般情況下,與理論計算和風洞實驗結果相比,根據(jù)飛行試驗數(shù)據(jù)辨識的氣動力系數(shù)精度最高,與實際狀況最為接近[2]。
國內外目前在氣動力系數(shù)辨識領域已開展了較多研究,使用的方法也是多種多樣。文獻[3-4]使用了遺傳算法及自適應混沌變異粒子群算法對旋轉穩(wěn)定彈進行了氣動系數(shù)辨識;文獻[5]提出一種采用雙加速度計測量的辨識方法,辨識出了自振角頻率、攻角和升力系數(shù)等。文獻[6]中提出了一種基于修正質點彈道模型的簡化模型,相比于以往的辨識方法,減小了彈道方程的復雜程度,并使用了最小二乘法進行辨識,但是辨識結果在精度上仍有提升的空間。
本文借鑒文獻[6]的思路,針對最小二乘法辨識精度較低的問題,提出采用一種改進粒子群算法構建氣動力系數(shù)辨識的模型,利用彈丸飛行的速度和轉速數(shù)據(jù),對旋轉穩(wěn)定彈的零升阻力系數(shù)、升力系數(shù)導數(shù)、馬格努斯力系數(shù)導數(shù)及極阻尼力矩系數(shù)導數(shù)進行辨識,以期為旋轉彈丸氣動力系數(shù)辨識提供一種新的思路和模型。
根據(jù)外彈道理論,六自由度剛體彈道模型可以全面地描述彈丸的飛行[7],但是該六自由度模型較為復雜,需要考慮的因素過多[8],從精度和計算量角度考慮,本文擬使用四自由度修正質點彈道模型,并針對該模型做出以下假設:
①彈丸關于彈軸對稱且飛行過程中無形變;
②彈丸飛行過程中保持穩(wěn)定;
③各個氣動力系數(shù)對攻角的導數(shù)可以認為是線性關系,即:
C=Cαα
(1)
式中:C為各氣動力系數(shù),α為攻角,Cα為各氣動力系數(shù)相對于攻角的導數(shù)。
④忽略對飛行影響較小的力和力矩;
(2)
(3)
(4)
(5)
(6)
(7)
式中:ρ為大氣密度;D為彈丸直徑;l為彈丸特征長度;m為彈丸質量;g為重力加速度;Fc為科氏力;u為彈丸的絕對速度矢量;ux、uy和uz為彈丸的絕對速度u在x,y,z方向上的分量;v為彈丸相對于風的速度矢量;|v|為相對速度矢量v的模長;I1為極轉動慣量;ωc為彈丸的轉速;αe為彈丸的攻角;Cx0為彈丸的零升阻力系數(shù);CL,α為彈丸的升力系數(shù)導數(shù),Cm,α為彈丸的靜力矩系數(shù)導數(shù),Cz,α為彈丸的馬格努斯力系數(shù)導數(shù),CMxz,ω為彈丸的極阻尼力矩系數(shù)導數(shù)。
由粒子群算法的特性可知[9],維度數(shù)量直接影響著適應度函數(shù)的收斂速度。現(xiàn)有研究表明[10],利用彈丸紙靶章動試驗可單獨獲得高精度的靜力矩系數(shù)導數(shù),因此,為了提高氣動系數(shù)辨識的精度和效率,本文建立的氣動系數(shù)辨識模型中將靜力矩系數(shù)作為已知量,故可對上述四自由度修正質點彈道方程做出一定的簡化。
將攻角方程(7)代入速度方程(2)中,得到簡化后的速度方程(8)為
(8)
從上式可以看出,適應度函數(shù)就減少了攻角所帶來的3個維度,在精度不變的基礎上有可能提高收斂的速度。
粒子群算法(簡稱PSO算法),與蟻群算法、遺傳算法和神經網絡算法等稱為智能算法。與傳統(tǒng)算法不同,這種智能算法在全局搜索上具有較強的能力,而且這類算法不依賴于模型本身的好壞,其目的在于在解空間中找出使得適應度函數(shù)最大或最小的解。這與Chapman-Kirk法[10]類似,但是不會出現(xiàn)Chapman-Kirk法中偶爾發(fā)生的正規(guī)方程系數(shù)矩陣成為準奇異矩陣的現(xiàn)象。
PSO算法主要分為設計變量、適應度函數(shù)和數(shù)學模型等幾個部分,其計算原理如圖1所示。
①設計變量θ。 本文中的設計變量即為待辨識的各個系數(shù)。即:
②適應度函數(shù)J。 常用的適應度函數(shù)有最小二乘準則、極大似然準則等。這里選用準則如下式:
式中:Ui為當前粒子的計算值,Ut0為時間區(qū)間[t0,t1]上t0時刻的測量值,Ut1為t1時刻的測量值。(Ui-Ut1)和(Ut1-Ut0)為2個向量之間的標準歐幾里得距離。由于ux,uy、uz和ωc具有不同的數(shù)量級和量綱,標準歐幾里德距離可以修正數(shù)量級和量綱帶來的差異,使結果更加精確。
③數(shù)學模型。 數(shù)學模型即為簡化后的彈道方程,即式(3)~式(6)和式(8)。
④終止條件設置。
利用轉速辨識:J<0.000 01或者k>200。
利用速度辨識:J<0.000 4或者k>1 000。
⑤粒子群算法主要參數(shù)。 速度vp與位置n主要由下式決定。
其中初始速度與位置在解空間內隨機產生。式中:kn為迭代次數(shù),w為慣性權重,c1和c2為學習因子,gi為當前粒子個體最優(yōu)點,p為全局最優(yōu)點,ξ1和ξ2為[-1,1]之間的一個隨機數(shù)。
標準PSO算法多用于規(guī)劃、金融、設計等多目標領域,算法本身具有運算簡單、容易實現(xiàn)等優(yōu)點。但是,在參數(shù)辨識領域,單純的標準PSO算法會有收斂速度慢、迭代方向性差、收斂步數(shù)過長等缺點。因此,針對彈丸氣動力系數(shù)辨識的特點,對PSO算法進行一定的改進與優(yōu)化。
①出界反射。 一般對待辨識的系數(shù)進行范圍上的估計,如阻力系數(shù)的取值一般在[0.1,0.5]內,當某個粒子超出該范圍時,需要對該粒子的位置進行處理。本文使用出界反射原理,即當某一維度超出范圍時,將該系數(shù)反射回來,這個反射僅在該維度上進行,不涉及其他的維度。反射公式為
θk=2Mk-θk-1
式中:M為待辨識系數(shù)的范圍,即PSO算法的搜索區(qū)間。
②分步辨識。 由于對轉速方程進行了簡化,忽略了赤道阻尼力矩系數(shù),使得轉速ωc僅是極阻尼力矩系數(shù)導數(shù)CMxz,ω的函數(shù)。故可先用轉速方程(6)對極阻尼力矩系數(shù)導數(shù)進行辨識,再將辨識后的極阻尼力矩系數(shù)導數(shù)代回方程(6)中,與其他方程一起構成完整的彈道方程,對其余的3個系數(shù)進行辨識。這樣的辨識順序會提升辨識的速度,提高辨識結果的精確度。
③變速度與區(qū)間收斂。 由PSO算法的特性可知,過大的速度有可能導致一次次地“越過”全局最優(yōu)解,所以隨著迭代次數(shù)的增多和適應度函數(shù)的收斂,速度應逐漸減小。一般情況下,可以使用自適應學習因子和慣性權重來放慢速度。而在本文中,因總迭代步數(shù)較少,故可依據(jù)適應度函數(shù)的值選擇學習因子和慣性權重。同時由于PSO算法是全局的優(yōu)化算法,其搜索的區(qū)間較大,若在迭代過程中隨著全局最優(yōu)解p的變化將區(qū)間半徑R逐漸變小,則能夠有效地提高收斂速度,區(qū)間中心為當前全局最優(yōu)解。本文中學習因子、慣性權重的變化和區(qū)間設置如表1。
表1 PSO參數(shù)設置
④初始區(qū)間的選擇。 由于彈丸的飛行過程是一個連續(xù)的過程,其氣動系數(shù)的變化也是連續(xù)的,在計算下一個區(qū)間時,可以以上一個區(qū)間的解θk-1為中心,區(qū)間半徑長度R為
將初始區(qū)間Mk取為
Mk=[θk-1-Rθk-1+R]
β的取值與整體數(shù)據(jù)密度有關,密度越大,β越小;密度越小,β越大。本文中選取β=10。
⑤粒子變異。 為避免陷入局部收斂,根據(jù)迭代的步數(shù)與當前J的值對粒子進行變異。由于氣動系數(shù)辨識的解空間中存在唯一解,或在某些條件下(如測量不準確)僅存在少數(shù)解,所以本文的變異選擇在區(qū)間內為重置當前粒子的位置,即隨機產生位置n與速度vp。
使用改進PSO算法與標準PSO算法對同一段數(shù)據(jù)進行辨識,其收斂速度如圖2所示,終止條件為適應度函數(shù)J<0.000 4。為了讓結果更加明顯,將J取對數(shù),即lgJ<-3.4。
可以看出,改進PSO算法在迭代次數(shù)達到36次時就已經達到了收斂要求,而標準PSO算法在迭代次數(shù)達到100次時還與迭代終止條件有著較大的差距??紤]到所采用的坐標是對數(shù)坐標,這個差距將使迭代次數(shù)成數(shù)十倍地提升。這說明改進PSO算法的收斂速度要遠遠高于標準PSO算法。
利用某155 mm榴彈的仿真彈道數(shù)據(jù)進行氣動力系數(shù)辨識,并將辨識的結果與仿真使用的氣動力系數(shù)理論值進行比較,以驗證改進PSO算法用于旋轉彈氣動系數(shù)辨識的可行性和有效性。在彈道仿真和氣動力系數(shù)辨識的過程中,氣象條件均選擇標準炮兵氣象條件,初速選擇930 m/s,射角ψ為35°,45°和55°,不考慮科氏力,重力加速度為9.81 m/s2。辨識模型的輸入量為:[t0,t1]時間區(qū)間上的3個速度ux,uy,uz和轉速ωc;辨識模型的輸入量為氣動力系數(shù)的辨識值。在實際應用中,彈丸速度可由坐標雷達測得,彈丸轉速可由彈載地磁傳感器測得。
將辨識后的結果與理論值對比,并繪制出各個系數(shù)隨馬赫數(shù)的變化曲線。
零升阻力系數(shù)的辨識結果與理論值對比如圖3所示??梢钥闯?零升阻力系數(shù)的辨識結果與理論值基本一致,最大相對誤差約為0.1%,平均相對誤差為0.036%。
CL,α/Cm,α的辨識結果如圖4所示。其中,圖4中的辨識值乘以對應馬赫數(shù)下的靜力矩系數(shù)為升力系數(shù)導數(shù)。平均誤差約為0.05%,最大誤差約為0.15%。多次辨識的誤差并不完全一致,這是因為PSO算法位置區(qū)間的初始化與前一個點的辨識值有著直接的關系,辨識的結果在辨識過程中產生了一步步的偏移,不同的辨識過程將在總體上產生不同的偏移,但是都在一個較小的范圍內。
極阻尼力矩系數(shù)導數(shù)的辨識結果如圖5所示。辨識結果的平均誤差為0.041%,最大誤差為0.06%。由于極阻尼力矩系數(shù)是由轉速方程(6)進行辨識,沒有其他未知量的干擾,所以辨識的精度要比零升阻力系數(shù)和CL,α/Cm,α更高。
Cz,α/Cm,α的辨識結果如圖6所示。與其他3個系數(shù)的辨識結果不同,馬格努斯力系數(shù)導數(shù)的辨識結果誤差較大,但是在不同馬赫段上誤差是不同的。總體上,馬赫數(shù)在0.8~1.25范圍內,辨識值與理論值的相對誤差在30%~50%之間;馬赫數(shù)在1.25~2.7之間,辨識值的相對誤差很大,結果已不可用。此外,從圖6可知,不同射角下的誤差也有所不同,55°射角下的辨識值要好于其他射角下的辨識值,這可能是由于在55°射角下全彈道上的動力平衡角相對較大,有利于馬格努斯力系數(shù)導數(shù)的辨識。
以上辨識結果表明,本文所采用的辨識方法對零升阻力系數(shù)Cx0、升力系數(shù)導數(shù)CL,α和極阻尼力矩系數(shù)導數(shù)CMxz,ω具有較高的辨識精度,而馬格努斯力系數(shù)導數(shù)的可辨識性較低[11],辨識值與理論值相差較大。3個射角下Cx0,CL,α和CMxz,ω沒有明顯的差別,這說明發(fā)射條件的改變對改進PSO算法在氣動力系數(shù)上的辨識沒有影響。
將辨識結果代回方程重新計算彈道,并與辨識所用的仿真值進行對比,對比結果如圖7所示。
從圖7可以看出,由氣動系數(shù)辨識值計算的彈道軌跡和速度變化均與測量結果較為吻合,其中彈道落點的測量值為(26 386,0,918),而根據(jù)氣動系數(shù)辨識值計算的彈道落點是(26 378,0,918),落點差為8 m。這個結果說明了辨識的結果具有較高的準確性。
為了更加貼近工程實際,將速度ux,uy同時加入相同的高斯白噪聲,共選擇2個工況:①期望E=0,標準差σ=0.01 m/s;②期望E=0,標準差σ=0.03 m/s。由于數(shù)據(jù)密度較密,數(shù)據(jù)初段2個采集點之間的差值在0.8 m/s左右,數(shù)據(jù)末端2個采集點之間的差值在0.04 m/s左右,而在實際工程運用中,需要對數(shù)據(jù)進行降噪處理,處理后的數(shù)據(jù)噪聲不會很大,所以這樣的高斯白噪聲是相對符合實際的。辨識出的零升阻力系數(shù)如圖8所示。
由圖8可以看出,該算法對帶噪聲的飛行數(shù)據(jù)具有較好的辨識能力,其中對噪聲較小數(shù)據(jù)的辨識精度要比噪聲較大數(shù)據(jù)的辨識精度要高。在高馬赫數(shù)下的辨識精度要比低馬赫數(shù)下的辨識精度要高,這是因為高馬赫數(shù)為彈道初段,速度變化較大,所以噪聲影響較小;而低馬赫數(shù)為彈道頂端與彈道末段,速度變化較小,所以噪聲影響較大。而由于算法中依區(qū)間收斂,在辨識過程中減小了噪聲對飛行數(shù)據(jù)的干擾。
本文以四自由度修正質點彈道模型為基礎,利用改進粒子群算法對零升阻力系數(shù)、升力系數(shù)導數(shù)、極阻尼力矩系數(shù)導數(shù)和馬格努斯力系數(shù)導數(shù)進行了辨識,取得了較好的結果,并得出以下結論:
①改進粒子群算法比標準粒子群算法具有更快的收斂速度。
②通過該方法得到的零升阻力系數(shù)、升力系數(shù)導數(shù)和極阻尼力矩系數(shù)導數(shù)具有較高的精度,滿足實際工程的需要。
③對于帶有噪聲的數(shù)據(jù)能夠進行有效地辨識,且能夠在一定程度上減小噪聲對辨識結果的影響。