水下爆炸沖擊波作用下船體剛體動態響應特性

王軍，郭君，孫豐，楊棣

(哈爾濱工程大學 船舶工程學院，黑龍江 哈爾濱 150001)

艦船在水下爆炸環境中主要發生2種動態響應:剛體運動及自身結構彈塑性變形破壞。若船體的剛體運動不能進行正確分析，則后續艦船的彈塑性動態響應計算及強度評估等工作便無法準確進行［1-2］。關于水下爆炸沖擊波作用下船體動態響應的研究開展了很多，但大部分是結構的彈塑性響應分析［3-5］，船體的剛體運動響應的影響一般被忽略。李海濤等［6］提出了一種計算船體梁在球面沖擊波作用下發生整體運動響應的理論方法，是針對試驗模型的船體梁開展的研究，炸藥在船體結構中部正下方，將球面波看作是移動的平面波，認為船體梁橫向截面受到的沖擊波壓力為均勻分布。張弩等［7］研究了船體在水下爆炸氣泡作用下的剛體運動對其鞭狀運動的影響，并未考慮沖擊波階段的影響。

針對目前研究進展，考慮到實船可能受到任意位置水下武器的攻擊，同時沖擊波是以球面沖擊波的形式作用在船體上，船體不同位置受到作用力的時間不同。為此在文中，首先分析船體在水下爆炸沖擊波作用下的剛體運動，即船體在水中爆炸發生后的起伏運動，然后將剛體運動和彈塑性振動分開研究。建立的理論模型以球面沖擊波的形式進行加載，并考慮沖擊波的反射、流場和船體的流固耦合作用，分析藥包在水下任意位置船體的剛體動態響應。之后通過數值模擬的方法計算出沖擊波入射能，得出船體剛體運動的動能及船體結構彈塑性變形等吸收的能量。此過程應用動量定理進行求解，從而間接求出船體振動變形等吸收的能量，進而對船體的彈塑性振動變形等進行分析。

1 船體垂向剛體運動計算模型

在計算艦船受到水下爆炸沖擊波的作用時，由于不同的船體外殼其計算模型不同，為理論計算帶來很大困難。文獻［9］中給出了以流固耦合的Taylor平板理論計算船體動態響應的合理性，因在已有理論基礎上，應用Taylor平板來代替復雜形狀的船體外殼，分析船體在沖擊波作用下的整體垂向運動。

對于鑄裝TNT球形藥包，空間任意一點爆炸沖擊波的壓力隨時間和位置的變化為

式中:Pm為峰值壓力;θ為衰減時間常數，對于板上不同的位置而言，其距藥包的距離不同，所以爆距不同。設藥包位置為O(xW，yW，zW)，水中聲速為c，板上離藥包最近點為A(x0，y0，z0)，如圖1所示，則沖擊波到達平板的時間:

以沖擊波開始作用平板時記為初始時刻，以圖1所示坐標系進行計算，則對于板上任意一點B(x，y，z)，相對于A點受到壓力有一時間延遲:

在td時刻之前，B點并未受到水下爆炸沖擊波壓力，在此定義一個便于計算的輔助函數:

此時平板任意一點受到的球面沖擊波入射壓力可以表示為

得到沖擊波入射壓力后，計算沖擊波反射壓力時采用流固耦合的Taylor平板理論，假定入射波全部反射，并且這種假定相當近似于緊接沖擊后的實際情況［9］。考慮波以與板面法向成α角度沖到鋼板上，如圖2所示，對于空氣背襯板，流固耦合界面處的壓力可表示為

式中:Pr(x，y，t)為反射壓力，應用Taylor平板理論求解流固耦合界面處的壓力為

圖1 球面沖擊波作用于剛性平板Fig.1 The spherical shock wave acting on a rigid plate

圖2 平板和流體相互作用Fig.2 Interaction between plate and fluid

整個平板受到球面沖擊波的作用力隨時間的變化可表示為

入射角α(x，y)也是隨平板上作用點變化而發生變化的量，用坐標可表示為

Taylor平板理論沒有正確模擬流體慣性(水動力質量)的影響，在此研究平板的垂向運動，需考慮附加質量的影響。為得到較準確的理論計算結果，不使用空氣背襯單面浸水的平板附加質量計算公式，而直接利用流固耦合動力學中的邊界元法［10］計算船體的垂向平動附加質量，計算公式為

式中:φ為速度勢，S為所有邊界面，ρf為流場密度。則由動量定理，平板所受垂向力與自身垂向速度關系為

式中:m為船體自身質量。聯立式(8)～(11)，則平板剛體運動速度隨時間變化為

上式含有隨位置和時間變化無積分性質的輔助函數δ(x，y，t)，同時等式兩邊都有要求解的未知函數v(t)，因此無法得到解析解。但將平板面積進行微分，計算時間離散，通過迭代方法，可以方便進行求解。因平板在受到爆炸沖擊波之前，垂向速度為零，初始沖量為零，所以求解思路可表示為

式中:大括號內表示平板整體在某時刻受到的壓力，Ii表示沖量對時間的積分值，Δx、Δy、Δt分別為板長度、寬度和計算時間的微元，根據求解精度及計算時間的需要可進行設置。

2 理論模型的試驗驗證

2.1 艦船實尺度艙段試驗驗證

為驗證本文理論計算方法的合理性，通過某艦船艙段水下爆炸試驗的高速攝影試驗資料，得到試驗艙段在水下爆炸沖擊波作用下的最大速度，將試驗值與理論計算值進行比較。試驗艙段模型長度方向取某艦船的5倍跨長L，寬度方向取船半寬B，高度方向取四甲板高度H，具體結構如圖3所示。

圖4為試驗時艙段在水下爆炸響應時的圖片，試驗過程中設置3個工況，如圖5所示，都采用8 kg的TNT炸藥，藥包位于艙段中部下方，具體位置如表1所示。通過高速攝影技術，可測出艙段在整個水下爆炸過程中的剛體位移響應，將位移微分便可得到艙段的速度時歷響應。

圖3 艙段結構尺寸Fig.3 Structural size of section model

圖4 水下爆炸試驗Fig.4 The test of underwater exp losion

圖5 工況設置Fig.5 Structural size of section model

表1 試驗工況Table 1 Test conditions

不同工況得到的試驗值和理論值的速度最大值如表2所示，爆距越小，艙段的剛體運動速度越大，總體上理論值和試驗值比較接近，誤差控制在10%以內，充分說明了本文理論計算方法在預報船體水下爆炸剛性運動的有效性。其中工況1的誤差最大，可能由于爆距較近時在試驗過程中模型發生了一定程度的旋轉，增大了試驗值測量的誤差。

以工況1為例，計算出的艙段整體受到的沖擊波垂向壓力如圖6，在求解空氣背襯板的反射壓力時應用了Taylor平板理論，此理論沒有考慮氣泡脈動、氣泡射流等影響，但已經能成功模擬經受水下爆炸時船體結構在早期運動的動態響應［9］。本文在此不考慮氣泡脈動對艦船剛體運動的影響，以此求出艙段剛體運動速度如圖7所示，出現最大剛體運動速度的時間大約為1.5 ms，在8 ms時速度趨近于0。

表2 試驗值與理論值對比Table 2 Comparison of experimental and theoretical results

圖6 艙段受到沖擊波壓力隨時間變化曲線Fig.6 Variation curve of the section receiving the shock wave pressure with time

圖7 艙段剛體運動速度隨時間變化曲線Fig.7 Variation curve of the section velocity of rigid movem ent with time

2.2 實船試驗驗證

為驗證本文方法應用于實船在水下爆炸沖擊環境中計算剛體運動速度的有效性，利用文獻［11］中的實船水下爆炸試驗數據進行對比。試驗船體長54 m，寬8.5 m，排水量570 t。應用本文方法計算得到的理論值為2.07 m/s，而實船試驗值為2.17 m/s，相對誤差5%。

實船的受力和剛體速度響應變化趨勢和艙段類似，但實船受到的沖擊波壓力比艙段大很多，作用時間也長，同時最大速度出現的時間比艙段晚，以沖擊波開始作用船體為初始時刻，艙段試驗的工況1出現最大剛體運動速度的時間大約為1.5 ms，而實船約為4 ms。

由艙段和實船的水下爆炸試驗可知，不論是艦船艙段還是整船試驗的剛性運動速度，以上理論分析得到的理論值與試驗值的相對誤差，都控制在10%以內，充分說明了本文理論計算方法在預報船體水下爆炸剛性運動的有效性。

3 剛體運動對艦船總體響應影響

上文實尺度艙段試驗的模型取材于某艦船，擬通過局部結構試驗來考核全船的局部強度，研究艦船在水下爆炸沖擊環境中的響應。下面以此艦船為研究對象，分析在不同工況的水下爆炸環境中艦船剛體速度響應對艦船總體響應的影響，從而為局部結構試驗設計提供理論指導。

3.1 剛體運動速度對艦船總體速度影響

艦船水下爆炸節點速度響應可認為是剛體運動速度和彈塑性變形響應速度的疊加，為把剛體運動速度從節點速度響應的結果分離出去，需要得到在這個過程中剛體運動速度。

現根據艦船各站的剖面慣性矩和分布質量，把整個艦船等效成變截面梁，船體共分20站，船體中部較重，兩端質量相對較小，將水下爆炸的沖擊波載荷根據不同分段受到不同的沖擊波作用力分別加載在變截面梁對應節點上，這樣分析的問題就成為變截面梁在沖擊載荷作用下的動態響應求解問題。此時可用多體系統的傳遞矩陣法或有限元法進行分析，傳遞矩陣法編程較復雜，且求解效率上不如有限元方法，同時有限元在求解變截面梁的動響應問題已相當成熟，求解精度也較高，所以在此選用有限元方法對變截面梁在沖擊載荷作用下的動態響應問題進行分析。在對實船進行考核時，采用1 t的TNT炸藥，使實船的殼體沖擊因子與艙段試驗的3個工況對應相等，藥包的位置設置在艦船試驗艙段中部的正下方，相應的爆距設置如表3所示。

表3 艦船試驗工況Table 3 Ship test conditions

在將爆炸沖擊波壓力轉換為節點力輸入時，不同位置力的作用時間不同，相對于船體中部正下方爆炸工況，船體中部先受到力的作用，之后隨著沖擊波的傳播，兩側開始受力。上面已驗證了理論方法在實船計算中的有效性，以工況Ⅲ為例，應用式(8)可編程得到船體受到的沖擊作用力，如圖8，第10站為船中位置，其相對于一側的2站在受力時，約提前18 ms，而位于船長1/4位置的第5站，受力時間則在第10站和第2站之間，約比第10站延遲8 ms。從圖中還可看出離藥包越近的位置，受力越大，但作用時間也相對于爆距遠的位置作用時間短。在船體受力達到峰值時，爆距越遠，則力變化的越緩慢。

圖8 沿船長不同位置受力隨時間變化Fig.8 Variation curve of the different locations along the ship receiving the force with time

沿船長方向取不同節點，通過有限元計算得到的不同工況節點最大速度如圖9所示。工況Ⅰ的速度響應最大，工況Ⅲ最小。3種工況下船體兩端的速度響應較其他位置偏大，這與船體首尾重量較輕有關。同時船體中間部分沿船長有較小的起伏變化，但整體在某一穩定值附近變化。通過本文計算船體剛體運動速度的方法可得到艦船的剛體運動速度，用有限算法得到節點整體速度響應減去剛體運動速度則可認為是船體結構彈塑性變形的速度。選取船體中部節點，以工況Ⅰ為例，繪制出整體變形速度、剛體運動速度及彈塑性變形速度隨時間的變化如圖10所示，在沖擊波階段，3種速度響應都在8 ms左右出現最大值，之后由于剛體速度逐漸趨于零，整體速度和彈塑性變形速度變化曲線基本重合。

船體在3種工況下沿船長分布的各個節點最大速度響應的平均值如表4所示，從表中可以看出對于沖擊因子最大的工況Ⅰ整體速度響應最大，剛體運動速度也最大，同時彈塑性變形也是最大的，根據剛體速度占整體速度比重，得出全船實際得到的節點速度響應，包含總速度1/4左右的剛體運動速度，如不考慮剛體運動速度的響應，會對整體結構動響應問題分析帶來很大的誤差。

圖9 不同工況船體節點速度沿船長變化規律Fig.9 The change law of hull nodes'speed along the ship length in different cases

圖10 船體中部節點速度響應Fig.10 Velocity response of the hull central node

表4 船體最大速度響應Tab le 4 M axim um velocity of the hull

3.2 剛體運動動能對艦船響應影響

通過上文分析，剛體運動速度在船體整體速度響應分析中不可忽略，由剛體運動產生的剛體運動動能會削弱沖擊波入射能量作用給船體結構彈塑性變形的變形能，下面將討論剛體運動動能對船體結構吸收沖擊波入射能量的影響。

沖擊波對結構的入射能量Es可通過結構遮擋的沖擊波能量求得［12］，在無限域流場中可認為沖擊波能量均布于整個球形沖擊波波陣面，沖擊波的入射能量為

式中:E為沖擊波總能量，R為爆距，Se為結構在垂直于沖擊波波陣面上的投影面積，可通過文獻［12］中的方法進行求解，圖11為數值模擬沖擊波波陣面求解全船在波陣面上投影面積Se示意圖。

圖11 全船數值模擬Fig.11 The schematic diagram of ship numerical simulation

在一次爆炸中，沖擊波攜帶的能量為

式中:Ef為能流密度，指沖擊波在傳播方向上經過某一固定的單位面積內的能流，其表達式［1］為

由式(15)、(16)即可求得結構受到的沖擊波能量Es，通過式(12)可得到船體剛體運動最大速度，即可求得船體剛體運動消耗的能量Er，則船體彈塑性變形等吸收的變形能為(包括船體結構的彈塑性變形能和流體的耗散能，在此不研究流體耗散能量，將變形能定義為沖擊波入射能減去剛體運動動能):

在對船體艙段進行試驗考核全船的水下抗爆性能時，只是按照沖擊因子相等來設置試驗，而實際要達到的目的是要彈塑性變形的能量相同，艙段和實船實現相同的彈塑性變形響應，從而考核實船的抗爆性能。艙段試驗的工況1、2、3分別對應全船水下爆炸工況的Ⅰ、Ⅱ、Ⅲ，他們的沖擊因子分別相等。但根據艙段和藥包的相對位置及艙段試驗設計時大量的數值仿真計算結果，艙段試驗的3個工況是按照要考核的藥包在攻角30°時全船艙段對應位置響應進行設置的。工況Ⅳ、Ⅴ、Ⅵ分別對應實船工況Ⅰ、Ⅱ、Ⅲ，藥量和爆距均對應相同，但設置攻角為30°，如圖12 所示。

根據試驗艙段在全船水下爆炸工況中的投影面積，可求出艙段的沖擊波入射能;已知全船剛體運動速度，根據艙段占全船重量的比重，可得到艙段在全船工況中的剛體運動動能。以此得到的全船上局部艙段受到的結構變形能，大小應與艙段試驗的變形能相等。在按照艙段和全船沖擊因子相等設置工況時，艙段試驗和全船不同工況的沖擊波能量如表5～6所示。

圖12 全船工況設置Fig.12 Test conditions of ship

表5 艙段試驗能量計算Table 5 Energy calculation of cabin test

表6 艙段在全船工況下的能量計算Table 6 Energy calculation of cabin in ship conditions

從表5和表6中數據可以看出，無論是艙段試驗還是全船水下爆炸工況，船體剛體運動的動能相對于沖擊波入射能不能忽略，尤其艙段試驗的剛體運動動能占沖擊波入射能的比重更大，在工況2時高達24%。

最終計算出的艙段試驗各工況的變形能與全船各個工況的變形能并不相等，盡管艙段和全船各工況的沖擊因子對應相等。艙段試驗工況1的變形能為1 230 kJ，而全船上在對應沖擊因子相同的工況Ⅳ時對應艙段位置的變形能為1 720 kJ，比艙段變形能高出40%。工況Ⅴ實船上對應艙段位置變形能和艙段試驗工況2的變形能比較接近，認為此試驗艙段工況可以考核對應實船沖擊因子的工況，而實船工況Ⅵ的對應艙段位置變形能比艙段試驗工況3的變形能小，可見在攻角30°時，大沖擊因子的實船工況，艙段試驗需要通過適當加大沖擊因子來達到實船相等的變形能，在沖擊因子較小時艙段變形能可能大于實船的變形能。

通過上面分析，若要考核全船某一工況的爆炸沖擊動響應問題，只是根據沖擊因子相同考核局部艙段的動響應，并不能反映實船的真實響應。較為準確的方法是按照本文思路計算全船上對應艙段位置的變形能，根據此變形能和本文計算模型剛體運動速度的方法得出艙段剛體運動動能，得到試驗的沖擊波入射能，最終確定試驗的藥量和爆距。具體求解流程如圖13所示。

圖13 局部結構試驗工況設置求解框圖Fig.13 The solving block diagram of the test conditions set of local structure

通過此求解思路，計算出本文艙段試驗的工況只能考核如表7所示實船攻角為30°的工況。

表7 艙段試驗能考核的實船工況Table 7 Real ship conditions of the cabin test assessment

若要考核實船要求的工況Ⅳ、Ⅴ、Ⅵ沖擊因子，則根據同樣的算法，求出艙段試驗需要設置的工況沖擊因子如表8所示。

表8 要考核實船工況需設置的艙段試驗工況Table 8 Cabin test conditions to be set for assessment of ship conditions

4 結論

本文將復雜船體外殼等效為剛性平板，研究了船體在水下爆炸環境中球面沖擊波作用下的剛體動態響應特性，計算出船體受到的沖擊波作用力及剛體運動速度，研究了船體整體速度響應和剛體運動速度之間及沖擊波入射能、剛體運動動能之間的關系，得出如下結論:

1)船體在球面沖擊波作用下，計及沖擊波的反射及船體與流場流固耦合的附加質量，對于艙段模型試驗及實船試驗，當藥包在水下任意位置時，應用本方法可得出與試驗結果較符合的剛體運動速度;

2)在本文設定的實船工況中，船體整體速度響應、剛體運動速度及彈塑性變形速度的最大值基本出現在同一時刻，剛體速度響應占總速度響應的1/4左右，在水下爆炸沖擊波作用下進行結構動響應問題分析時，需考慮結構剛體運動對其的影響;

3)艙段模型試驗中剛體運動動能占沖擊波入射能的比重大于全船中艙段剛體運動動能占沖擊波入射能的比重，為使艙段模型試驗局部結構吸收的變形能與全船對應位置相等，可通過改變沖擊因子的方法來設置局部結構試驗的工況，考核實船對應沖擊因子的結構動態響應。

［1］COLE R H.Underwater explosion［M］.Princeton:Princeton University Press，1948:270-352.

［2］李良軍，蔣偉康，艾艷輝，等.水下爆炸載荷下小型目標變形及沖擊損傷試驗［J］.上海交通大學學報，2010，44(10):1450-1455.LILiangjun，JIANG Weikang，AI Yanhui，et al.Experimental study on deformation and shock damage of small target subjected to underwater explosions［J］.Journal of Shanghai Jiaotong University，2010，44(10):1450-1455.

［3］LIJ，RONG JL.Experimental and numerical investigation of the dynamic response of structures subjected to underwater explosion［J］.European Journal of Mechanics，2012，32:59-69.

［4］JIN Q K，DING G Y.A finite element analysis of ship sections subjected to underwater explosion［J］.International Journal of Impact Engineering，2011，38(7):558-566.

［5］JEN C Y，TAI Y S.Deformation behavior of a stiffened panel subjected to underwater shock loading using the non-linear finite elementmethod［J］.Materials ＆ Design，2010，31(1):325-335.

［6］李海濤，朱錫，張振華.水下爆炸球面沖擊波作用下船體梁的剛塑性動響應特性［J］.工程力學，2010，27(10):202-207.LIHaitao，ZHUXi，ZHANG Zhenhua.Dynamic rigid-plastic response of ship-like beam subjected to underwater spherical shockwaves［J］.Engineering Mechanics，2010，27(10):202-207.

［7］ZHANG N，ZONG Z.The effect of rigid-body motions on the whipping response of a ship hull subjected to an underwater bubble［J］.Journal of Fluids and Structures，2011，27:1326-1336.

［8］姚熊亮.艦船結構振動沖擊與噪聲［M］.北京:國防工業出版社，2007:144-146.

［9］汪玉，華宏星.艦船現代沖擊理論及應用［M］.北京:科學出版社，2005:220-223.

［10］張阿漫，戴紹仕.流固耦合動力學［M］.北京:國防工業出版社，2011:35-47.

［11］張文鵬，宗智，汪玉，等.沖擊波作用下艦船剛體運動響應［J］. 噪聲與振動控制，2012，32(6):134-136.ZHANGWenpeng，ZONG Zhi，WANG Yu，et al.Motion response of rigid ship under shock wave［J］.Noise and Vibration Control，2012，32(6):134-136.

［12］曹宇.水下爆炸沖擊環境相似性研究［D］.哈爾濱:哈爾濱工程大學，2006:27-30.CAO Yu.Research on environment comparability of underwater explosion［D］.Harbin:Harbin Engineering University，2006:27-30.