水下爆炸氣泡射流作用下固支方板動態響應的數值仿真研究

，，，

(1 海軍航空工程學院 研究生管理大隊，山東 煙臺 264001；2海軍工程大學 船舶與動力學院， 武漢 430033)

傳統上將水下爆炸載荷分為兩個階段：沖擊波階段和氣泡階段。沖擊波載荷能夠對艦船結構造成嚴重的局部毀傷，而氣泡載荷主要引起艦船結構的整體毀傷，當氣泡脈動周期與艦船低階振動頻率接近時可能引起艦船結構的“鞭狀響應”(Whipping)[1]，當氣泡尺寸較大時可能造成艦船浮力損失，進而造成整體破壞[2]。藥包在結構近場爆炸時，氣泡受到結構等邊界影響容易形成射流，射流載荷也可能對結構產生毀傷作用。對于水下爆炸氣泡運動及其塌陷形成射流的過程已有很多研究[3-4]，在射流形成的條件、射流的影響因素等問題上取得了很多研究成果。但是國內關于由射流載荷引起的結構毀傷研究尚不多見。研究射流載荷引起的固支方板毀傷的特點、程度等問題對艦船防護設計和水中兵器戰斗部設計具有重要的軍事意義及參考價值。為此，模擬近水面處水下爆炸氣泡的脈動及產生射流的過程，并得到相應的固支方板結構的動態響應。

1 基本理論

1.1 歐拉算法簡介

對于水中爆炸問題，采用歐拉方法求解比較合適。當采用歐拉方法時，節點固定在空間中，由相關節點連接而成的單元僅僅是空間的劃分，分析對象的材料在網格中流動，材料的質量、動量和能量從一個單元流向另一個單元。見圖1。

圖1 Euler求解器中的網格

歐拉材料在單元之間的流動主要用質量守恒、動量守恒和能量守恒來描述。

(1)

(2)

(3)

式中：V——體積；

A——該體積區域的邊界；

n——邊界上的法向量；

u——區域內速度向量。

物質從一個單元傳遞到另一個單元，輸出單元和輸入單元都可以包含多種歐拉材料，每個歐拉單元最多可以有五種材料。

1.2 材料狀態方程

在水下爆炸問題的仿真研究中涉及空氣、水、炸藥及金屬。各種材料本構關系選取是否恰當，直接關系到計算結果的精度。

仿真計算中認為空氣是理想氣體，并且滿足

p=(γ-1)ρ·e

(4)

式中：p、ρ、e——氣體的壓力、密度和比能；

γ——氣體比熱比。

認為炸藥是高度壓縮高壓氣體，也適用式(4)。

假設水是無粘無旋的可壓縮流體。其壓力狀態方程用多項式表示為

(5)

式中：k——體積模量；

ρ、ρ0——總密度和參考密度。

金屬材料的模擬采用雙線性彈塑性本構模型，材料的屈服極限σy描述為

εp

(6)

式中：σ0——靜態屈服強度；

Eh——應變硬化模量；

εp——有效塑性應變。

材料在強沖擊載荷作用下，一般會產生比較明顯的應變率強化效應，特別是對于低碳鋼這種應變率敏感材料。計算中應變率強化模型采用Cowper-Symonds模型，動態屈服強度σd為

(7)

C、P——應變率系數，對于低碳鋼,一般取C=40，P=5。

材料失效模型采用最大塑性應變失效。

2 數值仿真模型

模型區域尺寸為3.0 m×3.0 m×2.4 m，見圖2。其中上面0.4 m為空氣，下面2.0 m為水。在模型中，結構與炸藥有兩種相對關系，炸藥在結構底部和側面，可以很好地代表船底和舷側兩種典型艦船局部結構。炸藥正對固支方板的中心點，底部爆炸爆距分別為0.2 m和0.1 m，側面爆炸爆距分別為0.20、0.15 m。

圖2 炸藥、水面與結構位置關系

水和空氣采用歐拉正六面體單元劃分，單元尺寸為0.05 m，共計172 800個歐拉單元。固支方板指方板在邊界處沿x、y、z三個方向上的位移和彎矩都為零，方板邊長為1.0 m，采用拉格朗日四邊形單元劃分，單元尺寸為0.016 7 m，共計3 600個拉格朗日單元。為了使耦合面成為封閉曲面，且計算過程中耦合面內體積不能為負，建立了9 360個“啞單元”。“啞單元”是虛擬單元，只起到封閉耦合面的作用，不影響結構響應。計算中采用一般耦合算法，耦合面取結構和“啞單元”組成的封閉曲面。

歐拉區域上表面壓力為101.33 kPa，下表面為121.38 kPa。水的密度取1 000 kg/m3，體積模量為2.2 GPa。空氣密度1.21 kg/m3，比熱比取為1.4，比能取210 kJ/kg，初始壓力為100 GPa。TNT炸藥密度為1 580 kg/m3，炸藥重5 g，比能為4 190 kJ/kg。在MSC.Dytran中炸藥可以用JWL或者高壓氣體方程來模擬，若假設藥包為球形，藥包半徑不足0.01 m，比歐拉網格尺寸小太多，在計算中會影響精度，因此在計算中取半徑為0.02 m，按照總能量守恒原則，其密度相應調整為149 kg/m3。按照理想氣體狀態方程計算，炸藥初始壓力為250 MPa。

方板結構模型材料計算采用雙線性強化模型,材料參數參照文獻[5]。彈性模量取為210 GPa，泊松比為0.3，密度為7 800 kg/m3，強化模量為250 MPa，靜態屈服極限為250 MPa，最大失效應變為0.2。

3 仿真結果及分析

3.1 氣泡運動仿真結果

整個計算時間為100 ms，底部爆炸仿真和側面爆炸仿真清晰地再現了氣泡運動并發生射流的現象。圖3和圖4顯示了炸藥爆炸時氣泡膨脹、收縮并發生沖向結構的射流的全過程。底部爆炸和側面爆炸在氣泡膨脹階段基本相似，但是在氣泡收縮階段側面爆炸比底部爆炸更早發生射流。側面爆炸氣泡脈動周期與自由場相比相差不大，底部爆炸時氣泡脈動周期變長。

圖3 底部0.1 m爆炸氣泡運動過程

圖4 側面0.15 m爆炸氣泡運動過程

3.2 射流速度計算結果

仿真發現，氣泡在塌陷形成射流時，流體質點速度很大，最高可以達到200 m/s，但是到達結構處流體質點剩余速度一般在40～50 m/s這一量級。圖5給出了射流發生時，結構附近流場質點速度矢量。從方向看，射流發生時周圍流體都向著板中間方向運動。對于底部和側面爆炸，發生射流時質點運動情況基本類似。

3.3 結構響應計算結果

從結構響應角度來看，沖擊波載荷過后，結構整體有比較大的塑性變形，射流載荷過后，板的中心塑性變形撓度增加較大，對周圍區域塑性變形影響不大。這與以上分析射流發生時流體質點運動特點相吻合。射流時，質點都被氣泡聚集在一起向結構中心運動，引起較大局部變形。圖6是底部0.20 m爆炸時結構等效塑性應變云圖。可以看出，沖擊波載荷過后，中間很大區域形成一個平臺，而射流載荷過后，應變和變形梯度增大。側面爆炸時，情況基本一致。

圖5 射流發生時結構附近流體節點速度矢量

圖6 等效塑性應變云圖

為準確研究結構上不同位置點的變形特點，在結構上取4個典型點A、B、C、D，見圖7a)，其中A為整個結構中點，其余3個點沿中心線等間距分布。

圖7 典型位置變形撓度隨時間變化曲線

圖7顯示結構的變形明顯分為兩個階段，前一階段是由沖擊波載荷引起的，后一階段是由氣泡射流引起的。還可以看出，沖擊波載荷引起中間很大區域形成一個平臺，而射流載荷引起變形梯度增大。各點的塑性變形撓度見表1中。

從總變形來看，底部爆炸要比側面爆炸變形大，但是射流引起的塑性變形在總變形中所占的比例基本相同，在射流發生的部位該比例超過30%；中間點射流引起的變形所占比例比較大，靠近固支邊的位置受射流影響較小，這說明射流的影響區域很小，主要引起射流區域附近結構的變形，這與沖擊波和氣泡脈動壓力載荷引起的塑性變形是不同的，壓力波引起的變形區域比較大。

表1 典型部位塑性變形撓度 %

注：wt、ws、wj分別為總塑性變形、沖擊波引起的塑性變形和射流引起的塑性變形。

4 結論

1) 炸藥在結構底部爆炸氣泡發生射流的時間要晚于同工況下側面爆炸時間，但是引起的結構變形大于側面爆炸。

2) 在射流發生時，周圍的水流都會向結構高速運動，算例中水流速度可以高達200 m/s，但是在向結構運動過程中速度會降低，到達結構位置的速度在40～50 m/s左右。

3) 在近場水下爆炸氣泡發生射流時，射流載荷會對結構產生較強的局部作用，對射流區域外的影響不大，因此射流載荷是很強的局部載荷，可能造成結構的嚴重局部毀傷。

4) 底部和側面爆炸中，在射流發生位置由射流引起的塑性變形都達到總塑性變形的30%以上，所以在近場水下爆炸研究中應該重視射流載荷的作用。

[1] 李玉節，潘建強，李國華.水下爆炸氣泡激起的船體鞭狀運動的試驗研究[J].船舶力學，2001,5(6):75-83.

[2] 朱 錫，方 斌.艦船靜置爆炸氣泡時總縱強度計算方法研究[J].海軍工程大學學報，2007,19(6):6-11.

[3] 張阿漫，姚熊亮.水深和藥量的變化對水下爆炸氣泡射流的影響研究[J].工程力學，2008,25(3):222-229.

[4] 牟金磊，朱 錫，張振華.近自由面水下爆炸氣泡現象的數值仿真研究[J].艦船科學技術，2008,30(4):113-116.

[5] 張振華，朱 錫，馮 剛，等.水下爆炸沖擊波作用下自由環肋圓柱殼動態響應的數值仿真研究[J].振動與沖擊，2005,24(1):45-48.