張 超,王學德,王 勇,曹友銓
(南京理工大學 能源與動力工程學院,江蘇 南京 210094)
低速旋轉尾翼式彈箭氣動特性數值研究
張 超,王學德,王 勇,曹友銓
(南京理工大學 能源與動力工程學院,江蘇 南京 210094)
為了研究低速旋轉對尾翼式彈箭氣動特性的影響,采用三維非定常N-S方程并結合滑移網格技術,在小攻角和全馬赫數下,對某尾翼彈在低轉速狀態下的繞流流場進行了數值模擬。以美國陸-海軍動導數計算標模驗證該文算法的有效性。結果表明該方法有較高的精確度。由不同馬赫數、轉速和滾轉角條件下的計算結果發現:縱向氣動特性(即升力、阻力、俯仰力矩)不隨轉速而變化,平均滾轉力矩系數和轉速為定比例關系,平均馬格努斯力系數隨轉速呈非線性變化,瞬時馬格努斯力系數隨滾轉角呈正弦變化。
彈箭;馬格努斯效應;滑移網格;氣動特性;數值計算
彈箭通常采用繞彈體縱軸旋轉的飛行方式以提高射擊密集度以及飛行穩定性,然而,旋轉的同時也會帶來諸如馬格努斯效應、陀螺效應等問題。相比于無旋彈箭,旋轉彈箭的空氣動力特性非常復雜。因此,研究旋轉彈箭相關氣動特性對彈箭設計,尤其是無控彈箭[1],有著極其重要的影響。
在21世紀以前,人們主要通過實驗和工程估算方法研究自旋彈箭的氣動特性。伴隨著計算機技術性能的大幅度提高以及相關數值計算軟件的日臻完善,加之具有計算周期短、精度高、成本低廉等優勢,這使得利用成熟的計算軟件對彈箭進行數值仿真成為主流趨勢。美國航天宇航局[2]分別采用風洞實驗、工程估算和數值計算3種方式對自旋彈丸的氣動特性進行了專門的研究。Sturek W B與Pechier M[3-4]對尖拱-圓柱型旋轉彈丸同時進行了風洞試驗和數值仿真,并將兩者的結果做了對比,發現兩者吻合程度較高,尤其是在小攻角范圍內。在國內,吳甲生[5]對旋轉彈箭空氣動力學效應進行了系統深入的分析,并指出旋轉彈箭未來的研究重點以及發展方向。雷娟綿[1]對大攻角、高轉速條件下的SOCBT(Secant-Ogive-Cylinder-Boattail)進行了數值模擬,并從流場結構闡述了馬格努斯效應的形成機理。劉周[6]采用不同的湍流模型對SOCBT彈丸進行數值模擬,發現不同的湍流模型對計算精度有較大影響。
從公開的文獻來看,國內外對旋轉彈丸進行了廣泛的研究,在高轉速、大攻角條件下,對不同形狀的無翼彈丸采用不同的數值方法和湍流模型進行了大量研究,并對馬格努斯效應作出了相應的機理分析,但缺少在低轉速且小攻角條件下,對尾翼式彈箭氣動特性隨馬赫數、轉速以及滾轉角進行定量分析與研究。本文采用滑移網格技術,先對ANF標模進行數值方法的有效性驗證,而后對某自旋式尾翼彈進行模擬,從給出的壓力云圖對馬格努斯效應產生原因進行分析與解釋,并且定量地給出了馬格努斯力、升力系數等氣動特性隨轉速、滾轉角以及馬赫數的變化關系,為飛行控制系統設計提供了相關依據。
滑移網格技術是上世紀90年代興起的一種基于非穩態思想的數值模擬方法,近年來得到了廣泛的應用,尤其是在彈箭數值仿真方面。滑移網格原理如圖1所示,整個流動區域被分為包含123單元和456單元的2個不同的區域,其中一個區域相對另一個區域做相對滑移運動來模擬真實的流動過程。ABC和EFG是各子區域與滑移面重合的邊界面。在數值模擬過程中,123單元和456單元在運動過程中會形成B′、C′、D′、E′、F′、G′節點,并沿著交界面做相對滑移運動。其中子塊2的信息通過面C′D′和面D′E′傳遞到子塊4和子塊5,其他單元之間的信息傳遞與此類似。
滑移網格技術作為一種特殊的動網格[7],適用于具有周期性旋轉運動物體的流場計算,而且具有計算精度高、占用內存少等優勢[1]。當模擬彈箭旋轉運動時,只需要一個非旋轉區域和一個包圍彈體的旋轉區域,兩者之間存在一個交界面,交界面的網格節點數不必相同,僅僅在滑移網格交界面進行數值插值,即可保證2個區域之間通量守恒。圖2(a)和圖2(b)給出了網格區域劃分示意圖,圖中d為彈體直徑。
采用基于滑移網格下的三維積分守恒型非定常N-S方程組,由質量守恒方程、動量守恒方程和能量守恒方程組成[1]:
式中:ρ,U,e分別為氣流的密度、速度和單位體積的能量,vm為彈體的運動速度,p為靜壓力,I為單位張量,τ為粘性張量,q為熱通量,S和V分別為控制體的面積分和體積分的積分區域,r為S的外法向單位向量。
采用有限體積法對空間進行離散,對流項二階通量差分分裂(FDS)AUSM+格式,粘性項采用中心差分格式[8]。
Realizablek-ε模型是對RNGk-ε湍流模型做了修正而產生的,能更為準確地反映流體真實的流動過程,尤其是旋轉流動計算中氣流的流動狀態[9]。
對自旋彈箭非定常的研究,采用雙時間步長技術,即采用物理時間步長和內循環時間步長。在進行計算時,一般在一個時間步內彈體旋轉0.1°~0.3°比較合適,在這個范圍內既能夠保證流場空間數據的連續性,又能使得滾轉物體有足夠的旋轉角度,這樣可以大大節省收斂所需的時間[10]。
選取ANF標模進行數值算法的有效性驗證,ANF模型經過大量風洞實驗和飛行試驗驗證,被廣泛采用為動導數的計算標模。表1給出計算條件,其中馬赫數Ma、雷諾數Re、總壓p0和總溫T0與風洞實驗條件一致,彈箭旋轉速度Ω=385.8 rad/s,Ω*為對應的無量綱轉速,Ω*=dΩ/c,d為彈體直徑,c為當地聲速,N為外循環中每轉動一周的總計算的迭代次數,Δt為外循環的計算步長,i為每一次外循環計算中內循環的迭代次數。

表1 計算條件
利用ICEMCFD,采用多塊對接網格生成技術對整個計算區域生成六面體網格,計算區域共生成900多萬網格。縱向采用C型網格,橫向采用O型網格。第一層距離物面為1.0×10-5m,網格伸展比為1.1,控制第一層網格到物面的距離和網格伸展比的目的是為了保證y+<1。圖3(a)和圖3(b)分別是彈體頭部和尾部附近網格。
彈箭表面設置為無滑移壁面條件,彈體跟隨旋轉區域以相同的轉速同步運動。非旋轉區域和旋轉區域的交界面設置為滑移邊界條件;非旋轉區域的外邊界設置為壓力遠場條件。

從圖4(a)和圖4(b)可知,俯仰力矩和法向力系數與實驗值以及參考文獻中計算值高度吻合,誤差均在5%以內;從圖4(c)和圖4(d)可知,在攻角α<30°時,馬格努斯力和力矩動導數與實驗值以及參考文獻中的計算值吻合程度較好。當攻角為60°時,本文數值計算結果與實驗值相差為20%,這是由于當攻角過大時,旋轉彈體附近流場變得更為復雜。綜上所述,在小攻角范圍內,本文采用的數值計算方法具有較高的可信度。
計算模型為某尾翼彈,由錐形部、圓柱段和尾翼組成,其尺寸外形如圖5所示。選取的計算狀態為:來流馬赫數Ma=0.9,1.1,1.4,1.6,2.0,來流溫度T0=288 K。該彈旋轉速度Ω=π,2π,4π rad/s,此時無量綱轉速分別為Ω*=0.003 45,0.006 9,0.013 8,來流攻角α=2°。一個計算周期采用1 440個物理時間步,即Δt=0.000 69 s(3種轉速采用相同的物理時間步和步長)。氣動力和力矩計算的參考點位于彈體質心位置,參考面積為彈體的最大橫截面積,參考長度為彈體直徑d。
為了研究網格的收斂性,在保證y+≤1的情況下,對來流條件為Ma=1.1,攻角α=2°,轉速Ω*=0.013 8,分別以網格數M=6×106,7×106,8×106,9×106共4種不同網格數對該彈箭流場進行數值模擬。圖6給出了瞬時馬格努斯力系數Cmg隨滾轉角Φ的變化規律,可知:隨著網格密度的進一步增大,馬格努斯力系數的最大幅值差也不斷變小,6×106個網格與9×106個網格的幅值差異最大為20.2%,7×106個網格與9×106個網格幅值最大相差7.8%,8×106個網格與9×106個網格幅值最大相差2.3%,可見隨著網格的加密,瞬時馬格努斯力系數隨滾轉角的變化趨于收斂。為了進一步驗證網格的收斂性,表2給出Ma=1.1,攻角α=2°,轉速Ω*=0.013 8,不同網格數下氣動特性計算結果。可知:平均升力系數CL、平均阻力系數CD、平均馬格努斯力系數Cmg隨網格數M的增大逐漸收斂。因此,后序計算以900萬網格為準。

MCLCDCmg6×1060.652190.336020.0033017×1060.618140.317040.0031418×1060.614060.304890.0030779×1060.600740.300390.003001
圖7和圖8分別給出了彈箭圓柱段和尾翼段在不同轉速下的壓力分布云圖,其中觀察方向為彈體尾部向前看。由于攻角和彈體自轉的存在,會產生一個繞彈體縱軸x旋轉的環流(環流方向即為彈體旋轉方向)和一個指向y軸正向的橫流。由圖7(a)可知,當彈體旋轉時,彈體左側橫流方向與環流方向相反,這使得左側氣流速度降低,而右側環流與橫流方向相同,則氣流速度增大。最終使得氣流對彈體左側壓力大于其對彈體右側的壓力,即彈體左右存在壓力差,所以產生了指向彈體右側的側向力,即馬格努斯力。綜合對比圖7(a)~圖7(c)發現,隨著轉速的增大,左側高壓區明顯擴大,而右側低壓區也在逐漸擴大,彈體左右兩側壓差不斷增大,所以馬格努斯力也在增大。通過比較不同轉速下的壓力云圖,可知彈體轉速是影響馬格努斯力大小的一個重要因素。從圖8可知,當彈箭處于低速旋轉狀態時,在迎風區域,由于橫流的影響遠遠大于環流的影響,因而該區域的尾翼兩側區域壓力分布較為對稱,所以產生的馬格努斯力也較小。而背風區由于只受到環流的影響,所以尾翼兩側出現較為明顯的壓力分布不對稱現象,產生較大的馬格努斯力。而且隨著轉速的增大,對背風區的影響也要遠遠大于對迎風區的影響。通過以上分析,并且綜合比較圖7和圖8可知:轉速是馬格努斯效應的影響因素之一,并且彈體處馬格努斯力方向與環流方向一致,而尾翼處馬格努斯力方向與環流方向相反。
通過數值模擬獲得了尾翼彈在各種飛行條件下的平均阻力系數、平均升力系數、平均滾轉力矩系數等氣動特性數據,如圖9所示。
圖9(a)~圖9(c)是不同轉速下尾翼彈的平均阻力系數CD、平均升力系數CL以及平均俯仰力矩系數Cm隨馬赫數的變化曲線,可知:在全聲速范圍內,縱向氣動特性(阻力、升力、俯仰力矩)隨馬赫數變大呈先增大后減小的趨勢,且縱向氣動特性不隨轉速的變化而變化。
圖9(d)是不同轉速下平均滾轉力矩系數Crm隨馬赫數的變化曲線,并且每個馬赫數下的滾轉力矩均為負值,這是由于該導彈模型不存在尾翼斜置角,導彈無法獲得使其主動滾轉的力矩(即導轉力矩為0)而只存在阻礙其旋轉的滾轉阻尼力矩,故而平均滾轉力矩系數均為負值;隨轉速增加,平均滾轉力矩系數增大,與轉速成定比例關系;在同一轉速下,隨馬赫數增大,平均滾轉力矩系數的絕對值呈先增加后減小的變化趨勢,與常規滾轉力矩系數變化規律一致。
圖9(e)給出了不同轉速下平均馬格努斯力Cmg隨馬赫數的變化曲線,可知:平均馬格努斯力系數隨馬赫數的增加呈先增大后減小的趨勢,各個馬赫數下的平均馬格努斯力系數與轉速為非線性增大關系。圖9(f)是轉速Ω*=0.013 8、馬赫數Ma=1.1的條件下,瞬時馬格努斯力隨滾轉角的變化曲線,可知:在一個旋轉周期內有4個小周期,每個小周期的馬格努斯力系數隨滾轉角近似呈現正弦曲線變化。在滾轉角從0°到90°的過程中,馬格努斯力系數先增大,當轉到Φ=27°左右時達到峰值,然后減小,在Φ=69°左右時達到波谷,然后再次增大。瞬時馬格努斯力隨滾轉角變化會產生4個相同周期,這是由于該彈箭的尾翼布局為“十字型”,彈箭每旋轉90°時,此時彈箭再一次處于“十字型”布局,彈箭周圍的流場性質與前一次“十字型”布局流場性質相同。自旋彈箭在飛行過程中,馬格努斯力可以分解為粘性分量(由于流體粘性的存在而產生的那一部分馬格努斯力)和壓差分量(彈體兩側的壓力差產生的那一部分馬格努斯力)。表3給出了在轉速Ω*=0.013 8、馬赫數為1.1且不同滾轉角下,各個部件的粘性分量Cmg,S(通過物面摩擦應力積分得到的馬格努斯力)和壓差分量Cmg,P(通過物面壓力積分得到的馬格努斯力)對馬格努斯力系數Cmg的貢獻值。可知,粘性分量對馬格努斯力系數的貢獻始終為正值(對馬格努斯力系數的最大貢獻比例不超過5%)。當滾轉角Φ=0°,50°,90°時,由于彈體和尾翼產生的馬格努斯力始終是一個正值,而另一個為負值,所以總的馬格努斯力在0值上下浮動。當滾轉角Φ=27°時,彈體和尾翼產生的馬格努斯力均為正值,這使得總的馬格努斯力最大(波峰)。而當滾轉角Φ=68°時,彈體和尾翼產生的馬格努斯力均為負值,這使得總的馬格努斯力最小(波谷),即在一個小周期內馬格努斯力系數隨滾轉角近似呈現正弦曲線變化。
圖9(g)給出轉速Ω*=0.013 8下不同部件馬格努斯力系數隨馬赫數的變化曲線。可知:彈體、尾翼、全彈的平均馬格努斯力系數均隨馬赫數的增大呈先增大后減小的趨勢。并且馬格努斯力主要是由尾翼產生的,尤其是在Ma=0.9,2時,彈體的馬格努斯力幾乎可以忽略。各部件平均馬格努斯力系數如表4所示。表中Δ表示部件提供平均馬格努斯系數百分比。
對比圖9(a)~圖9(c),可知轉速對平均阻力、升力、俯仰力矩系數沒有影響,而非定常計算所需的時間遠遠大于定常計算時間,所以可以采用定常方法來得到升力、阻力、俯仰力矩系數。對比圖9(a)~圖9(e)發現:平均滾轉力矩系數與平均馬格努斯力系數的數值大小約為平均升力的1%~10%,這是由于該彈箭的尾翼屬于“十字型”布局,并且當彈箭不旋轉時,彈箭攻角平面兩側流場對稱分布,不會產生馬格努斯力和滾轉力矩。而當彈體逆時針低速旋轉時,由于旋轉產生環流的線速度遠遠小于攻角產生的橫流速度,所以旋轉引起的氣動系數屬于小量。

表3 不同滾轉角下的馬格努斯力系數

表4 不同部件馬格努斯力系數
本文選取ANF標模進行數值算法的有效性驗證。數值計算結果表明,本文選取的計算方法具有較高的計算精度。在此基礎上,研究了某低速旋轉尾翼彈箭氣動特性隨馬赫數、轉速以及滾轉角的變化規律,得到以下結論:
①通過比較不同馬赫數和轉速下的縱向氣動數據,可知縱向氣動特性不隨轉速變化而變化。隨馬赫數增加,縱向氣動特性均呈現先增加后減小的變化規律,與常規靜態氣動變化規律相同。
②隨馬赫數增加,平均滾轉力矩系數絕對值呈現先增加后減小的變化趨勢,與常規滾轉力矩變化規律一致。隨轉速增加,平均滾轉力矩系數增大,且平均滾轉力矩系數與轉速呈現典型的線性變化規律。
③隨馬赫數增加,平均馬格努斯力系數呈現先增加后減小的變化趨勢;隨轉速增加,平均馬格努斯力系數增大,且平均馬格努斯力系數與轉速為非線性變化關系。在同一轉速和馬赫數下,瞬時馬格努斯力系數與滾轉角為近似正弦曲線變化規律,在一個滾轉周期內,瞬時馬格努斯力有4個自然變化周期,每一個周期為90°,即尾翼從十字型布局變化到下一個十字型布局為一個變化周期。瞬時馬格努斯力隨滾轉角的變化在一個較大正值和較小負值間振動。
[1] 雷娟綿,李田田,黃燦.高速旋轉彈丸馬格努斯效應數值研究[J].兵工學報,2013,34(6):718-725.
LEI Juan-mian,LI Tian-tian,HUANG Can.A numerical investigation of magnus effect for high-speed spinning projectile[J].Acta Armamentarii,2013,34(6):718-725.(in Chinese)
[2] DESPIRITO J.Lateral jet interaction on finned projectile in supersonic flow:AIAA-2012-0413[R].Reston:AIAA,2012.
[3] STUREK W B,KAYSER L D,NIETUBICZ C J.Computations of Magnus effects for a yawed,spinning body of revolution[J].AIAA,1978,7(16):687-692.
[4] PECHIER M,GUILLEN P,CAYZAC R.A combined theoretical-experimental investigation of Magnus effects[C]//16th AIAA Applied Aerodynamics Conference.Barges:AIAA,1998:778-788.
[5] 苗瑞生,吳甲生.旋轉彈空氣動力學[J].力學進展,1987,17(4):479-488.
MIAO Rui-sheng,WU Jia-sheng.Aerodynamics of spinning projectiles[J].Advances in Mechanics,1987,17(4):479-488.(in Chinese)
[6] 劉周,謝立軍,楊云軍,等.彈丸旋轉空氣動力效應非定常數值模擬[J].航空學報,2016,37(5):1 401-1 410.
LIU Zhou,XIE Li-jun,YANG Yun-jun,et al.Unsteady numerical simulation of aerodynamic effect of a spinning projectile[J].Acta Aeronautica et Astronautica Sinica,2016,37(5):1 401-1 410.(in Chinese)
[7] 米百剛,詹浩,朱軍.基于CFD數值仿真技術的飛行器動導數計算[J].空氣動力學學報,2014,32(6):834-839.
MI Bai-gang,ZHAN Hao,ZHU Jun.Calculation of dynamic derivatives for aircraft based on CFD technique[J].Acta Aerodynamica Sinica,2014,32(6):834-839.(in Chinese)
[8] 趙博博,劉榮忠,郭銳,等.扭曲尾翼彈箭的馬格努斯數值研究[J].固體火箭技術,2015,38(4):465-471.
ZHAO Bo-bo,LIU Rong-zhong,GUO Rui.Numerical prediction of the Magnus effect for twist fin swept flight projectile[J].Journal of Solid Rocket Technology,2015,38(4):465-471.(in Chinese)
[9] 吳萍.二維彈道修正彈氣動布局分析[D].南京:南京理工大學,2014.
WU Ping.Aerodynamic configuration analysis of two-dimension trajectory correction projectile[D].Nanjing:Nanjing University of Science & Technology,2014.(in Chinese)
[10] 宋琦,楊樹興,徐勇,等.滾轉狀態下卷弧翼火箭彈氣動特性的數值模擬[J].固體火箭技術,2008,31(6):552-554.
SONG Qi,YANG Shu-xing,XU Yong.Numerical simulation on aerodynamic characteristic of rolling rocket with wrap around fins[J].Journal of Solid Rocket Technology,2008,31(6):552-554.(in Chinese)
[11] BHAGWANDIN V A.Numerical prediction of roll damping and magnus dynamic derivatives for finned projectiles at angle of attack[C]//30th AIAA Applied Aerodynamics Conference.New Orleans:AIAA,2012:1-17.
ResearchonAerodynamicCharacteristicsofMissileWithFinsUnderLow-speedRotation
ZHANG Chao,WANG Xue-de,WANG Yong,CAO You-quan
(School of Energy and Power Engineering,Nanjing University of Science and Technology,Nanjing 210094,China)
In order to research the effect of low-speed rotation on aerodynamic characteristics of missile with fins,the flow-field around the missile with fins under small attack-angle and total Mach number was simulated by applying 3D unsteady N-S equations and slide mesh technology.The Army-Navy Basic Finner Missile was used for verifying the validity of the numerical simulation method.The results show that this method has high accuracy.The results under different Mach numbers,rotation speed and roll angle show that the longitudinal aerodynamic characteristics(lift,drag,pitching moment)is independent of spin.The average rolling-moment-coefficient is contantly proportionate to rotation speed,and the average Magnus-force-coefficient exhibits a nonlinear change with the rotation speed.
missile;Magnus effect;sliding mesh;aerodynamic characteristics;numerical prediction
TJ303.4
A
1004-499X(2017)04-0040-08
2017-05-02
江蘇省青藍工程資助項目;國家留學基金委出國留學資助項目
張超(1991- )男,碩士研究生,研究方向為飛行器設計。E-mail:15996298561@163.com。
王學德(1977- ),男,副教授,研究方向為航空/航天飛行器外流場的數值模擬算法。E-mail:wangxuede2000@njust.edu.cn。