陳健偉,王良明,李子杰,王 琦
(1.南京理工大學(xué) 能源與動(dòng)力工程學(xué)院,江蘇 南京 210094;2.北京遙感設(shè)備研究所,北京 100854)
相關(guān)研究表明,氣象因素在影響火炮射擊精度的所有因素中占70%以上[1-2],其中,風(fēng)是影響彈箭飛行的最主要?dú)庀笠蛩?因此,風(fēng)對(duì)彈箭飛行的影響一直是彈道學(xué)中的研究熱點(diǎn)之一。文獻(xiàn)[3-4]基于隨機(jī)風(fēng)場(chǎng)的統(tǒng)計(jì)特性,建立了空間簡(jiǎn)化風(fēng)速模型;文獻(xiàn)[5]基于火箭探空數(shù)據(jù)和熱成風(fēng)公式,建立了100 km以?xún)?nèi)大氣風(fēng)場(chǎng)分布模型,并分析了高空風(fēng)場(chǎng)對(duì)遠(yuǎn)程彈箭彈道的影響;文獻(xiàn)[6-7]根據(jù)外彈道氣象學(xué)理論,分析了幾種不同類(lèi)型風(fēng)場(chǎng)對(duì)彈丸彈道特性的影響,可以看出,目前對(duì)于彈道風(fēng)場(chǎng)的建模主要以簡(jiǎn)化平均風(fēng)模型和隨機(jī)擾動(dòng)風(fēng)模型為主,并沒(méi)有考慮實(shí)際情況下常見(jiàn)的典型風(fēng)場(chǎng)現(xiàn)象。
低空風(fēng)切變是指海拔600 m高度內(nèi)風(fēng)速與風(fēng)向的突變[8],其形成與天氣狀況、地形環(huán)境有著密切的關(guān)系。微下?lián)舯┝魇堑涂诊L(fēng)切變的一種,其形成概率高、影響范圍大,當(dāng)飛機(jī)彈箭在飛行過(guò)程中遭遇風(fēng)切變時(shí),突然變化的氣流速度和方向會(huì)改變其受力狀態(tài),從而影響其飛行軌跡和穩(wěn)定性[9-11]。目前,對(duì)于微下?lián)舯┝黠L(fēng)場(chǎng)的研究主要圍繞飛機(jī)展開(kāi),如飛行軌跡控制[12]、探測(cè)預(yù)警[13]等,而對(duì)于微下?lián)舯┝鲗?duì)外彈道過(guò)程的影響鮮有研究。基于此,本文基于流體力學(xué)基本原理,建立微下?lián)舯┝黠L(fēng)場(chǎng)的非對(duì)稱(chēng)多渦環(huán)模型,并與六自由度剛體彈道模型相結(jié)合,通過(guò)算例仿真,研究了微下?lián)舯┝黠L(fēng)場(chǎng)模型在彈道仿真中的應(yīng)用,以期為惡劣氣候環(huán)境下的外彈道研究提供思路與參考。
在空間范圍內(nèi),微下?lián)舯┝魍ǔ1憩F(xiàn)為如下的流動(dòng)形式:高度300~800 m范圍內(nèi),高壓氣旋(云層)中部形成一股向下的氣流,氣流撞擊地面后,呈輻射狀水平散開(kāi),并在外圍形成尾流旋渦,實(shí)景圖和示意圖分別如圖1和圖2所示。
圖1 微下?lián)舯┝鲗?shí)景圖
圖2 微下?lián)舯┝餍纬墒疽鈭D
在地面坐標(biāo)系Oxyz中,以地面上方一點(diǎn)OP為中心設(shè)置主渦環(huán),其半徑為R,如圖3所示,其在地面系內(nèi)的曲線方程可表示為
(1)
圖3 渦環(huán)設(shè)置示意圖
圖3中,下標(biāo)P表示主渦環(huán),I表示鏡像渦環(huán)。根據(jù)流體力學(xué)基本原理可知,主渦環(huán)流線方程為
(2)
式中:Γ為渦環(huán)強(qiáng)度,其大小取決于由渦環(huán)中心向下的誘導(dǎo)速度vz(0)及渦環(huán)半徑R。
Γ=2Rvz(0)
(3)
對(duì)于地面坐標(biāo)系內(nèi)的任意點(diǎn)OM(xM,yM,zM),rmax為該點(diǎn)到主渦環(huán)的最大距離,rmin為該點(diǎn)到主渦環(huán)的最小距離,F(k)為橢圓積分函數(shù),有:
(4)
在0≤k≤1的情況下,F(k)有近似表達(dá)式:
(5)
通過(guò)對(duì)渦環(huán)流線方程進(jìn)行偏導(dǎo)計(jì)算,可得到主渦環(huán)徑向誘導(dǎo)速度和軸向誘導(dǎo)速度:
(6)
(7)
實(shí)際風(fēng)場(chǎng)環(huán)境中,微下?lián)舯┝髦行南蛳碌臍饬髯矒舻孛嫦蛩闹苌㈤_(kāi)形成外流時(shí),地面處的風(fēng)速垂直分量應(yīng)該為0,因此,通過(guò)設(shè)置一個(gè)與主渦環(huán)參數(shù)相同,位置與主渦環(huán)關(guān)于水平面對(duì)稱(chēng)的鏡像渦環(huán)來(lái)滿(mǎn)足這一條件,鏡像渦環(huán)中心為OI(xP,yP,-zP),且有ψI=-ψP,因此,主渦環(huán)與鏡像渦環(huán)的誘導(dǎo)速度正好在水平面處相互抵消為0。通過(guò)與主渦環(huán)類(lèi)似的計(jì)算方法,求得鏡像渦環(huán)在地面坐標(biāo)系下的速度分量,最后,將兩渦環(huán)產(chǎn)生的誘導(dǎo)速度分量疊加,得到點(diǎn)OM處的合速度為
(8)
點(diǎn)OM的在地面系內(nèi)的流線方程為
(9)
對(duì)于2種特殊情況:
①rP=0時(shí),由流線方程計(jì)算出的誘導(dǎo)速度趨近正無(wú)窮,不合理,因此該點(diǎn)速度由渦環(huán)位函數(shù)偏導(dǎo)計(jì)算:
(10)
②渦環(huán)線上,rP=R,zM=zP處的點(diǎn),以渦環(huán)線為中心,構(gòu)造半徑為r的封閉環(huán)形圓柱,由于渦量在圓柱內(nèi)均勻分布,因此圓柱內(nèi)沿半徑方向誘導(dǎo)速度線性分布,如圖4所示。
圖4 渦核內(nèi)風(fēng)速分布示意圖
當(dāng)點(diǎn)OM位于環(huán)形圓柱內(nèi)時(shí),由渦環(huán)曲線方程聯(lián)立直線OPOM在水平面內(nèi)的投影直線方程,可計(jì)算出點(diǎn)M的坐標(biāo)和點(diǎn)N的坐標(biāo),而后根據(jù)點(diǎn)N的坐標(biāo),結(jié)合式(6)和式(7)即可獲得此時(shí)點(diǎn)OM處的風(fēng)矢量。
實(shí)際情況下,微下?lián)舯┝黠L(fēng)場(chǎng)往往不是只由單個(gè)風(fēng)暴中心(渦環(huán))組成的,且水平輻射的外流也是非對(duì)稱(chēng)的,甚至可能還存在上升氣流,因此,引入傾斜變換矩陣,通過(guò)多渦環(huán)疊加來(lái)模擬區(qū)域范圍內(nèi)的微下?lián)舯┝黠L(fēng)場(chǎng)。
如圖5所示,令主渦環(huán)坐標(biāo)系與地面坐標(biāo)系的夾角向量為(φθ0)T,對(duì)應(yīng)的旋轉(zhuǎn)矩陣為M(φ),則鏡像渦環(huán)與地面系夾角為(-φ-θ0)T,對(duì)應(yīng)的旋轉(zhuǎn)矩陣為M(θ)。通過(guò)矩陣變換M(φ)M(θ)完成主渦環(huán)坐標(biāo)系到地面坐標(biāo)系的轉(zhuǎn)換,類(lèi)似地,由矩陣變換M(-φ)M(-θ)實(shí)現(xiàn)鏡像渦環(huán)坐標(biāo)系到地面坐標(biāo)系的轉(zhuǎn)換。其中:
(11)
(12)
結(jié)合式(8),最終地面坐標(biāo)系下三維風(fēng)矢量為
(13)
圖5 傾斜渦環(huán)示意圖
美國(guó)JAWS(聯(lián)合機(jī)場(chǎng)天氣計(jì)劃)美聯(lián)邦空管署FAA(federal aviation administration,FAA)的氣象研究資料[14]顯示,復(fù)雜微下?lián)舯┝鳝h(huán)境時(shí),多風(fēng)暴核在空間的位置分布間隔較為接近,且往往同時(shí)存在上升氣流中心和下沉氣流中心。參考文獻(xiàn)[15]的飛行仿真實(shí)驗(yàn)數(shù)據(jù),進(jìn)行多渦環(huán)模型參數(shù)設(shè)置。
在地面坐標(biāo)系一定區(qū)域內(nèi)設(shè)定3個(gè)渦環(huán),風(fēng)場(chǎng)模型主要參數(shù)如表1所示,其中渦環(huán)1、渦環(huán)2為下沖氣流,渦環(huán)3為上升氣流;除渦環(huán)中心坐標(biāo)外,渦環(huán)2和渦環(huán)3的其余模型參數(shù)同渦環(huán)1。表1中,P1,P2,P3分別為渦環(huán)1、渦環(huán)2、渦環(huán)3在地面系下的中心點(diǎn)。
表1 多渦環(huán)風(fēng)場(chǎng)模型參數(shù)
圖6和圖7分別給出了微下?lián)舯┝鞫鄿u環(huán)模型計(jì)算得到的水平截面和垂直截面內(nèi)風(fēng)矢量分布情況。由圖6可以看出,氣流在每個(gè)渦環(huán)的中心位置向四周呈輻射狀流動(dòng),風(fēng)速沿徑向逐漸衰減;由圖7可以看出,渦環(huán)中心附近的下沖氣流具有較高的風(fēng)速,氣流運(yùn)動(dòng)較為復(fù)雜,且由于非對(duì)稱(chēng)傾斜角的存在,風(fēng)速向一側(cè)偏移,到達(dá)地面后產(chǎn)生回流。
圖6 z=100 m時(shí)水平面風(fēng)矢量分布
圖7 y=0時(shí)垂直面風(fēng)矢量分布
綜合圖6和圖7可以看出,所建立的非對(duì)稱(chēng)多渦環(huán)微下?lián)舯┝黠L(fēng)場(chǎng)模型具有良好的空間三維特性,能夠較為合理有效地描述實(shí)際微下?lián)舯┝黠L(fēng)場(chǎng)的風(fēng)速分布特性。
根據(jù)外彈道學(xué)理論[16],將運(yùn)動(dòng)中的火箭彈視為常質(zhì)量的自由剛體進(jìn)行受力分析,忽略火箭燃?xì)鈶T性作用力和力矩,同時(shí)不考慮發(fā)動(dòng)機(jī)燃料燃燒對(duì)質(zhì)心移動(dòng)加速度以及轉(zhuǎn)動(dòng)慣量的影響。通過(guò)聯(lián)立彈道坐標(biāo)系及第一彈軸坐標(biāo)系內(nèi)火箭彈的質(zhì)心運(yùn)動(dòng)和繞心運(yùn)動(dòng)方程,即可得到野戰(zhàn)火箭彈六自由度剛體彈道模型:
(14)
(15)
(16)
(17)
dm/dt=mb
(18)
式中:t為彈丸飛行時(shí)間;v為彈丸質(zhì)心速度大小;m為彈丸質(zhì)量;θ1和ψ2分別為彈丸速度高低角和速度方位角;Fx2、Fy2、Fz2為彈道坐標(biāo)系下作用在彈丸質(zhì)心的合力分量;X,Y,Z分別為彈丸在地面坐標(biāo)系下的位置分量;ωξ、ωη、ωζ分別為彈丸轉(zhuǎn)動(dòng)角速度在彈軸坐標(biāo)系下的分量;Mξ、Mη、Mζ分別為彈丸所受合外力矩在彈軸坐標(biāo)系下的分量;φa和φ2分別為彈軸高低角和彈軸方位角;JC和JA分別為彈丸極轉(zhuǎn)動(dòng)慣量和赤道轉(zhuǎn)動(dòng)慣量;γ為彈丸自轉(zhuǎn)角;mb為火箭發(fā)動(dòng)機(jī)燃料質(zhì)量燃燒速率,當(dāng)獲得火箭彈彈體參數(shù)、發(fā)射條件、氣象條件后,通過(guò)數(shù)值積分算法可計(jì)算野戰(zhàn)火箭彈的完整彈道??紤]風(fēng)的影響時(shí),速度坐標(biāo)系下作用在彈上的力可以表示為彈體相對(duì)速度的函數(shù):
F2=f(v-ω2)
(19)
式中:ω2為速度坐標(biāo)系內(nèi)的風(fēng)矢量,當(dāng)獲得地面坐標(biāo)系下風(fēng)速分量時(shí),可將該分量投影至速度坐標(biāo)系進(jìn)行彈道解算。
基于外彈道學(xué)基本理論,將地面系下的三維風(fēng)矢量計(jì)入六自由度剛體彈道模型進(jìn)行計(jì)算,通過(guò)數(shù)值仿真來(lái)驗(yàn)證風(fēng)場(chǎng)模型計(jì)算效果。
以某122 mm口徑尾翼火箭彈為例,其彈體參數(shù)和發(fā)射條件如表2所示。表中,d為彈徑,l為彈長(zhǎng),v0為火箭彈發(fā)射初速,ε為射角,φ為射向。在主動(dòng)段彈道計(jì)算時(shí),認(rèn)為發(fā)動(dòng)機(jī)燃料質(zhì)量勻速率減小,且轉(zhuǎn)動(dòng)慣量和質(zhì)心位置均相應(yīng)隨時(shí)間勻速率變化。仿真中,除微下?lián)舯┝黠L(fēng)場(chǎng)參數(shù)外,其余氣象條件以炮兵標(biāo)準(zhǔn)氣象條件為準(zhǔn)。以火箭彈發(fā)射位置為原點(diǎn),射向?yàn)镺x軸建立地面坐標(biāo)系,微下?lián)舯┝鳒u環(huán)中心位置同表1。
表2 火箭彈彈體參數(shù)與發(fā)射條件
根據(jù)風(fēng)切變的定義,風(fēng)速大小是風(fēng)切變強(qiáng)度的決定性因素。表3給出了JAWS的多普勒雷達(dá)測(cè)量統(tǒng)計(jì)結(jié)果和美軍飛行數(shù)據(jù)記錄器FDR提取的微下?lián)舯┝黠L(fēng)切變數(shù)據(jù)比較[14]。
由表3微下?lián)舯┝黠L(fēng)切變中心風(fēng)速變化換算得,渦環(huán)中心處的風(fēng)速變化范圍為7.8~27.2 m/s,因此,在該范圍內(nèi)設(shè)置微下?lián)舯┝鞫鄿u環(huán)模型的中心誘導(dǎo)風(fēng)速分別為10 m/s,15 m/s,20 m/s,25 m/s,進(jìn)行風(fēng)場(chǎng)干擾下的彈道仿真,并與彈道無(wú)風(fēng)情況下進(jìn)行比較。表4給出了彈道特征參數(shù)的仿真結(jié)果,表中,vT為火箭彈末速度。圖8為三維彈道軌跡圖。圖9和圖10分別為不同中心誘導(dǎo)風(fēng)速下火箭彈高低攻角δ1和方向攻角δ2隨時(shí)間的變化曲線。
表3 微下?lián)舯┝黠L(fēng)切變多普勒雷達(dá)觀測(cè)數(shù)據(jù)和FDR記錄數(shù)據(jù)比較
表4 彈道特征參數(shù)仿真結(jié)果
由表4和圖8可以看出,受到微下?lián)舯┝黠L(fēng)切變風(fēng)場(chǎng)的影響,火箭彈的彈道特征量發(fā)生了較大變化,具體表現(xiàn)在:與無(wú)風(fēng)情況下相比,火箭彈在空中飛行時(shí)間變長(zhǎng),落點(diǎn)側(cè)偏增大,射程減小,彈丸末速度降低;隨著渦環(huán)中心誘導(dǎo)風(fēng)速的增加,相關(guān)彈道特征量的變化愈加明顯,可見(jiàn)風(fēng)場(chǎng)強(qiáng)度是決定風(fēng)切變對(duì)彈箭飛行過(guò)程影響程度的重要因素。
圖8 彈道軌跡圖
圖9 高低攻角曲線
圖10 方向攻角曲線
由圖9和圖10可以看出,在主動(dòng)段,尾翼火箭彈初速較低,抗擾動(dòng)能力弱,當(dāng)穿越氣流變化復(fù)雜的風(fēng)場(chǎng)影響區(qū)時(shí),高低攻角和方向攻角都發(fā)生了較大波動(dòng),高低方向的攻角波動(dòng)出現(xiàn)了多個(gè)峰值,是因?yàn)槭艿搅硕鄠€(gè)微下?lián)舯┝鳒u環(huán)中心垂直氣流的影響。在推力作用下,火箭彈很快穿越高度幾百米范圍之內(nèi)的風(fēng)切變影響區(qū)域,而后在彈體靜穩(wěn)定力矩作用下,攻角運(yùn)動(dòng)逐漸收斂趨于平衡。
基于渦環(huán)法建立了非對(duì)稱(chēng)多核微下?lián)舯┝黠L(fēng)切變風(fēng)場(chǎng)模型,通過(guò)彈道仿真計(jì)算,得到如下結(jié)論:①非對(duì)稱(chēng)多渦環(huán)微下?lián)舯┝黠L(fēng)場(chǎng)模型具有良好的空間三維特性,能夠在一定程度上較好地描述微下?lián)舯┝黠L(fēng)切變的風(fēng)場(chǎng)流動(dòng)特征;②低初速?gòu)椉?如部分火箭彈)在飛行過(guò)程中受到微下?lián)舯┝黠L(fēng)場(chǎng)干擾,飛行時(shí)間、末速度、落點(diǎn)位置和角運(yùn)動(dòng)過(guò)程均會(huì)發(fā)生變化,從而產(chǎn)生射擊精度偏差。如何探究風(fēng)切變環(huán)境下彈箭飛行的控制方法,是后續(xù)進(jìn)一步研究的內(nèi)容。