余 俊，王海坤，劉國振，郝 軼，張得揚(yáng)

（1.中國船舶科學(xué)研究中心，江蘇無錫 214082；2.深海技術(shù)科學(xué)太湖實(shí)驗(yàn)室，江蘇無錫 214082）

0 引 言

隨著現(xiàn)代海戰(zhàn)武器的高速發(fā)展，實(shí)戰(zhàn)條件下艦船遭受的威脅越來越大，對艦船水下沖擊與防護(hù)性能的研究受到了高度的重視。水面艦艇除了會遭受船底爆炸沖擊損傷之外，近水面爆炸也會產(chǎn)生重要的威脅[1]。近水面爆炸現(xiàn)象的本質(zhì)特征屬于多相流體（氣-液-汽）相互作用及其與結(jié)構(gòu)的耦合運(yùn)動過程，其非線性耦合效應(yīng)主要體現(xiàn)在兩個重要方面：首先是爆轟氣體產(chǎn)物、水、空氣等多組分流體之間的對流運(yùn)動,其可歸納為多相可壓縮流體界面捕捉問題；其次是自由面的存在使得流場中會出現(xiàn)稀疏波的傳播，水中會出現(xiàn)典型空化現(xiàn)象,可以歸結(jié)為相變問題。要想準(zhǔn)確理解和掌握近水面爆炸過程中的多流體非線性耦合效應(yīng)，需要同時重點(diǎn)處理好多相可壓縮流體界面捕捉和相變問題。

目前國內(nèi)外已經(jīng)對近水面爆炸現(xiàn)象進(jìn)行了大量的研究[2-4]。在實(shí)驗(yàn)研究方面，Marcus等[5]開展了水下近自由面爆炸試驗(yàn)研究，發(fā)現(xiàn)水下爆炸產(chǎn)生的片空化區(qū)域在閉合潰滅時會產(chǎn)生較大的壓力，獲得了近水面爆炸過程中包括空化效應(yīng)的沖擊波和氣泡階段的壓力曲線;Brett等[6]開展了一系列水下圓柱殼近距離爆炸試驗(yàn)，對沖擊波和空化聯(lián)合加載下圓柱殼的變形進(jìn)行了測量分析，發(fā)現(xiàn)空化載荷的幅值與沖擊波幅值相當(dāng)，而且空化潰滅載荷造成的殼體變形要超過整個沖擊波階段產(chǎn)生的總變形；Kleine等[7]利用紋影法拍攝到了近水面爆炸試驗(yàn)過程中液體空化的密度變化過程圖象，較為清晰地展示了爆炸氣泡膨脹以及空化域內(nèi)部演化圖像；Cui 等[8]開展了一系列不同邊界條件下近水面小藥量藥包爆炸試驗(yàn)，獲得了爆炸氣泡與自由面運(yùn)動的完整過程。由于實(shí)驗(yàn)研究獲得的數(shù)據(jù)有限，在穩(wěn)定性和重復(fù)性上難以滿足要求，難以獲得流場的精細(xì)特征。數(shù)值模擬能夠較為全面地展示流場的精細(xì)特征，使得其在近水面爆炸研究中具有一定的優(yōu)勢。這里著重介紹基于多相可壓縮流體計(jì)算模擬的研究情況，對于基于聲學(xué)單元和邊界單元的數(shù)值方法不作討論。自從Aanhold 等[9]提出cut-off 空化模型之后，由于其使用簡單方便，被廣泛應(yīng)用在水下近自由面或近結(jié)構(gòu)爆炸分析中。Wardlaw等[10]利用cut-off模型分析了充水圓筒中心藥包爆炸過程中沖擊波在圓筒內(nèi)部產(chǎn)生的局部空化現(xiàn)象，獲得了空化潰滅載荷的典型特征；Shukla[11]、Yu[12]等采用cut-off 空化模型與五方程相結(jié)合的方法對帶有空化效應(yīng)的近水面爆炸多流體運(yùn)動問題進(jìn)行了研究；Wu 等[13]采用cut-off模型與流體局部間斷迦流金（LDG）法相結(jié)合分析了水下遠(yuǎn)場沖擊波與自由面作用產(chǎn)生空化的過程。為了克服cut-off模型無法描述空化域內(nèi)液相和汽相組分的問題，Liu[14]、Schmidt[15]和Xie[16]等相繼提出了isentropic 模型、Schmidt 模型和Modified Schmidt模型。這些模型本質(zhì)上都屬于單流體空化模型，在流體壓力降低到飽和態(tài)壓力之前均視為純流體，當(dāng)降到飽和態(tài)壓力之后，各種模型的具體差別體現(xiàn)在空化域內(nèi)汽液混合流體狀態(tài)方程的構(gòu)造方法不同。除了上述單流體空化模型之外，雙流體模型也在不斷發(fā)展，主要以Chiapolino[17-18]、Saurel[19]、Pelanti[20-21]等提出的四方程、六方程等模型為主，該模型通過持續(xù)追蹤空化質(zhì)量分?jǐn)?shù)的演化來捕捉空化的產(chǎn)生、發(fā)展與潰滅過程。雙流體模型則將初始流體視為液相和汽相的某種混合，在流場發(fā)生改變時（如沖擊波傳播和氣泡運(yùn)動等包含間斷時）液-汽兩相之間會發(fā)生對流運(yùn)動，同時還伴隨有質(zhì)量和熱量交換的相變過程，空化被視為流體區(qū)域內(nèi)部蒸汽含量不斷累積的過程。因此雙流體模型描述的空化演化過程更符合物理本質(zhì)和思維邏輯，但雙流體模型相比單流體模型，雙流體計(jì)算過程非常復(fù)雜繁瑣，導(dǎo)致目前水下爆炸領(lǐng)域內(nèi)絕大部分的空化研究都采用單流體模型。

本文擬采用Chiapolino 等[18]提出的基于多相可壓縮流體四方程的雙流體空化模型來重點(diǎn)研究近水面爆炸過程中的多相流運(yùn)動。首先對計(jì)算模型的控制方程及其數(shù)值方法進(jìn)行了簡要介紹，然后采用一維激波管問題對計(jì)算模型中的相變轉(zhuǎn)換以及激波捕捉功能進(jìn)行了考核，最后針對近水面爆炸現(xiàn)象，將計(jì)算結(jié)果與試驗(yàn)結(jié)果進(jìn)行對比，并重點(diǎn)分析空化域內(nèi)部的典型物理特性。

1 計(jì)算模型簡介

1.1 控制方程

不考慮粘性和熱傳導(dǎo)效應(yīng)的多相可壓縮流體運(yùn)動方程可表示為[18]

式中，ρ、u、p和E分別代表混合流體的密度、速度矢量、壓力和單位質(zhì)量總能，E=e+0.5u·u，這里e為單位質(zhì)量內(nèi)能；Yk代表混合流體中各相介質(zhì)的質(zhì)量分?jǐn)?shù)?？刂品匠蹋?）中參數(shù)k代表了混合流體中各相成分，為了計(jì)算模型中表達(dá)簡便起見，這里約定如下：k=1代表液相成分；k=2代表蒸汽相成分；k=3,…，N代表其它不發(fā)生相變的氣相或液相成分。本文中令k=3 代表空氣（不可冷凝），k=4 代表爆轟氣體產(chǎn)物。

在流體計(jì)算模型中普遍采用的是SG（Stiffed-Gas）狀態(tài)方程，該方程能夠描述多種液相或者氣相流體的運(yùn)動特性。近年來NASG（Nobel-Able Stiffed-Gas）狀態(tài)方程逐漸被采納，該方程修正了SG方程對液相介質(zhì)中分子排斥力效應(yīng)描述的不足[22]。該方程形式為

式中，?k=1/ρk、Tk、gk、ck分別代表k相流體介質(zhì)的比體積、溫度、吉布斯（Gibbs）自由能與聲速，參數(shù)γk、、Ck、qk，qk'、bk等可以通過流體熱動力學(xué)試驗(yàn)中與溫度相關(guān)的飽和態(tài)壓力、各相比體積、各相的焓等數(shù)據(jù)來擬合得到[22]。在飽和平衡態(tài)下液相與汽相的Gibbs 自由能相等，即g1=g2，化簡可得飽和態(tài)下溫度與壓力的關(guān)系式：

式中下標(biāo)“sat”代表飽和態(tài)，參數(shù)As、Bs、Cs、Ds、Es分別為

式中Cp,k表示第k相流體定壓比熱容。液態(tài)水、水蒸汽以及空氣的NASG 狀態(tài)方程參數(shù)如表1所示，其中W為物質(zhì)的摩爾質(zhì)量。本文算例涉及到水中空化問題所采用的計(jì)算參數(shù)均采用表1中的參數(shù)[22]。

表1 液態(tài)水、水蒸氣和空氣的NASG狀態(tài)方程參數(shù)Tab.1 NASG coefficients for liquid,vapor and air

1.2 相變轉(zhuǎn)換模型

控制方程（1）并未考慮到各相之間發(fā)生的相變過程，對于流體中發(fā)生液相與其對應(yīng)汽相之間的物質(zhì)與熱量轉(zhuǎn)換，可以參考化學(xué)反應(yīng)過程中判斷反應(yīng)方向的判據(jù)，認(rèn)為在達(dá)到平衡態(tài)之后液相及其蒸汽相之間的吉布斯自由能相等。同時結(jié)合控制方程（1）的假定，認(rèn)為系統(tǒng)達(dá)到平衡態(tài)后滿足如下關(guān)系式[18]：

將公式（2）和（3）代入，得

這里ppartial代表氣體混合物中蒸汽相的分壓，由公式可知其與蒸汽相的摩爾質(zhì)量分?jǐn)?shù)成正比[18]。非線性方程（6）～（8）可以采用迭代方法進(jìn)行求解，如Newton-Raphson方法。

1.3 數(shù)值離散

控制方程（1）與狀態(tài)方程（2）以及空化相變方程（6）～（8）共同構(gòu)成了考慮相變效應(yīng)的多相可壓縮流體控制方程組，可以利用相變松弛法則進(jìn)行分步求解[23]。第一步求解齊次雙曲型方程組（1），需要結(jié)合狀態(tài)方程（2），利用二階MUSCL-Hancock方法以及HLLC近似黎曼求解器。通過該步獲得中間變量（v,e,p,T,Yk）；第二步采用Newton-Raphson 迭代方法求得相變轉(zhuǎn)換平衡態(tài)時的狀態(tài)量（p*,T*,Yk=1,2）。

2 計(jì)算與討論

2.1 一維激波管測試與驗(yàn)證

本算例主要用來與Chiapolino 等相同計(jì)算工況[18]進(jìn)行對比驗(yàn)證。其中整個流場的蒸汽相初始質(zhì)量分?jǐn)?shù)設(shè)置為0.2，液相與氣相的質(zhì)量分?jǐn)?shù)分別為0.1和0.7，計(jì)算域?yàn)閇0,1]m。流場中密度和壓力的初始間斷面位于x=0.5 m，左、右側(cè)壓力分別為2×105Pa和105Pa，兩側(cè)密度分別為1.94 kg/m3、1.02 kg/m3。

圖1 為t=1 ms 時刻流場中密度、壓力、速度、溫度、液相與蒸汽相質(zhì)量分?jǐn)?shù)分布曲線。其中P-T 和no P-T 曲線分別代表是否考慮相變轉(zhuǎn)換時的計(jì)算工況，initial代表初始條件，Chiapolino 代表對比的數(shù)據(jù)（網(wǎng)格數(shù)為100）。同時為了考察計(jì)算模型的網(wǎng)格收斂性，分別采用100、200 和400 個網(wǎng)格單元均勻劃分計(jì)算域。由密度和壓力分布曲線可知，考慮液相和汽相之間的相變轉(zhuǎn)換之后，向右傳播的壓縮波和向左傳播的稀疏波的波速都稍有下降。由溫度曲線可知，向右傳播的壓縮波波后溫度在兩種工況下都有所升高，但波后穩(wěn)定區(qū)的溫度在考慮相變轉(zhuǎn)換時為346.3K,而不考慮相變時為362.6 K，差異較大。同理，向左傳播的稀疏波波后穩(wěn)定區(qū)的溫度都有所下降，但是考慮相變轉(zhuǎn)換后為344.7 K，而不考慮時則為329.4 K。這是由于在壓縮波傳播的波陣面以及稀疏波傳播過程中液相與汽相發(fā)生了質(zhì)量和熱量交換，在達(dá)到平衡態(tài)的過程中溫度存在變換。而汽-液兩相的質(zhì)量轉(zhuǎn)換可以參考質(zhì)量分?jǐn)?shù)曲線得知，壓縮波過后流體中液相發(fā)生了蒸發(fā)，使得汽相質(zhì)量分?jǐn)?shù)增加。該計(jì)算結(jié)果與通常理解的沖擊波作用后由于壓力升高容易引起液相發(fā)生冷凝的認(rèn)識相矛盾，這里需要結(jié)合溫度變化進(jìn)行聯(lián)合分析。由于考慮了相變轉(zhuǎn)換，沖擊波作用后波后溫度由初始的337.5 K 上升到346.3 K,雖然波后壓力也由初始的105Pa上升到1.4×105Pa，但是從飽和態(tài)曲線P-T（公式（3））來分析會發(fā)現(xiàn)溫度升高引起的蒸發(fā)效應(yīng)要大于壓力升高帶來的冷凝效應(yīng)，最終效果就是壓縮波后處于蒸發(fā)狀態(tài)。同理，向左傳播稀疏波的最終效果是波后處于冷凝狀態(tài)。由圖1 可知，隨著網(wǎng)格數(shù)的增加，計(jì)算結(jié)果的收斂性較好，本文的計(jì)算結(jié)果與Chiapolino的結(jié)果吻合較好。

圖1 1 ms時刻流場各物理量分布曲線圖Fig.1 Diagram of different variables in fluid at t=1 ms

2.2 水下近自由面爆炸現(xiàn)象

本節(jié)利用上面介紹的計(jì)算模型與Cui 等[8]的試驗(yàn)觀察結(jié)果進(jìn)行對比。該試驗(yàn)在邊長2 m 的立方體水槽內(nèi)開展，采用等效于5.2 g TNT 的藥包起爆，藥包中心位于水下0.13 m 處。藥包采用等效爆轟模0.0型09近m似，爆處轟理氣，根體據(jù)產(chǎn)G物e內(nèi)e r-部H密u n度ter和理壓論力模分型別[24]以為及16數(shù)06值kg經(jīng)/m驗(yàn)3和，初10始9P球a，形采爆用轟理氣想體氣產(chǎn)體物狀半態(tài)徑方設(shè)程置，γ為=1R.80=，CV=695，W=290?？諝?、水和爆炸氣泡內(nèi)的初始蒸汽質(zhì)量分?jǐn)?shù)分別為10-9、10-8和10-6?？諝?、水和爆炸氣泡的初始溫度分別設(shè)置為295 K、295 K 和1120 K,其中空氣和爆轟氣體產(chǎn)物處于過熱狀態(tài)，水處于飽和狀態(tài)?？諝狻⑺捅馀輧?nèi)的初始蒸汽體積分?jǐn)?shù)分別為1.56×10-9、1.39×10-5和8.07×10-7。為了節(jié)省計(jì)算資源，采用二維軸對稱模型進(jìn)行模擬，如圖2 所示。其中y軸為對稱旋轉(zhuǎn)軸，y=0 為初始水平面位置，藥包位于(0,-0.13)m處。計(jì)算域采用非均勻結(jié)構(gòu)性網(wǎng)格劃分，如圖2 所示（圖中顯示的網(wǎng)格密度為實(shí)際的1/2）。其中區(qū)域[0, 0.4]×[-0.6, 0.2] m2內(nèi)采用均勻網(wǎng)格，dx=dy=0.001 m，區(qū)域外采用等比數(shù)列格式的遞增步長網(wǎng)格，等比為1.02。整個計(jì)算域網(wǎng)格為650×1200，計(jì)算域大小為[0,7.34]×[-7.5,1.1]m2。y軸邊界設(shè)置為對稱條件，其它如空氣層上方、水域下方以及右側(cè)邊界均設(shè)置為無反射條件，CFL設(shè)置為0.5。

圖2 近水面爆炸二維軸對稱計(jì)算模型與網(wǎng)格劃分示意圖Fig.2 Schematic of 2D axisymmetric model and mesh size of underwater explosion near free surface

圖3 顯示了在考慮相變條件下t=0.084 ms、0.168 ms、0.334 ms、0.5 ms、0.666 ms、1.0 ms、4.176 ms 和6.176 ms 時刻流場中的密度、壓力和蒸汽體積分?jǐn)?shù)分布云圖，其中所有圖片采用統(tǒng)一大小的觀察窗，尺寸為[-0.6,0.6]×[-0.7,0.35]m2。由于圖3 中的各物理量的內(nèi)部分布跨度非常大，未單獨(dú)列出各云圖的刻度線，壓力云圖中的黑色虛線表示對應(yīng)時刻爆炸氣泡界面的位置。0.084 ms 時刻，初始沖擊波剛到達(dá)自由面不久，其產(chǎn)生的反射稀疏波使得自由面下方附近出現(xiàn)蒸汽聚集，圖中對應(yīng)時刻蒸汽體積分?jǐn)?shù)云圖顯示的局部高亮區(qū)域內(nèi)數(shù)值范圍為0.4‰～2.5‰。隨著稀疏波不斷向下傳播，0.168 ms 時刻稀疏波已經(jīng)達(dá)到爆炸氣泡上壁面，并反射壓縮波，可以發(fā)現(xiàn)此時爆炸氣泡表面附近的局部范圍處于相對高壓區(qū)。隨著稀疏波在爆炸氣泡表面反射成壓縮波以及自由面不斷反射稀疏波，旋轉(zhuǎn)軸上的蒸汽含量不斷降低（對應(yīng)空泡潰滅過程），而旋轉(zhuǎn)軸外側(cè)的蒸汽含量不斷增加（對應(yīng)空化產(chǎn)生過程），如圖3中的0.334 ms和0.500 ms所示。隨著空泡的不斷產(chǎn)生和潰滅，空泡中心區(qū)域逐漸變薄，并向外擴(kuò)張，如圖3中的0.666 ms和1.0 ms所示。在4.167 ms和6.167 ms時刻空泡基本上完全消失，但爆炸氣泡一直處于向外膨脹過程，并不斷抬高水面位置。在整個過程中空化由開始的單連通域不斷演化為多連通域，顯示出渦環(huán)形狀，渦環(huán)逐漸向外擴(kuò)展，并逐漸變薄直至最終完全消失。

圖3 典型時刻流場的密度、壓力和蒸汽體積分?jǐn)?shù)分布云圖Fig.3 Contour plots of the density,pressure and vapor volume fractions at typical times

在上述整個時間段內(nèi)，自由面上方附近的空氣中蒸汽含量一直比較高，并且包含蒸汽的區(qū)域不斷擴(kuò)展。這是由于水中初始沖擊波作用自由面后向空氣透射壓縮波、水中的液相和汽相同時向空氣中擴(kuò)散的結(jié)果，在擴(kuò)散過程中仍然存在液-汽兩相之間的相變轉(zhuǎn)換。這是因?yàn)樗嫔戏匠跏伎諝庵惺遣话合嗨?，而且水蒸汽含量極低（10-9），而圖3中蒸汽體積分?jǐn)?shù)云圖中紅色區(qū)域代表自由面上方，其內(nèi)部的液相水的質(zhì)量分?jǐn)?shù)范圍在0.6～0.8，水蒸汽質(zhì)量分?jǐn)?shù)范圍在0.1～0.3，剩下的為空氣。

文獻(xiàn)[8]中只記錄了通過觀察窗獲得的部分流場圖像，如圖4 中第一行所示。試驗(yàn)獲得的圖片的中間縱剖面尺寸與數(shù)值模擬中的[-0.4,0.4]×[-0.68,0.2]m2基本相當(dāng)。在第一行的試驗(yàn)圖片中左下角兩位數(shù)值代表圖片序號，右下角為記錄的時刻（單位為ms），其中0.167 ms代表裝藥起爆時刻。由圖4可知，小藥量裝藥在近水面爆炸時在水面下方產(chǎn)生的空化域形態(tài)與“云空化”非常相似，水中的空泡體積非常小，形成云帶狀，其發(fā)展過程中能夠觀察到明顯的渦環(huán)形狀。為了與試驗(yàn)進(jìn)行對比，需要對“云空化”區(qū)域內(nèi)的蒸汽含量進(jìn)行判斷，選取蒸汽體積分?jǐn)?shù)是一個比較合適的對象。這里假定當(dāng)流體單元中的蒸汽體積分?jǐn)?shù)大于0.5‰時就標(biāo)記為空化單元，單元所在的區(qū)域即為空化域。有關(guān)空化域內(nèi)蒸汽含量的判據(jù)目前仍然沒有統(tǒng)一的標(biāo)準(zhǔn)，這里選取的0.5‰只是數(shù)值計(jì)算上的一個參考。將圖3 中蒸汽體積分?jǐn)?shù)云圖轉(zhuǎn)換為空化域紋影圖，如圖4 第二行所示。其中紅色區(qū)域?yàn)樗驴栈?，黃色區(qū)域?yàn)樽杂擅嫔戏降母吆空羝袇^(qū)域，綠色為空氣，藍(lán)色為水域，淺藍(lán)色為爆炸氣泡。經(jīng)過對比可知，計(jì)算獲得的水中空化域與試驗(yàn)觀察到的“云空化”范圍基本相當(dāng)，兩者都清晰地展示了空化域形成渦環(huán)，并逐漸向外擴(kuò)張和變薄的過程，最終完全潰滅。

圖4 典型時刻的試驗(yàn)與計(jì)算結(jié)果對比Fig.4 Comparison of experimental and numerical results at typical times

圖5 顯示了空化域內(nèi)的最大、最小和平均溫度變化以及空化域總體積時程曲線。由圖5 可知空泡總體積的膨脹過程與收縮過程的時間基本相當(dāng)，整個空化運(yùn)動周期約為1.14 ms。空化域膨脹和收縮過程中內(nèi)部最高溫度變化范圍為294.6～309.7 K，最低溫度變化范圍為283.7～294.6 K。膨脹到最大體積時空化域內(nèi)部的最高和最低溫度分別為309.7 K 和283.7 K，但在空化域變化過程中其內(nèi)部體積平均溫度基本恒定為294.6 K，與室溫相當(dāng)。

圖6 為流場中三個坐標(biāo)點(diǎn)處壓力時程曲線圖，分別為水下0.05 m 處距離爆心不同水平距離處的測點(diǎn)。由圖可知，三個測點(diǎn)處的空化潰滅過程中最大壓力分別為2.15 MPa、0.98 MPa 和0.92 MPa。圖7 為空化域內(nèi)最大、最小和平均壓力變化時程曲線。由圖可知，空化域內(nèi)最大壓力在體積膨脹階段的大部分時刻保持在0.018 MPa 左右，而在收縮階段則主要在0.018 MPa 上方振蕩，最大值達(dá)到0.029 MPa。空化域內(nèi)最小壓力在中間的絕大部分時間內(nèi)變化較為平緩，在3200～4800 Pa 范圍內(nèi)變動，只是在空化剛開始產(chǎn)生和最后潰滅階段變化梯度較大。空化域內(nèi)體積平均壓力變化范圍為9000～18 000 Pa，其變化趨勢與最小壓力變化趨勢基本保持一致。由此可見，空化域內(nèi)壓力并不是維持在飽和態(tài)壓力不變，其具有分布不均勻性、變化梯度大等特點(diǎn)。同時結(jié)合圖5可知其內(nèi)部溫度的變化范圍在283.7～309.7 K，空化域內(nèi)部發(fā)生較為明顯的相變轉(zhuǎn)換現(xiàn)象。

圖5 空化域內(nèi)最大、最小和平均溫度以及空化域體積的時程曲線Fig.5 Time history of the maximun,minimun,average temperatures and cavitation domain volumes in cavitation domain

圖6 三個測點(diǎn)處壓力時程曲線圖Fig.6 Time history of pressure at three measuring points

圖7 空化域內(nèi)最大、最小和平均壓力時程曲線Fig.7 Time history of the maximum,minimum and average pressures in cavitation domain

圖8 為空化域內(nèi)的蒸汽相體積分?jǐn)?shù)的最大和平均值變化時程曲線，空化域內(nèi)蒸汽相體積分?jǐn)?shù)的最小值為前面設(shè)定參考值的0.5‰。由圖8 可知，空化域內(nèi)的蒸汽相體積分?jǐn)?shù)在0.5‰～173‰，其峰值并未出現(xiàn)在空化體積膨脹到最大的時刻。盡管空化整個運(yùn)動過程中蒸汽相最大體積分?jǐn)?shù)并不低，但是空化域內(nèi)部整體體積平均值變化區(qū)間在0.5‰～10‰，說明蒸汽的平均含量比較低，這與前面介紹的實(shí)驗(yàn)中觀察到的“云空化”較為一致，同時也說明了本文計(jì)算模型具有一定的合理性和精確性。圖9為空化域內(nèi)密度的最大、最小和平均值變化時程曲線。由圖9可知，空化域內(nèi)的最大密度變化范圍為1051～1061 kg/m3，體積平均密度變化范圍為1039～1051 kg/m3，最小密度變化范圍860～1061 kg/m3，其峰值也未出現(xiàn)在空化體積膨脹到最大的時刻。結(jié)合圖8可知，空化域內(nèi)蒸汽相的含量總體上非常小，對于流體密度降低的程度有限，空化域內(nèi)總體體積平均密度變化較小。

圖8 空化域內(nèi)蒸汽體積分?jǐn)?shù)?2最大值、平均值變化曲線Fig.8 Time history of the maximum and average vapor volume fractions ?2 in cavitation domain

圖9 空化域內(nèi)密度的最大值、最小值和平均值的時程曲線Fig.9 Time history of the maximum,minimum and average densities in cavitation domain

3 結(jié) 論

針對近水面爆炸過程中出現(xiàn)的爆轟氣體產(chǎn)物、水以及空氣等多流體復(fù)雜運(yùn)動，以及各種流體內(nèi)部可能出現(xiàn)的液相水和汽相水之間的相變轉(zhuǎn)換等現(xiàn)象，本文采用了基于相變效應(yīng)的多流體運(yùn)動計(jì)算模型進(jìn)行處理，獲得了如下結(jié)論：

（1）利用一維激波管問題對計(jì)算模型精確性和收斂性測試表明，本文的模型計(jì)算結(jié)果與Chiapolino的結(jié)果吻合較好，能夠精確地捕捉?jīng)_擊波與稀疏波的傳播過程，展示了沖擊波作用之后流體蒸發(fā)以及稀疏波作用之后流體冷凝的過程，最終的狀態(tài)往往是兩種過程的綜合作用的結(jié)果；

（2）近水面爆炸過程中爆炸氣泡運(yùn)動過程的計(jì)算結(jié)果與試驗(yàn)結(jié)果符合較好，說明本文建立的多流體耦合計(jì)算模型能夠捕捉多種界面之間的復(fù)雜運(yùn)動過程。

（3）多流體計(jì)算模型中引入了相變轉(zhuǎn)換，可獲得近水面爆炸過程中水面下方的液體內(nèi)部發(fā)生的空化現(xiàn)象；采用蒸汽體積分?jǐn)?shù)0.5‰作為空化域識別的判據(jù)，數(shù)值模擬獲得的空化域范圍與試驗(yàn)結(jié)果較為一致，空化域都是由開始的單連通域演化為渦環(huán)形態(tài)，并向外擴(kuò)展直至完全潰滅。

（4）數(shù)值模擬結(jié)果發(fā)現(xiàn)，空化域在運(yùn)動過程中其內(nèi)部蒸汽的體積分?jǐn)?shù)范圍為0.5‰～173‰，但內(nèi)部平均值的范圍僅為0.5‰～10‰，說明蒸汽的平均體積含量比較低，這與實(shí)驗(yàn)中觀察到的“云空化”現(xiàn)象具有一致性。同時該現(xiàn)象也可以從空化域內(nèi)部密度的平均體積含量來得到印證，其變化范圍為1039～1051 kg/m3，說明空化域內(nèi)部平均密度與常態(tài)下液相水的密度非常接近。另外，空化域內(nèi)的壓力變化范圍為3200～29 000 Pa，溫度變化范圍為283.7～309.7 K，說明空化域內(nèi)部的壓力和溫度一直處于動態(tài)變化之中，并不是維持在某個恒定的飽和態(tài)不變。

近水面爆炸過程中除了存在多流體之間的非線性耦合效應(yīng)，還會與周圍結(jié)構(gòu)發(fā)生流固耦合作用，其過程更加復(fù)雜，將在后續(xù)計(jì)劃中研究。