波浪影響下的潛射導彈出水姿態擾動分析

楊繼鋒，陳 捷，劉丙杰，李淑媛

(1.海軍潛艇學院， 山東 青島 266099； 2.火箭軍駐北京地區軍事代表局， 北京 100085)

潛射導彈通常在水下40 m以內深度組織發射，彈體運動過程經歷了復雜物理環境的交替變化，先后穿越了固海和海空界面，期間受到平臺運動、海流、海浪、海面風場等水文氣象因素以及“筒口效應”和“空化效應”等特殊物理現象的影響[1]。導彈水下發射技術極其復雜，涉及到多項流、非線性，尺度效應、負雜耦合等問題。其中波浪作為表層海水最重要的擾動形式，導彈穿行的水域是波浪擾動最為明顯的區域，波浪力的沖擊作用是導彈姿態擾動的重要影響因素[2]。

對于潛射導彈水下發射的數值仿真，傳統的數學模型很難對非線性特征進行合理的特征捕捉和仿真計算，造成仿真結論與試驗結果有較大差異。為了能夠更為真實的研究全濕流狀態下導彈發射過程中的相關規律，需要借助專業的流體仿真軟件，解決發射環境的數值模擬、空化與多項流模型的構建，進行動邊界的設計、出水自由面的處理以及流固耦合負雜計算等多方面問題。

本文對全濕流狀態的潛射導彈，以ANSYS 14.0軟件為仿真平臺，進行發射環境的構造和導彈水中運動規律的數值仿真，研究波浪對導彈運動姿態的擾動規律。

1 基于ANSYS軟件的彈道仿真分析

水下發射流場的構建，其數學方程分別包含代表質量守恒的連續性方程，代表動量守恒的N-S方程和代表能量守恒的能量方程[3]。考慮粘性影響，需要根據湍流模型，補充湍動能及湍流耗散率方程。自由界面和空泡界面的捕捉需要對空氣、水、蒸汽等多項體積積分求得，因此需引入相應變量，體積分數滿足對流方程[4]。

1.1 控制方程

控制方程的引用由所研究的問題所決定，主要功能是控制所研究問題的發展變化。本文將多項流介質作為可變密度的單一流體，各項同享壓力場、速度場，采用控制流動的基本方程作為連續性方程和動量方程，在此基礎上補充了蒸汽量連續方程。對于流體力的加載，利用ANSYS自帶的網格節點內的積分算法進行解算[5]。具體的控制方程如下。

連續方程：

(1)

動量方程和能量方程分別為：

(2)

(3)

其中，混合物有效傳遞系數：

(4)

蒸汽相和氣相連續性方程分別為：

(5)

(6)

1.2 湍流模型

導彈變速穿行于流體中，運動過程引起彈體周圍流體的復雜變化，具有明顯的湍流特征。湍流引起的波動具有小尺度和高頻率的特征，介質間伴隨著動能的相互轉換和液體濃度的交織變化。k-ε湍流模型可以對平板和圓柱射流的發散比率進行精確的預測，而且對于旋轉流動、強逆壓梯度的邊界層流動、流動分離和二次流有很好的表現，適宜描述導彈運動過程中的流體特性。本文選用該模型模擬湍流環境，脈動動能方程和動能耗散率方程分別為：

Gk+Gb-ρε-YM

(7)

(8)

式中，Gk、Gb分別表示由均速度梯度和浮升力引起的湍流動能；YM表示壓縮性的影響；C2和C1ε為常數；σe和σk分別表示k和ε的湍流Prandtl數。其中：

模型中取值：C1ε=1.44,C2=1.9,σe=1.2，σk=1.0。

1.3 空化模型

對于自然空泡流動，由于相變導致密度不為常數，造成基本流動方程不封閉，需要加載空化模型，建立密度與其他物理量的關系。對空泡的數值計算有兩類方法，一種是基于氣液交界面跟蹤的邊界元法，另一種是忽略泡內流動，將流場視為整體，通過求解Navier-Stokes方程的多項流法。本文以后者為方向，采用預處理的準雙時間步長法，建立混合介質體積連續性方程和蒸汽相體積分數的運輸方程，模型如下[7]：

(9)

(10)

(11)

2 基于UDF宏命令的仿真環境構造

2.1 彈體運動模型UDF設計

彈體運動模型采用ANSYS軟件進行二次開發，通過DEFINE_CG_MOTION宏命令進行運動過程的干預和控制。為保證彈體周圍網格質量，將彈體及周圍流體計算域共同作為運動體。這種方法雖然保證了網格質量，但卻犧牲掉了 ANSYS 六自由度運動參數計算功能。為此，彈體的運動參數需要通過流體力進行計算和讀取并賦值給運動區域。編譯后的UDF程序可以確保以庫文件的形式動態鏈接到運動區域。此時，該區域整體將保持剛性運動，而區域外的網格仍然按照重構的方式進行處理。

2.2 基于層變的環境模型網格構造

考慮導彈質心可能的運動范圍和波浪模型的完整性，設計選取計算域范圍為40 m×20 m×70 m。40 m的長度可以保證計算域中能夠形成至少1個完整的波形，盡量避免出口回流對計算的影響。如圖1所示，計算域按網格形式設計為2個部分。在彈道周圍主要計算區域采用非結構化網格，以適應導彈的不規則幾何外形和運動過程中多自由度的平移和旋轉運動。為保證計算過程中的網格質量，上下變形區各層網格保持同步生長和潰滅，彈體區網格相對獨立，不受影響。在距離導彈運行軌跡較遠處以結構化網格為主，降低網格密度，減少計算量。兩種網格采用對接方式連接[8]。在導彈出筒的初始階段，受發射筒的徑向制約，導彈基本做一維運動，此時彈體周圍的柱形區域網格隨彈體一同運動，該區域網格采用交接面與外圍網格銜接。隨著彈體出筒，新的網格在彈體下方及筒口平面生成。

圖1 彈體網格和計算域網格構造示意圖

邊界條件設置：彈體設置為壁面邊界條件(wall)，左側面為速度入口(velocity-inlet)，右側面為自由出流邊界(outflow)，前后面對稱設置(symmetry)。由于自由出流邊界條件不能與壓力出口邊界條件同時使用，所以模型將上下兩面都定義為壁面[9]。雖然這與實際情況略有不符，但由于計算域較大，且在豎直方向上沒有大范圍的流動出現，故此邊界設置可以滿足計算要求。選用標準k-ε模型作為湍流模型，壓力速度耦合方式采用PISO格式。其中水氣交界采用VOF模型設置水體區域，自由液面初始位置根據發射深度設定為40 m。

對于不同的問題，初始條件的影響會有所不同，如果初始條件給的不好，計算過程會逐漸擺脫初始條件的影響，產生較大的計算誤差。由于導彈彈射出筒，速度很高，如果計算開始時刻即進行運動的模擬，很容易造成計算結果發散。因此，在模擬飛行過程之前，先賦予流場一定的初值，在靜默條件下進行一定步數的靜態求解，使流場中的壓力、速度、密度等元素的分布達到接近于真實的水平[10]。然后將迭代后的流場作初始流場開始仿真計算，從而消除初始條件誤差帶來的影響。

計算過程中，前置一段導彈運動空間，設定導彈初速為0，尾部施加推力，使導彈加速運動，在流場穩定時加速到出筒速度，位置達到模擬位置，取消推力。這樣做使導彈在正式開始計算時，彈體周圍流場已經有了一定的分布形態，和實際中導彈出筒時的流場較為接近。

2.3 波浪模型構造

在發射環境的構造中，對于海浪的模擬是核心工作，合理的方案設計以及正確的邊界條件設置是模擬高品質海浪的基礎。數值波浪模型的建立包括造波及消波方式的選擇、對造波邊界的設置、VOF模型的應用以及計算步長的設置等諸多問題。

目前，數值造波技術分為基于勢流理論的邊界元法、基于解不可壓縮流體的粘性數值造波法和理想流體源函數造波法三類。具體方法主要有以下5種:空間周期波、模擬造波機造波法、入射速度法造波、預先描述入射波法和引入離散奇點作為造波源造波的方法。另外，通過對水波方程展開形式的不同，可以得出不同的波浪理論，如線性波理論、二階Stokes波理論、三階及高階Stokes波理論等。其中，大連理工大學齊鵬等采用有限差分法對對粘性不可壓流體NS方程進行離散，結合VOF方法和k-ε湍流模型，采用速度邊界造波的方式，對三維數值波浪水池中規則波進行了良好的模擬。

本文將通過入射速度法造波，模擬生成線性波。此外，計算域邊界條件設置如下:左側為速度入口,右側和上側為壓力出口邊界,下側為壁面條件,同時設定VOF模型為明渠流。由于基于DEFINE_PROFILE宏命令設定的波浪擾動無法滿足在速度入口邊界上完全一致，需要植入函數對計算域邊界網格進行掃描，確定網格中心坐標并進行函數賦值，即可實現均一速度入口的開發編制[11]。

計算過程中，造波中涉及到氣液交界面采用VOF模型處理方法，并考慮流體的重力，流場計算物性參數見表1。綜合考慮計算效率和圖形效果，設置計算域中單元格面積在0.01～0.02 m2，1、3、5級波浪模型參數如表2所示，網格總數控制在24 100以下，時間步長選擇在0.001 s。利用CFD軟件經模擬計算的三級海浪的壓力和速度云圖見圖2。

表1 流場計算物性參數

表2 波浪參數

圖2 三級波浪壓力云圖(a)和速度云圖(b)

3 仿真分析

按上述原理和方法構建仿真環境，利用FLUENT軟件，模擬平臺速度為1 m/s條件下的水下發射，研究導彈運動過程壓力場、速度場的特性以及波浪對彈體姿態的擾動情況。

如圖3和圖4所示，導彈在水下運動過程中，彈頭對流體的沖擊提升了該區域流體的能量，形成了局部的高能流體團，該區域的動能和壓力都高于其他區域，流體在高壓作用下向四周加速流動。而沿彈體方向的流體受彈體幾何形狀的限制，流域呈收縮狀，速度進一步增加，在肩部形成壓力最低的區域，并將進一步演化成肩部空化區。隨著影響區域的增大，影響效果逐漸降低，此過程形成壓力場的劇烈變化，此變化隨著距彈體距離的增加而逐漸衰減。導彈尾部受推力控制形成高壓區，推力隨著彈體的升高逐漸降低，高壓區影響范圍也逐漸減小。首尾高壓區不沿軸線對稱，形成驅動彈體轉動的力矩，使導彈姿態發生偏轉。在接近彈體表面區域存在與導彈具有相當速度的流層，為粘性作用引起的流體運動。

圖3 彈體頭部壓力云圖(a)和速度云圖(b)

圖4 彈體尾部壓力云圖(a)和速度云圖(b)

圖5和圖6為三級海浪下不同時刻導彈運動時的截面相圖和截面壓力云圖，可以較為清晰的反應導彈的運動狀態。

圖5 T=0.812 s時刻，彈體截面相圖(a)和截面壓力云圖(b)

圖6 T=1.224 s時刻，彈體截面相圖(a)和截面壓力云圖(b)

取導彈軸線上彈頭、彈尾兩個端點A和B，如圖7所示，兩個端點在X方向上的位移分別為XA、XB，在Z方向的位移分別為ZA、ZB，位移的大小是隨著時間變化的。導彈在彈體坐標系內俯仰姿態的計算方法如下：

俯仰角偏差Δφ：

(12)

圖7 導彈出水俯仰角偏差示意圖

導彈水中運動過程姿態變化曲線如圖8和圖9。波浪對彈體姿態影響十分明顯，俯仰角速度在波浪運動區域具有明顯的增大趨勢，隨著彈體出水，角速度變化得以逐漸抑制，呈下降趨勢。導彈俯仰姿態呈逐漸增高的趨勢，并且隨浪級的增高而增大。在1級海浪條件下，導彈出水時的俯仰角偏差可達為7.86°，當浪級增至5級時，出水俯仰角偏差增至13.2°，變化十分明顯。

圖8 導彈俯仰角變化曲線

圖9 導彈俯仰角速度變化曲線

4 結論

基于FLUENT軟件的潛射彈道導彈出水過程仿真，能夠較為客觀的反映導彈的運動規律，特別是在流體作用影響以及非線性特征提取方面具有明顯優勢。該仿真結論對彈體水中運動姿態變化的體現細致具體，導彈運動姿態受波浪的影響程度和變化規律與波浪力的脈動規律基本一致，計算結果可以為深入研究海洋環境因素對導彈發射安全影響提供參考，豐富該研究領域的仿真技術經驗，為實現潛射彈道導彈大深度安全發射提供理論支撐。