蔡禮港,常秋英,楊 超,段曉亮
(1.北京交通大學 機械與電子控制工程學院, 北京 100044; 2.中國兵器工業導航與控制技術研究所, 北京 100089)
火箭彈主要用于遠程作戰,適宜攻擊面積大、定位不太精確的集群目標。高超聲速火箭彈一般是指在大氣層內或跨大氣層飛行速度大于5馬赫的火箭彈[1],與常規火箭彈相比,由于氣動加熱現象的存在,彈體在飛行過程中壁面附近的溫度急劇升高,達到一定溫度后發生燒蝕[2]。對舵翼而言,燒蝕使其氣動外形發生明顯改變,導致彈道軌跡變化甚至飛行失穩,高超聲速飛行對火箭彈的氣動性能與防熱保護措施提出更高的要求[3-4]。
舵翼一般采用高強度鋼制成,在高超音速飛行過程中,舵翼內部伴隨熱傳導、燒蝕等現象,而舵翼結構內部熱力場的改變又會影響外圍流場流動特性與氣動加熱特性,產生氣動熱環境與材料熱響應的耦合效應[5]。耦合計算的求解方式有兩種,分別為緊耦合與松耦合[6-7]。緊耦合方式是將流體區域與固體區域的控制方程統一求解[8],由于兩區域控制方程的離散、求解方法均不相同,統一求解難度大,難以用于解決實際復雜工程問題。松耦合方法對其進行了簡化處理,把流體區域與固體區域分開求解[9-11],并在流場和固體場之間傳遞數據,更新每一次迭代后的固體場與流體場,這種簡化方法對于實際工程應用是可以接受的[7,12]。
為了獲得某型號高超聲速火箭彈舵翼的氣動熱環境和燒蝕過程,并以此作為舵翼熱防護措施的設計和優化依據,以D6AC鋼舵翼為例,基于真實飛行工況,采用有限元方法進行氣動熱與燒蝕仿真,獲得動態燒蝕過程與燒蝕形貌,并與電弧加熱風洞實驗結果進行對比。
直角坐標系下的三維非定常可壓N-S方程是連續性方程、能量方程、氣體狀態方程,可表示為:
(1)
式中,x,y,z分別為三個直角坐標方向變量,t為時間變量,Q為流場守恒變量,Ec,Fc,Gc分別為三個坐標方向上的無粘對流通量,Ev,Fv,Gv分別為三個坐標方向上的粘性耗散通量,具體表達式參見文獻[13]。
湍流流動的數值模擬主要是求解Reynolds平均N-S方程,也就是對Reynolds平均方程中的Reynolds應力項引入湍流模型,使方程組成封閉的可解方程。方程中的粘性系數由層流粘性系數和湍流粘性系數組成:
μ=μl+μt
(2)
其中μl為層流粘度,由Sutherland公式近似給出:
(13)
式中:μ0為粘性系數,T為環境溫度,TS為蘇士南常數,TC為273.16 K。
μt為湍流粘度,由湍流模型給出。采用Spalart-Allmaras湍流模型。
對舵翼計算域進行結構化網格劃分,在外部流場的前端、后部與上下區域進行局部加密,如圖1所示。
仿真計算中,舵翼主要由D6AC鋼制成,其外形如圖2所示。舵翼體網格劃分為六面體網格,以保證燒蝕邊界的規整,便于實現流固區域數據映射,如圖3所示。
燒蝕過程需要使用不同的網格調整策略:當舵翼表面燒蝕速率較低時,采用網格自適應方法進行調整;當燒蝕速率較高時,需要采用生死單元法調整網格。具體而言,使舵翼燒蝕區域的固體網格失效實現固體邊界后退,相應地激活在燒蝕區域的流體網格實現流體邊界的前進。需要注意的是,流體區域被激活的網格,在劃分網格階段就已經生成了,只是在開始計算時將其設置為失效。
以熱力學侵蝕為基礎,不考慮內部氣體溢出和復雜化學反應現象,結合舵翼飛行環境,采用的舵翼燒蝕模型如下:
(4)
由于氣動熱仿真與燒蝕仿真采用不同的網格。數據映射的目的在于把氣動熱仿真所得到的各點溫度、壓強和切應力等參數信息傳遞到燒蝕網格中。具體而言,針對燒蝕網格的每一節點,在氣動熱網格中搜索對應位置附近的8個節點,使用8個節點的數據插值計算燒蝕網格節點的溫度與壓強,直到整個流固交界面計算完畢。
首先,搜索找到附近8個象限中距離中心點最近點,相互聯結形成圍繞插值點的六面體單元。如圖4所示,六面體的中心點9為插值點。
確定六面體單元后,使用三維空間反距離加權插值方法,待估點的值等于各個已知點值的加權平均值。插值計算公式如下:
(5)
其中fi是搜索到的8個數據點的函數值,u為常數。
(6)
di(x,y,z)表示點(x,y,z)到點(xi,yi,zi)的距離。定義冪次u的大小可調節相鄰數據對插值點的影響。
對飛行器實際飛行工況的氣動加熱狀況進行分析,根據實際飛行速度、高度隨時間的變化擬合成曲線,實現非穩態過程的仿真計算。
圖5標示了用于定量化比較的四個監控點A、B、C、D和四條特征直線1~4的位置。圖6~圖8給出了四個監控點的溫度、壓強和切應力隨時間的變化曲線,各圖均對飛行時間進行無量綱處理。舵翼溫度、壓強與切應力的變化與飛行馬赫數變化趨勢相近:在加速階段,隨著馬赫數的增大,舵翼表面溫度、切應力和壓強等參數迅速增大;在20%飛行時間左右,各參數均達最大值;20%飛行時間之后,各參數快速減小。
根據氣動熱仿真結果,定義燒蝕網格中節點的溫度、壓強以及剪切力,建立各參數的變化規律,結合前面的燒蝕損傷模型,計算舵翼燒蝕表面各個節點的退化速率,求得燒蝕量。圖9給出了D6AC鋼舵翼特征直線的燒蝕比率變化情況。
在圖9中,前10%飛行時間為氣動加熱過程,舵翼在15%飛行時間左右開始發生燒蝕;在20%飛行時間左右燒蝕速率達到最大,此時舵翼前緣溫度、壓力與切應力也達到最大;在30%飛行時間燒蝕速率趨于平緩;到40%飛行時間,燒蝕比率保持穩定,燒蝕基本結束。燒蝕比率變化與氣動熱仿真中溫度、壓強與切應力的變化趨勢吻合。D6AC鋼舵翼的線燒蝕比率為19.16%,前緣區域被燒蝕,中后部區域保持完整。仿真得到的舵翼燒蝕形貌如圖10所示。
根據高超聲速火箭彈的飛行馬赫數和彈道高度變化,在仿真過程得到相應焓值。在中國航天空氣動力技術研究院電弧加熱風洞實驗室FD04風洞中進行實驗,模擬D6AC鋼舵翼在飛行過程中的燒蝕損傷,評估其耐燒蝕能力。由于電弧加熱風洞不能實現出口焓值連續變化,通過逐級調整逼近理論焓值,如圖11所示。
電弧加熱風洞實驗布局如圖12所示,被電弧加熱的高溫氣體通過噴管膨脹在試驗段形成高溫、高速氣流,在給定參數范圍內模擬舵翼在飛行過程中的氣動加熱環境。舵翼試樣使用D6AC鋼作為主要材料,樣件外形與實際舵翼保持一致,如圖13所示。實驗中,把火箭彈舵翼試件放在在矩形噴管出口處,并且確保舵翼與來流無偏角,前緣與噴管出口平行距離10 mm,模型上端點距離噴管出口上方10~15 mm。實驗采用了軌道模擬試驗技術,通過逐次調節流入電弧加熱器的氣體流量和電弧加熱器的電源功率,利用多個臺階近似模擬飛行過程中瞬態參數的變化。
實驗結果如圖14所示。根據試樣的燒蝕外形判斷,由于氣動現象產生的熱量,使得試件溫度臨近D6AC鋼的熔點,D6AC鋼舵翼上部前緣發生明顯燒蝕。大部分熔化的D6AC鋼被高速氣流帶離舵翼本體,進入高速流場的下游。而小部分熔融狀態的D6AC鋼在表面壓強及壁面切應力等因素作用下被擠壓,附著在舵翼下游表面并且冷卻后形成殘渣。另外,在進行風洞試驗時,固定支架導致附近流場的流動狀況與仿真相比發生變化,部分熱量沿著固定支架散失,導致在根部前緣位置舵翼相對保持完整。
最后,將燒蝕仿真結果(圖10)與實驗結果(圖14)對比可知,仿真得到的燒蝕外形與實驗燒蝕結果相差不大。通過稱重、測量表面積等手段,獲得舵翼質量燒蝕比率與面積燒蝕比率如表1所示。量化結果表明,仿真燒蝕結果的誤差在工程應用允許誤差范圍內。可使用此仿真方法預測不同材料舵翼的動態燒蝕過程與最終燒蝕形貌,作為設計、優化舵翼防熱保護措施的依據。

表1 D6AC鋼舵翼實驗與仿真結果
1) 氣動熱仿真結果表明,在飛行過程中,D6AC鋼舵翼溫度、壓強與切應力的變化先迅速增加,在20%飛行時間到達最大值,隨后一直減小直至飛行結束;
2) 燒蝕仿真結果表明,D6AC鋼舵翼在15%飛行時間開始發生燒蝕,在20%飛行時間燒蝕速率達到最大,在30%飛行時間燒蝕速率趨于平緩,到40%飛行時間,燒蝕基本結束。最終,前緣區域被燒蝕,中后部區域基本完整;
3) 電弧加熱風洞實驗結果表明,D6AC鋼舵翼根部前緣區域相對完整,其余前緣區域被燒蝕,中后部區域基本完整;
4) 對實驗和仿真結果的量化對比表明,仿真誤差在工程應用允許誤差范圍之內。未來可以使用該仿真方法獲得不同材料舵翼的動態燒蝕過程與燒蝕區域,為舵翼材料的選擇、熱防護的設計和優化提供依據。