周東輝，施紅輝，魯建華，薛明瑞

(1. 浙江工業(yè)職業(yè)技術(shù)學(xué)院，浙江 紹興 312000；2. 浙江理工大學(xué) 機(jī)械與自動(dòng)控制學(xué)院，浙江 杭州 310018)

0 引 言

并聯(lián)出水問(wèn)題來(lái)源于水下彈道導(dǎo)彈的齊射出水、水下槍炮多管并射形成彈幕攔截或攻擊低空飛行目標(biāo)等海戰(zhàn)背景，具有很強(qiáng)的工程應(yīng)用價(jià)值。當(dāng)采用并聯(lián)發(fā)射方式潛射導(dǎo)彈/射彈時(shí)，可以顯著地增加武器的突防概率和毀傷效果。并聯(lián)出水過(guò)程中不僅涉及了流體介質(zhì)的突變，氣相、固體、液相的三相耦合流動(dòng)，還存在航行體之間的流場(chǎng)彼此干擾，導(dǎo)致航行體的流場(chǎng)和運(yùn)動(dòng)特性非常復(fù)雜。

導(dǎo)彈/射彈以較高速度水下發(fā)射時(shí)，其表面發(fā)生空化形成空泡，空泡的演化過(guò)程及形態(tài)關(guān)系到水下航行體的受力情況，進(jìn)而影響到航行體的彈道特性。因此，空泡演化特性的研究一直是該領(lǐng)域的研究重點(diǎn)。在國(guó)外，LOGVINOVICH 等[1]提出的空泡截面獨(dú)立膨脹原理，可用于計(jì)算空泡的形狀。SAVCHENKO 等[2]提出了計(jì)算超空泡輪廓的半經(jīng)驗(yàn)公式。WAUGH 等[3]通過(guò)實(shí)驗(yàn)研究了發(fā)射角度和空泡對(duì)導(dǎo)彈出水姿態(tài)的擾動(dòng)影響，給出了出水空泡的形態(tài)。NGUYEN 等[4]利用數(shù)值模擬方法研究了射彈勻速出水過(guò)程的超空泡流動(dòng)，給出了射彈出水過(guò)程的阻力變化和超空泡形態(tài)。在國(guó)內(nèi)，賈會(huì)霞等[5,6]利用高速攝像技術(shù)開(kāi)展了超空泡射彈出水的實(shí)驗(yàn)研究，獲得了出水過(guò)程中空泡的演化過(guò)程和射彈的速度變化，給出了不同模型和參數(shù)下超空泡出水形態(tài)尺寸的變化規(guī)律。施紅輝等[7]采用實(shí)驗(yàn)和數(shù)值模擬相結(jié)合的方法研究了超空泡射彈出水的現(xiàn)象，分析了射彈頭部形狀對(duì)射彈運(yùn)動(dòng)、自由面破碎、超空泡形狀的影響。別慶等[8]基于Mixture 多相流模型，采用動(dòng)網(wǎng)格技術(shù)和6DOF 技術(shù)建立了潛射導(dǎo)彈水下運(yùn)動(dòng)模型，獲得了變速狀態(tài)下導(dǎo)彈的空泡特性。魯傳敬等[9]基于多相流模型，采用動(dòng)網(wǎng)格技術(shù)數(shù)值模擬了細(xì)長(zhǎng)體出水過(guò)程。魏海鵬等[10]比較了不同多相流模型在航行體出水過(guò)程中對(duì)流體界面捕獲、壓力計(jì)算方面的差異，給出了VOF 模型和Mixture 模型的適用范圍。顏開(kāi)等[11]對(duì)出水空泡的行為進(jìn)行了很多數(shù)值計(jì)算，給出了空泡周?chē)俣葓?chǎng)和壓力場(chǎng)的變化規(guī)律。ZHANG 等[12]數(shù)值模擬了出水空泡潰滅的過(guò)程，獲得了出水速度對(duì)空泡形態(tài)及潰滅壓力的影響。陳瑛等[13]利用大渦模擬的數(shù)值方法對(duì)導(dǎo)彈出水過(guò)程中空泡脫離和潰滅進(jìn)行了研究。

針對(duì)并聯(lián)出水問(wèn)題，MNASRI 等[14]利用VOF 模型和動(dòng)網(wǎng)格技術(shù)進(jìn)行了雙圓柱體低速并聯(lián)出水的二維數(shù)值模擬，分析了出水過(guò)程中自由面的變形以及流場(chǎng)的相互干擾。盧佳興等[15]進(jìn)行了回轉(zhuǎn)體齊射出水實(shí)驗(yàn)，分析了回轉(zhuǎn)體間的流體動(dòng)力干擾和艇速對(duì)其影響。畢鳳陽(yáng)等[16]建立了多細(xì)長(zhǎng)體水下齊射多相流動(dòng)與多體運(yùn)動(dòng)耦合數(shù)值模擬的計(jì)算模型，分析了典型工況下的齊射擾動(dòng)特性。

綜上所述，現(xiàn)有的研究主要集中在單航行體出水過(guò)程的空泡演化、流體動(dòng)力以及彈道特性等，但對(duì)高速航行體并聯(lián)出水的超空泡流動(dòng)特性研究還偏少。本文基于VOF 多相流模型，引入Schnerr-Sauer 空化模型、SST κ -ω湍流模型、6DOF 剛體運(yùn)動(dòng)模型，利用重疊網(wǎng)格技術(shù)建立高速航行體并聯(lián)出水的數(shù)值模型，模擬發(fā)射速度為200 m/s 時(shí)，單航行體出水以及不同并聯(lián)間距下雙航行體出水的自由運(yùn)動(dòng)過(guò)程，獲得了航行體的超空泡演化特性和運(yùn)動(dòng)特性，可為并聯(lián)超空泡射彈出水的工程應(yīng)用提供參考。

1 控制方程及數(shù)值方法

1.1 控制方程

采用VOF 多相流模型描述水、水蒸汽、空氣構(gòu)成的多相流動(dòng)系統(tǒng)，該模型將三相當(dāng)作密度可變的單一介質(zhì)的混合相來(lái)處理，各相共享同一壓力差和速度場(chǎng)?；旌舷嗟倪B續(xù)性方程和動(dòng)量方程分別為：

式中： ρm為混合相密度， ρm=αvρv+αgρg+α1ρ1， αv，αg， α1分 別為水蒸氣、空氣、水的體積分?jǐn)?shù)；ui為i軸方向的速度；P為壓力； μm為混合相動(dòng)力粘度，μm=αvμv+αgμg+α1μ1， μv、 μg、 μ1分別為水蒸汽、空氣、水的動(dòng)力粘度；uj為j軸方向的速度；gi為重力加速度的分量，F(xiàn)i為體積力的分量。

采用SST κ -ω湍流模型[17]對(duì)RANS 方程提供湍流封閉?？栈F(xiàn)象采用Schnerr-Sauer 空化模型描述[10]，其控制方程如下：

1.2 邊界條件、網(wǎng)格劃分及數(shù)值方法

數(shù)值模型中采用的航行體模型為截錐型射彈，如圖1 所示。該彈由1 個(gè)圓臺(tái)和1 個(gè)圓柱組成，質(zhì)量為2.98 g，圓柱段直徑D=6 mm，錐段長(zhǎng)L0=18 mm，射彈全長(zhǎng)L=48 mm，空化器直徑D0=3 mm。

圖1 幾何模型示意圖Fig. 1 Schematic diagram of geometric model

計(jì)算域和邊界條件設(shè)置如圖2 所示，計(jì)算域?yàn)殚L(zhǎng)方體，高780 mm，長(zhǎng)210 mm，寬120 mm。計(jì)算域的底部為壓力入口，計(jì)算域的側(cè)面和上部為壓力出口，采用用戶(hù)自定義場(chǎng)函數(shù)定義邊界面上的壓力?？諝庥蚋?00 mm，水域高480 mm，航行體頭部距離自由液面的高為400 mm，航行體A 和航行體B 的表面均設(shè)定為無(wú)滑移壁面。上述計(jì)算模型中壓力與速度耦合的求解采用C o u p l e d 算法，壓力場(chǎng)和空間離散采用PRESTO!格式，體積率離散采用Modified HRIC。定義兩航行體軸線(xiàn)之間的距離為并聯(lián)間距△d。航行體的質(zhì)心在y方向的位移為縱向位移Sy，在y方向的速度為縱向速度Vy，取豎直向上為正，定義航行體無(wú)量綱縱向位移

圖2 計(jì)算域及邊界條件設(shè)置Fig. 2 Computational domain and boundary condition setting

在數(shù)值模擬中重疊網(wǎng)格技術(shù)是研究流場(chǎng)中多物體運(yùn)動(dòng)最為直接有效的方法之一，因此本文選擇重疊網(wǎng)格技術(shù)實(shí)現(xiàn)2 個(gè)航行體并聯(lián)出水的數(shù)值模擬。計(jì)算區(qū)域分為背景網(wǎng)格區(qū)域和子網(wǎng)格區(qū)域，由于航行體A 和航行體B 分別獨(dú)立運(yùn)動(dòng)，因此需要?jiǎng)澐? 個(gè)子網(wǎng)格區(qū)域。圖3 為計(jì)算域的網(wǎng)格劃分示意圖，均采用結(jié)構(gòu)化網(wǎng)格進(jìn)行劃分，全計(jì)算域的網(wǎng)格總數(shù)約為216.5 萬(wàn)。背景網(wǎng)格區(qū)域中的航行體運(yùn)動(dòng)區(qū)域和自由液面附近進(jìn)行局部加密，以便精確捕捉空泡界面和自由液面的變形，外部背景區(qū)域網(wǎng)格較疏，加快收斂速度。2 個(gè)子網(wǎng)格區(qū)域大小相同，均為包裹航行體的圓柱，其長(zhǎng)度為10D，直徑為2D，并且采用相同的網(wǎng)格劃分方式，對(duì)航行體壁面附近區(qū)域進(jìn)行網(wǎng)格加密。計(jì)算中，航行體的運(yùn)動(dòng)通過(guò)6DOF 剛體運(yùn)動(dòng)模型求解。

圖3 網(wǎng)格劃分示意圖Fig. 3 Schematic diagram of meshing

1.3 數(shù)值方法驗(yàn)證

根據(jù)文獻(xiàn)[7]中射彈模型2 的實(shí)驗(yàn)數(shù)據(jù)，對(duì)水下超空泡射彈出水進(jìn)行數(shù)值方法驗(yàn)證。圖4 給出了出水超空泡演化過(guò)程的實(shí)驗(yàn)和數(shù)值模擬結(jié)果，相鄰圖片的時(shí)間間隔為Δt=0.25 ms。通過(guò)對(duì)比可以發(fā)現(xiàn)，數(shù)值模擬較好地模擬出了實(shí)驗(yàn)中的出水超空泡演化過(guò)程。為了進(jìn)一步驗(yàn)證數(shù)值模擬方法有效性，將射彈的無(wú)量綱縱向位移的數(shù)值結(jié)果與實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比，如圖5 所示。可以看出，射彈無(wú)量綱縱向位移的數(shù)值結(jié)果與實(shí)驗(yàn)數(shù)據(jù)具有較好的一致性，最大誤差約為3.9%。通過(guò)上述比較驗(yàn)證了本文數(shù)值模擬方法的有效性。

圖4 射彈出水超空泡的演化過(guò)程對(duì)比（Δt=0.25 ms）Fig. 4 Comparisons of water-exit supercavity evolutions of the projectile between numerical simulation and experiment(Δt=0.25 ms)

圖5 關(guān)于射彈無(wú)量綱縱向位移的數(shù)值模擬和實(shí)驗(yàn)結(jié)果[7]對(duì)比Fig. 5 Comparisons of dimensionless longitudinal displacements of the projectile between numerical simulation and experiment

2 結(jié)果與討論

本文以單航行體出水為參照工況，考慮并聯(lián)間距為2D、3D、4D、5D四種工況，分別對(duì)不同工況開(kāi)展數(shù)值模擬研究，計(jì)算工況如表1 所示。表中VA0，VB0，V0分別為航行體A、航行體B、單航行體出水時(shí)的初速度，方向?yàn)樨Q直向上。

表1 計(jì)算工況表Tab. 1 Numerical simulation cases

2.1 并聯(lián)出水過(guò)程的超空泡演化特性

圖6（a）給出了航行體單獨(dú)出水過(guò)程的超空泡演化過(guò)程?？梢钥闯?，航行體在水下運(yùn)動(dòng)至出水過(guò)程中經(jīng)歷了形成自然超空泡、超空泡隨動(dòng)、超空泡局部潰滅等。在水下運(yùn)動(dòng)階段，由于航行體的運(yùn)動(dòng)速度較高，導(dǎo)致航行體肩部位置的壓力降低到水的飽和蒸汽壓以下，此時(shí)液相水會(huì)不斷汽化形成空泡，隨著空泡不斷的發(fā)展形成超空泡，超空泡的形狀為細(xì)長(zhǎng)的橢球體，航行體被包裹在超空泡內(nèi)部，只有頭部沾濕，如圖6（a）中t=1 ms 所示。隨著航行體繼續(xù)向上運(yùn)動(dòng)，超空泡跟隨航行體向上運(yùn)動(dòng)。航行體穿越水面時(shí)，被航行體帶入空氣中的空泡在大氣壓力的作用下發(fā)生潰滅，如圖6（a）t=2.8 ms 所示。由于航行體出水速度較高，水面以下的空泡來(lái)不及隨航行體出水而被剝離在自由液面以下，如圖6（a）中t=3 ms 所示。當(dāng)航行體完全進(jìn)入空氣之后，出現(xiàn)了向上運(yùn)動(dòng)的水的噴濺，如圖6（a）中t=3.4 ms 所示。6（b）給出了并聯(lián)間距△d= 4D時(shí)，航行體并聯(lián)出水過(guò)程中超空泡演化過(guò)程。為了便于分析，以t=2 ms 的圖像為例，定義2 個(gè)航行體相隔較近的一側(cè)為內(nèi)側(cè)，相隔較遠(yuǎn)的一側(cè)為外側(cè)?？梢钥闯?，2 個(gè)航行體都各自形成了超空泡，航行體被包裹在超空泡內(nèi)部。由于受到相鄰航行體的限制，超空泡內(nèi)側(cè)的擴(kuò)張受到了抑制，而超空泡的外側(cè)自由擴(kuò)展，導(dǎo)致雙航行體超空泡外側(cè)輪廓的曲率要大于內(nèi)側(cè)輪廓，形成了一對(duì)非對(duì)稱(chēng)超空泡，但雙超空泡形態(tài)在空間上呈現(xiàn)出較好的鏡面對(duì)稱(chēng)特征，如圖6（b）中t= 1~1.8 ms 所示。隨著雙航行體繼續(xù)向上運(yùn)動(dòng)，超空泡繼續(xù)發(fā)展，其直徑和長(zhǎng)度增加。在穿越水面階段，并聯(lián)出水與單獨(dú)出水相似，出現(xiàn)水面以上的空泡潰滅、還未潰滅的空泡被自由液面阻隔等行為，如圖6（b）t=2.6~3 ms 所示。另外還可以看出，超空泡的后半段以并聯(lián)航行體的中軸線(xiàn)為基準(zhǔn)相互靠攏。在雙航行體完全出水后，噴濺的水花是不對(duì)稱(chēng)的，外側(cè)噴濺高度大于內(nèi)側(cè)，如圖6（b）中t= 3.4 ms 所示。

圖6 單獨(dú)出水和并聯(lián)出水的超空泡演化Fig. 6 Supercavity evolutions of the single vehicle exiting water and two vehicles exiting water in parallel

圖7為t=2 ms 時(shí)，不同并聯(lián)間距下的出水超空泡形態(tài)對(duì)比?？梢钥闯觯S著航行體并聯(lián)間距的減小，超空泡尾部向內(nèi)側(cè)擴(kuò)張的程度逐漸增加。當(dāng)并聯(lián)間距△d為2D和3D時(shí)，雙超空泡尾部壁面融合形成“U”形空泡，而△d為3D和4D的工況沒(méi)有出現(xiàn)雙超空泡的尾部融合。其原因主要是由于隨著雙航行體內(nèi)側(cè)流體域變小，該區(qū)域接受雙航行體傳遞動(dòng)能的水流體減少，引起每個(gè)流體質(zhì)點(diǎn)獲得的動(dòng)能平均值增加，根據(jù)伯努利定理，內(nèi)側(cè)區(qū)域的平均壓力必然降低，增強(qiáng)了超空泡向內(nèi)側(cè)擴(kuò)張的速度，從而導(dǎo)致雙空泡在內(nèi)側(cè)壁面發(fā)生了融合。

圖7 t=2 ms 時(shí)，不同并聯(lián)間距的出水超空泡形態(tài)對(duì)比Fig. 7 Comparisons of the shape of supercavitation with different parallel distances when t = 2 ms

2.2 航行體并聯(lián)出水過(guò)程中的運(yùn)動(dòng)特性

定義航行體的偏航角θ為當(dāng)前航行體軸線(xiàn)與y軸的夾角，逆時(shí)針旋轉(zhuǎn)方向?yàn)槠浇堑恼较?。圖8 給出了并聯(lián)間距△d= 4D時(shí)，雙航行體偏航角隨時(shí)間變化的曲線(xiàn)。從圖中可知，并聯(lián)出水時(shí)，航行體A 向外側(cè)進(jìn)行偏轉(zhuǎn)，航行體B 向另一外側(cè)進(jìn)行偏轉(zhuǎn)，即2 個(gè)航行體的運(yùn)動(dòng)向著兩者頭部遠(yuǎn)離、尾部靠近的方向偏轉(zhuǎn)，從整體上看航行體A 和航線(xiàn)體B 的偏轉(zhuǎn)運(yùn)動(dòng)表現(xiàn)為良好的對(duì)稱(chēng)性。航行體偏轉(zhuǎn)的原因主要是航行體誘導(dǎo)的超空泡沿軸線(xiàn)不對(duì)稱(chēng)，導(dǎo)致航行體受到了不對(duì)稱(chēng)的水動(dòng)力。在此種工況下，隨著航行體運(yùn)動(dòng)時(shí)間的增加，偏航角也越來(lái)越大，在時(shí)間為3.4 ms 時(shí)，偏航角達(dá)到了3.06○。

圖8 并聯(lián)間距△d = 4D 時(shí)，雙航行體的偏航角隨時(shí)間變化Fig. 8 Variations of yaw angle of two vehicles with time when parallel distance △d = 4D

鑒于兩航行體運(yùn)動(dòng)的對(duì)稱(chēng)性，將基于航行體A 進(jìn)行并聯(lián)出水過(guò)程的運(yùn)動(dòng)分析。圖9 給出了不同工況下航行體的縱向位移、縱向速度隨時(shí)間的變化?？梢钥闯觯? 種工況下航行體的速度變化趨勢(shì)相似，在運(yùn)動(dòng)初始階段，航行體的速度衰減很快，這是由于超空泡尚未形成，受到水的阻力較大；隨著航行體產(chǎn)生穩(wěn)定的超空泡，減小了受到的水的阻力，速度衰減趨勢(shì)逐漸變緩；2.6 ms 之后，航行體的速度衰減變得很小，原因?yàn)?.6 ms 時(shí)，航行體開(kāi)始穿越水面進(jìn)入空氣中，在穿越過(guò)程中，航行體頭部進(jìn)入空氣中，尾部在空泡中，頭部流體介質(zhì)由液相水變成氣相，密度突變，從而導(dǎo)致航行體受到的壓差阻力變得很小，當(dāng)航行體完全進(jìn)入空氣中，航行體只受到量級(jí)較小的空氣阻力。通過(guò)對(duì)比縱向速度和無(wú)量綱縱向位移的放大圖可以發(fā)現(xiàn)，雖然各工況下航行體縱向速度和無(wú)量綱縱向位移相差較小，但是仍有區(qū)別。對(duì)于并聯(lián)出水，隨著并聯(lián)間距的減小，航行體的縱向速度衰減略微增快，無(wú)量綱縱向位移變小。航行體單獨(dú)出水的縱向速度和縱向無(wú)量綱位移大于并聯(lián)出水的，這表明并聯(lián)出水對(duì)超空泡航行體的減阻產(chǎn)生不利影響。

圖9 不同工況下航行體無(wú)量綱縱向位移、縱向速度隨時(shí)間的變化Fig. 9 Variations of dimensionless longitudinal displacement and longitudinal velocity of vehicle with time under different cases

圖10 為不同工況下航行體出水過(guò)程的偏航角變化。從圖中可知，航行體單獨(dú)出水時(shí)，航行體幾乎沒(méi)有發(fā)生偏轉(zhuǎn)，偏航角的變化近似為 0 度附近的一條直線(xiàn)。對(duì)于并聯(lián)出水，并聯(lián)間距△d為2D和3D時(shí)，航行體的偏航角隨時(shí)間的增加先增加后減小，其原因?yàn)楹叫畜w向外側(cè)偏轉(zhuǎn)，當(dāng)偏航角達(dá)到某一臨界值時(shí)，本文工況條件下約為3.4°，航行體的尾部刺穿超空泡內(nèi)測(cè)壁面沾濕，產(chǎn)生了橫向沾濕面作用的均布力，均布力產(chǎn)生向內(nèi)側(cè)的偏轉(zhuǎn)力矩，此力矩大于航行體受到的向外側(cè)偏轉(zhuǎn)的力矩，因此航行體開(kāi)始向內(nèi)側(cè)偏轉(zhuǎn)，對(duì)應(yīng)偏航角減小?！鱠為4D和5D時(shí)，隨著并聯(lián)間距的增大，航行體受到相鄰航行體流場(chǎng)干擾的影響減弱，航行體向外側(cè)的偏航角減小，一直到航行體完全出水時(shí)偏航角還沒(méi)有達(dá)到臨界值，所以航行體只發(fā)生了向外側(cè)偏轉(zhuǎn)。因此為確保航行體出水彈道的穩(wěn)定性，需要合理控制并聯(lián)間距。

圖10 不同工況下航行體出水的偏航角變化Fig. 10 Variations of the yaw angle of vehicles exiting water under different cases

3 結(jié) 語(yǔ)

本文采用數(shù)值模擬方法對(duì)高速航行體并聯(lián)出水問(wèn)題進(jìn)行研究，獲得主要結(jié)論如下：

1）并聯(lián)出水過(guò)程中，航行體周?chē)鲌?chǎng)彼此發(fā)生干擾，出水超空泡內(nèi)側(cè)擴(kuò)張受到抑制，導(dǎo)致超空泡外側(cè)輪廓的曲率要大于內(nèi)側(cè)輪廓。而且當(dāng)并聯(lián)間距較小時(shí)，2 個(gè)出水超空泡的尾部發(fā)生融合。

2）相比航行體單獨(dú)出水，并聯(lián)出水過(guò)程中航行體的縱向速度衰減較快，使航行體的減阻性能小幅下降。

3）航行體并聯(lián)出水受到不對(duì)稱(chēng)的水動(dòng)力作用，航行體的彈道發(fā)生偏轉(zhuǎn)，當(dāng)并聯(lián)間距較小時(shí)，航行體先發(fā)生向外側(cè)偏轉(zhuǎn)之后轉(zhuǎn)向內(nèi)側(cè)偏轉(zhuǎn)，當(dāng)并聯(lián)間距較大時(shí)，航行體只發(fā)生向外側(cè)的偏轉(zhuǎn)。