劉 冰,劉 科,薛瑞娟,張凌雲(yún),樂貴高

(1.南京理工大學 機械工程學院，南京 210094；2.中國船舶集團公司第七一三研究所，鄭州 450000)

1 引言

在現(xiàn)代導彈潛射技術中，從彈體有無附加保護角度分類，可分為濕式發(fā)射和干式發(fā)射2種。干式發(fā)射借助專門的運載器包裹形式完成水下運動和出水跨介質(zhì)運動；在水面附近導彈與運載器實現(xiàn)分離動作，導彈轉入空中飛行彈道，完成整個導彈水下發(fā)射過程。運載器出水時介質(zhì)由水變?yōu)榭諝?，并受海洋中風浪的強烈干擾作用，在氣液交界面處作用時間短、變化快、變化因素多，飛行器外部力學環(huán)境急劇變化，容易產(chǎn)生劇烈的振動與沖擊，帶來出水載荷的特性不明確和飛行穩(wěn)定性問題。彈筒分離階段雖然時間歷程很短，卻賦予了導彈空中飛行的初始條件，直接影響導彈空中彈道的啟控條件和穩(wěn)定性，對潛射導彈成功發(fā)射至關重要。導彈的彈射點位置不同會造成導彈發(fā)射過程中最終離開運載器的姿態(tài)產(chǎn)生差異。若運載器出水初始姿態(tài)角小于90°，受力形態(tài)一般是非對稱的。這種非對稱力可引起飛行器姿態(tài)偏離設計狀態(tài)。若沖出水面的飛行器無法達到預定的姿態(tài)指標，即使點火仍將導致失穩(wěn)和發(fā)射失敗。所以深入開展浮動平臺水面分離特性影響研究非常必要。

劉曜計算分析了艇速對潛射導彈水中運行彈道的影響，對導彈與運載器水面分離彈道也進行了數(shù)值計算研究，但是沒有考慮波浪、潛流和陣風對分離運動的影響。彭正梁采用勢流理論計算水動力，對高速航行體與運載器水面分離過程進行了二維彈道仿真，但是所做的二維仿真分析研究沒有考慮到彈器分離過程中導彈姿態(tài)角的變化。段崇一等建立導彈水下發(fā)射二維流場網(wǎng)格模型，基于二維不可壓N-S方程，對導彈出管過程中不同時刻彈體的壓力載荷分布進行解算，但是采用二維的網(wǎng)格模型，難以完全考慮導彈所處的真實環(huán)境。李晶對水動力環(huán)境下的兩剛體相對運動過程進行數(shù)值仿真，計算得出潛射導彈發(fā)射的速度和姿態(tài)角變化，但是沒有考慮彈器分離過程中適配器對分離過程的影響，也沒有在彈器分離模型施加波浪的作用。顧媛媛采用切片法分析水面分離時浮筒和航行器的受力情況，利用MATLAB建立了浮筒與航行器分離的仿真模型，計算得到航行器的姿態(tài)變化，但建立的彈器分離模型不能準確計算分離過程中浮筒與流體的相互耦合作用，也沒有考慮水面波浪的影響。張紅軍等用Mixture等模型和動網(wǎng)格技術，對考慮潛艇平動影響的導彈冷彈射出筒過程進行了三維非定常數(shù)值模擬。陸辰昱使用VOF方法比較分析了橫向流速和單一的傾斜角度對運載器式潛射導彈分離的影響，沒有涉及不同的彈射時機和出水角度對運載器式潛射導彈分離特性的影響，也沒有考慮到適配器對彈射分離的影響。目前，有關彈射點位置和波浪對運載器式潛射導彈水面彈道彈射分離影響的研究未見報道，而且更多僅局限于二維的分離運動，沒有充分考慮導彈發(fā)射過程中姿態(tài)角的變化，未能考慮液體表面邊界條件，用流動的任何區(qū)域氣體和液體模糊的處理，給計算自由液面的運動特性帶來較大誤差，而且沒有考慮到導彈與運載器分離過程中適配器對分離過程產(chǎn)生的影響，缺少彈器動態(tài)分離過程中適配器的相互作用模型。

本文建立高精度VOF模型和六自由度耦合運動模型，對運載器出水過程中不同彈射起始高度、不同出水角度下運載器式潛射導彈水面分離的質(zhì)心位移，姿態(tài)角及軸向速度等進行了三維仿真分析。本文提出的仿真方法具有較高的精度，所形成的高精度數(shù)值仿真方法以及潛射導彈水面分離研究特性結果，對工程應用具有良好的參考價值。

2 分離過程及計算模型建立

運載器從水下發(fā)射到出筒過程主要包括：以水面平均高度作為基準水面高度，運載器在基準水面高度以下1.5 m開始運動，在水中依靠初速和浮力迅速上升，當運載器以一定速度沖出水面至一定高度，達到分離時機，頭罩首先分離，導彈助推發(fā)動機點火工作產(chǎn)生推力，運載器受到同樣大小的反向推力，導彈與運載器間的鎖定銷被推力剪斷，導彈沿著運載器筒體快速滑行，適配器相繼脫離，最后導彈離開運載器。建立固定坐標系，運載器固定坐標系，航行體固定坐標系，位于初始時刻運載器頂端中心，位于運載器質(zhì)心，位于航行體質(zhì)心。沿著運載器中心軸方向，沿著航行體中心軸方向。運載器動坐標系及導彈動坐標系與固定坐標系之間的歐拉角為滾轉角，俯仰角，偏航角。水面分離示意圖如圖1。運載器式潛射導彈采用單筒單彈結構。運載器的直徑為0.7 m，長度為6 m，材料密度為1 656.46 kg/m；導彈直徑為0.375 m，長度為3.9 m，材料密度為1 179.69 kg/m，航行體質(zhì)心與頭部頂點距離2.56 m，運載器質(zhì)心距頭部頂點3.3 m，適配器剛度為700 kN/m，阻尼系數(shù)為200 N·s/m。

圖1 水面分離示意圖Fig.1 Diagram of water surface separation

為研究復雜海情對運載器動平臺分離特性的影響，共設立4組計算工況，具體計算條件如表1所示，分別進行彈射起始高度、初始俯仰角對動平臺分離特性影響研究。工況1彈射起始高度為運載器頭部出水面2 m，工況3彈射起始高度為運載器頭部出水面4 m，與工況2彈射起始高度為運載器頭部出水面3 m進行對比；工況4的初始傾斜角度為10°，與工況2的初始傾斜角度為5°進行對比。

表1 計算工況Table 1 Computational conditions

3 數(shù)值模擬方法的建立

3.1 自由液面模型

基于分數(shù)容積障礙網(wǎng)格法，建立笛卡爾坐標系下的三維不可壓縮Navier-Stokes方程，形式如下：

(1)

式(1)中：是流體密度；是壓力；，，是流體在時刻，在點(,,)處的速度矢量；是單位體積流體受的外力，若只考慮重力，則=；常數(shù)是運動粘度。

與標準-模型相比，RNG-湍流模型考慮了平均應變率的影響，能更好地處理高應變率及流線彎曲較大的流動，所以可以更精確地計算運載器出水過程波浪破碎的狀況。RNG-模型的湍流產(chǎn)生項和耗散項分別表示為：

方程：

(2)

方程：

(3)

式(3)中：湍流黏性=RNG；模型常數(shù)RNG=0085,2RNG=168,RNG=RNG=0071 9。

3.2 VOF方法

對自由液面采用由Hirt提出的流體體積函數(shù)(VOF)進行求解。VOF方法利用計算網(wǎng)格單元中流體體積量的變化和網(wǎng)格單元本身體積的比值函數(shù)來確定自由面的位置和形狀。若在某時刻網(wǎng)格單元中=1，則說明該單元全部為指定相流體所占據(jù)，為流體單元，若=0，則該單元全部為另一相流體所占據(jù)，相對于前相流體，則稱為空單元；當0<<1時，則該單元為包含兩相物質(zhì)的交界面單元。采用VOF兩相流計算模型，定義水為第一相、空氣為第二相。在數(shù)學上通常假定任意函數(shù)(,)，其定義為:

(4)

若不考慮劇烈的相變，根據(jù)連續(xù)性介質(zhì)概念，函數(shù)隨質(zhì)點運動而保持不變，其隨體導數(shù)為零，流體體積傳輸方程的形式為:

(5)

VOF方法的數(shù)值差分格式見式(6)，在Δ,中對式(5)進行積分，其積分形式為:

(6)

3.3 FAVOR方法

FAVOR(fractional area/volume obstacle representation)網(wǎng)格處理技術能夠在結構化的網(wǎng)格內(nèi)定義獨立復雜的幾何體，即使模型很復雜也可以精確地描述外型，實現(xiàn)利用簡單的矩形網(wǎng)格表示任意復雜的幾何形狀，避免了以往有限差分難以精確擬合物體邊界的缺點。FAVOR方法可以計算網(wǎng)格表面的開放面積分數(shù)和開放體積分數(shù)，并基于這些參數(shù)重建幾何結構。這種方法提供了簡單而準確的方法,來表示域中的復雜曲面，而不需要貼體網(wǎng)格。

采用的FAVOR方法本質(zhì)上是一種浸沒邊界方法，網(wǎng)格生成的工作量和難度遠小于貼體網(wǎng)格方法，模型運動與背景網(wǎng)格無關，無論模型的運動幅度有多大，網(wǎng)格均無需劃分。采用的FAVOR方法比ALE動網(wǎng)格方法更適用于求解大幅度運動的流固耦合問題，并且具有較高的精度。與傳統(tǒng)的浸沒邊界法相比，F(xiàn)AVOR方法對邊界的刻畫精度更高，可用來精確模擬高雷諾條件的流動，在一定條件下彌補了貼體網(wǎng)格的不足。FAVOR方法能夠準確計算界面面積，增強數(shù)值穩(wěn)定性，計算沿固體邊界的平流和應力，可以精確地計算得到固液交界面模型的相對運動。

3.4 六自由度運動模型

作用于高速航行體和運載器的外力和外力矩為和，則運動方程為：

(7)

式(7)中：和分別為剛體的線動量和角動量；為旋轉角速度；為線速度。

剛體運動的三維彈道方程組為：

(8)

式(8)中：、、分別表示運載器或航行體的角速度矢量在其動坐標系軸向上的分量；、、為流體附加質(zhì)量，、為流體附加靜矩，、、為流體附加轉動慣量；右端項是作用在剛體上的合力與合力矩在各坐標軸上的分量。航行體受到的作用力包括：重力2,燃氣發(fā)生器推力2，適配器約束內(nèi)力2,下標為各道適配器的標號。運載器所受作用力包括：重力1，浮力1，流體慣性作用力1，流體粘性作用力1，適配器約束力1，燃氣發(fā)生器反推力1，波浪擾動力1，流體慣性作用力項已包含在彈道方程的左端。

采用隱式數(shù)值算法計算物體運動并處理流體和運動物體之間的相互作用。潛射導彈密度比水更大，在與水面進行接觸碰撞過程中產(chǎn)生很大的非線性力，附加質(zhì)量也很大，采用隱式算法，計算運動的穩(wěn)定性更好、收斂更快。動量對流采用顯示的數(shù)值逼近算法，設置了最大時間步尺寸來維持數(shù)值穩(wěn)定性，與隱式算法相比較具有更高的計算精度。

3.5 適配器作用力

在航行體發(fā)射出筒的過程中，適配器隨著航行體一起向外移動，適配器的作用力可以分為阻尼力和彈性力。在航行體發(fā)射過程中，適配器的變形是非常微小的，認為適配器的變形在彈性范圍之內(nèi)，適配器的彈性力正比于變形量。阻尼力與適配器的阻尼系數(shù)和適配器節(jié)點的相對運動速度有關。

適配器作用力計算公式為：

(9)

3.6 計算邊界條件

設置長12 m、寬7 m、高19.8 m的整體區(qū)域，液面在距網(wǎng)格區(qū)域底面高9.5 m處。采用六面體網(wǎng)格單元，使用包含2個網(wǎng)格塊組合的嵌套網(wǎng)格，網(wǎng)格塊1包含運載器及其周圍流體區(qū)域，長3 m、寬2 m、高17.7 m，設置網(wǎng)格邊長0.04 m，網(wǎng)格數(shù)量166.1萬個，網(wǎng)格塊2的范圍占據(jù)整體計算區(qū)域，網(wǎng)格邊長0.1 m，網(wǎng)格數(shù)量166.3萬個，整體網(wǎng)格總數(shù)為332.4萬個。在網(wǎng)格區(qū)域的左側和右側邊界設置對稱邊界，前側和后側邊界設置波浪來流和出流，下側設置壁面邊界，上側設置壓力出口邊界。針對包含運載器及其周圍流體區(qū)域的內(nèi)側網(wǎng)格塊1進行網(wǎng)格無關性驗證，網(wǎng)格無關性分析設置如表2所示，設置0.03 m、0.04 m、0.05 m等3種網(wǎng)格，得出3種網(wǎng)格工況初始傾斜角度為5°條件下筒彈分離時刻運載器的俯仰角，0.04 m網(wǎng)格尺寸與0.03 m工況運載器俯仰角比較接近并且網(wǎng)格數(shù)較少，所以選用0.04 m的網(wǎng)格尺寸，既能保證較高的計算精度而且計算工作量較小。邊界條件類型如圖2所示。網(wǎng)格無關性分析設置表如表2所示。

表2 網(wǎng)格無關性分析設置表Table 2 Grid independence analysis setting

圖2 邊界條件類型示意圖Fig.2 Boundary condition type

3.7 模型驗證

根據(jù)上述自由液面模型，對含旋轉角速度的標準臺球入水過程進行數(shù)值模擬，將數(shù)值結果與實驗對比，驗證數(shù)值方法精度和有效性。標準臺球直徑=0.057 2 m，質(zhì)量170 g，初始入水速度為=5.45 m/s，并且以=199 rad/s的角速度旋轉。圖3給出了計算速度云圖與實驗高速攝影的對比。臺球入水具有初始的旋轉角速度，既包含平動也包含轉動，可以說明6自由度模型的可靠性。

圖3 速度場計算結果與實驗高速攝影圖像Fig.3 Calculated results of velocity magnitude and high-speed photographs

圖4分別表示了臺球入水過程中質(zhì)心運動軌跡和豎直方向位移的計算值和試驗值。由圖4可以看出，質(zhì)心運動軌跡和豎直方向位移的計算值與試驗值較為吻合，所以仿真模型具有較高的精度。

圖4 質(zhì)心運動軌跡和豎直方向位移曲線Fig.4 Trajectory of mass center and vertical displacement

4 計算結果及分析

4.1 彈射起始高度對分離特性的影響研究

選取不同的彈射起始高度，開展跨介質(zhì)運動仿真，并且考慮導彈發(fā)射過程中適配器運動過程，得出彈射起始高度對海面分離的影響規(guī)律。這里選取彈射起始高度為運載器的頭部高出水面2 m、頭部高出水面3 m和頭部高出水面4 m，海浪的浪級都是5級，出水速度都是9 m/s,初始傾斜角度均為5°，且初始俯仰角速度也均為5°/s的3種不同工況(見圖5)，也就是工況1、工況2、工況3進行對比。得出3種發(fā)射狀態(tài)下運載器和導彈在、、軸的坐標隨時間的變化曲線。

圖5 工況1(a)，工況2(b)，工況3(c)分離狀態(tài)圖像Fig.5 The image of separating state in working condition 1 (a),2 (b)and 3 (c)

質(zhì)心位置

圖6和圖7表示運載器和導彈的位移沿軸、軸、軸的變化，給出了工況1、工況2、工況3的位移。

圖6 工況1、工況2、工況3運載器位移曲線Fig.6 Displacement of launch tube in working condition 1,2 and 3

圖7 工況1、工況2、工況3導彈位移曲線Fig.7 Displacement of missile in working condition 1,2 and 3

由圖6、圖7可知，工況3運載器在方向位移最大，工況1運載器在方向位移最小。工況3導彈在方向位移最大，工況1和工況2導彈在方向位移基本相同。這3種工況運載器和導彈在方向有輕微波動。工況3運載器在方向位移最大，工況1運載器在方向位移最小。而且3種工況運載器在發(fā)射過程中始終沿方向運動。工況3導彈在方向位移最大，工況1導彈在方向位移最小。

姿態(tài)角

圖8和圖9是出水過程運載器和導彈的姿態(tài)角曲線，給出了工況1、工況2、工況3的姿態(tài)角。

由圖8、圖9可知，在俯仰角、偏航角、滾轉角這3個姿態(tài)角里，俯仰角和偏航角對導彈能否成功發(fā)射的決定性更大。3種工況運載器和導彈俯仰角都由5°逐漸減小，工況1運載器最后的俯仰角為-3°；工況2運載器最后的俯仰角為-66°；工況3運載器最后的俯仰角為-11°，表明工況3運載器的俯仰角的變化程度大于工況2，工況2運載器的俯仰角的變化程度大于工況1。3種工況導彈俯仰角由5°逐漸減小，工況1導彈最后的俯仰角為-33°；工況2導彈最后的俯仰角為-44°；工況3導彈最后的俯仰角為-89°，表明工況3導彈的俯仰角的變化程度大于工況2，工況2導彈的俯仰角的變化程度大于工況1。3種工況運載器和導彈偏航角波動范圍不超過3°，表明這3種工況運載器和導彈偏航角的變化都比較小。3種工況運載器的滾轉角從0°逐漸增大，工況3運載器的滾轉角為254°，比工況1和工況2的滾轉角更大，工況1和工況2的滾轉角基本相同。工況3導彈的滾轉角變化程度大于工況2，工況2的導彈滾轉角變化程度大于工況1。工況3運載器和導彈姿態(tài)角變化幅度大于工況2，工況2運載器和導彈姿態(tài)角變化幅度大于工況1。彈射起始高度越高，運載器和導彈姿態(tài)角變化幅度越大。

圖8 工況1、工況2、工況3運載器姿態(tài)角Fig.8 Attitude angles of launch tube in working condition 1,2 and 3

圖9 工況1、工況2、工況3導彈姿態(tài)角Fig.9 Attitude angles of missile in working condition 1,2 and 3

軸向速度

圖10表示出水過程運載器和導彈軸向速度，給出了工況1、工況2、工況3的軸向速度。

由圖10可知，3種工況運載器的初始速度為9.5 m/s，從平均水面以下1.5 m開始運動，當運動到平均水面的出水速度為9 m/s。在整個發(fā)射過程中，運載器軸向速度逐漸減小。工況1運載器0.46 s開始受到彈射力的作用，工況2運載器0.59 s開始受到彈射力的作用，工況3運載器從0.77 s開始受到彈射力的作用，3種工況運載器軸向速度最后減小至-1.5 m/s。在整個發(fā)射過程中，3種工況導彈軸向速度先略微減小，在施加彈射力以后速度逐漸增大，并且最后增大到20 m/s左右滿足出筒速度要求。3種工況施加彈射力以后導彈的加速度基本相同。

圖10 工況1、工況2、工況3運載器和導彈軸向速度曲線Fig.10 Axial speed of launch tube and missilein working condition 1,2 and 3

4.2 動平臺初始傾斜角度對分離特性的影響研究

選取典型的發(fā)射傾斜角度，此次研究針對出水姿態(tài)角接近垂直的運載器，因此分別選取初始傾斜角度為5°和10°的2種不同工況(見圖11)，也就是工況2和工況4進行對比，海浪的浪級都是5級，出水速度都是9 m/s，且初始俯仰角速度也均為5°/s，得出初始發(fā)射傾斜角度對分離特性的影響規(guī)律。

圖11 工況2(a)，工況4(b)分離狀態(tài)圖像Fig.11 The image of separating state in working condition 2 (a)and 4 (b)

質(zhì)心位置

圖12和圖13分別表示運載器和導彈的位移沿軸、軸、軸的變化，給出了工況2和工況4的位移。

圖12 工況2、工況4運載器位移曲線Fig.12 Displacement of launch tube in working condition 2 and 4

圖13 工況2、工況4導彈位移曲線Fig.13 Displacement of missile in working condition 2 and 4

由圖12、圖13可知，工況4比工況2運載器在正方向有更大的位移，工況4比工況2導彈在正方向也有更大的位移，工況2和工況4運載器和導彈在方向有輕微波動，工況4比工況2運載器和導彈在正方向位移略微更大，工況2和工況4運載器和導彈在方向位移也基本相同。

姿態(tài)角

圖14和圖15表示出水過程運載器和導彈姿態(tài)角，給出了工況2和工況4的姿態(tài)角。

由圖14、圖15可知，在俯仰角、偏航角、滾轉角這3個姿態(tài)角里俯仰角和偏航角對導彈能否成功發(fā)射的決定性更大。工況2和工況4運載器和導彈偏航角波動范圍不超過25°，表明工況2和工況4運載器和導彈偏航角的變化都比較小。工況2運載器俯仰角從5°逐漸減小至-66°；工況4運載器俯仰角從10°減小至-106°，表明工況4比工況2發(fā)射過程中運載器俯仰角的變化幅度更大。工況2導彈俯仰角從5°逐漸減小到-44°；工況4導彈俯仰角從10°減小到-79°，表明工況4比工況2發(fā)射過程中導彈俯仰角的變化幅度更大。工況2運載器的滾轉角從0逐漸增大到16°，工況4運載器的滾轉角從0增大到12°，工況2導彈滾轉角從0增大到16°，工況4導彈滾轉角從0增大到152°，表明工況2和工況4發(fā)射過程中運載器和導彈滾轉角變化基本相同。整體來看工況4比工況2發(fā)射過程中運載器和導彈的姿態(tài)角變化更大。

圖14 工況2、工況4運載器姿態(tài)角曲線Fig.14 Attitude angles of launch tube in working condition 2 and 4

圖15 工況2、工況4導彈姿態(tài)角曲線Fig.15 Attitude angles of missile in workingcondition 2 and 4

軸向速度

圖16表示了工況2和工況4運載器和導彈軸向速度曲線。

由圖16可知，運載器軸向速度逐漸減小，不施加彈射力時2種工況運載器軸向速度先略微減小，在施加彈射力以后，2種工況運載器軸向速度最后減小至-0.5 m/s。當不施加彈射力時，2種工況導彈軸向速度先略微減小，在施加彈射力以后，2種工況導彈軸向速度在推力作用下最終都增加到20 m/s左右滿足出筒的速度要求。彈射過程中工況4運載器軸向速度減小趨勢慢于工況2，工況4導彈軸向速度增大趨勢慢于工況2。

圖16 工況2、工況4運載器和導彈軸向速度曲線Fig.16 Axial speed of launch tube and missilein working condition 2 and 4

出水傾斜角度優(yōu)化分析

進行傾斜角度對運載器式潛射導彈分離影響的優(yōu)化分析。設置出水時刻導彈的傾斜角度為3°，仿真計算得出相應的筒彈分離時刻導彈俯仰角為0.82°。由于出水時刻傾斜角度5°相應的筒彈分離時刻導彈俯仰角為-4.4°，所以筒彈分離時刻導彈俯仰角為0對應的出水時刻傾斜角度在3°～5°。設置出水時刻導彈的傾斜角度為3.2°仿真計算得出相應的筒彈分離時刻導彈俯仰角為0.1°接近為0，導彈接近垂直彈射離開運載器，又設置出水時刻導彈傾斜角度3.3°進行仿真計算得到筒彈分離時刻導彈俯仰角為-0.4°，所以最優(yōu)的出水傾斜角度應該是3.2°。

5 結論

本文對不同彈射位置、出水角度的運載器式潛射導彈分離特性進行了數(shù)值模擬分析，得到以下結論：

1)建立了一套復雜水動力環(huán)境下的運載器式潛射導彈出水分離仿真方法，計入運載器跨越水氣介質(zhì)交界面的影響，添加了潛射導彈出水過程的慣性力、粘性力和波浪擾動力的因素，綜合考慮適配器在筒彈分離過程的影響，能夠模擬運載器和導彈從水下發(fā)射到筒彈分離結束的完整過程。

2)各種工況下運載器和導彈在方向和方向有較大的位移，在方向有輕微的晃動,各種工況下運載器和潛射導彈的姿態(tài)角和軸向速度的變化穩(wěn)定性較好，證明了這種分離方式的合理性。

3)發(fā)射過程中比較工況1、工況2、工況3，彈射起始高度距離水面位置越大，運載器和導彈的位移和姿態(tài)變化更大，其中俯仰角增加比較明顯。

4)比較工況2、工況4，傾斜10°出水比傾斜5°出水在發(fā)射過程中運載器和導彈的位移和姿態(tài)角，發(fā)現(xiàn)俯仰角的變化更大，而且運載器和導彈軸向速度減小趨勢更慢，最優(yōu)的出水傾斜角度是3.2°。