艦船艙內爆炸載荷簡化載荷計算模型

陳鵬宇，侯海量，金 鍵，李 茂，朱 錫，焦立啟

(1. 海軍工程大學 艦船與海洋學院，湖北 武漢 430033；2. 海軍研究院，北京 100161)

0 引 言

隨著現代反艦武器的發展，各類大型艦船面對的威脅形式由船體外遠距離爆炸向精確制導式的穿入船體內爆，由于船體內部空間的限制，爆炸載荷在船體內部不斷反射匯聚，因此其對結構和人員的毀傷作用將比在船體外爆炸要更加顯著。為了更好地抵御戰斗部艙內爆炸載荷，研究艙內爆炸載荷的特性和規律就成為了防護領域的熱點問題。M.kurki[1]對艦船結構在艙內爆炸載荷作用的變形破壞模式進行了數值仿真研究。侯海量等[2-3]開展了典型艙室模型的艙內爆炸試驗，研究了艙內爆炸載荷的作用過程及艙室板架在內爆載荷下的失效模式。孔祥韶等[4]研究了3種不同的角隅連接結構型式對沖擊波在角隅匯聚情況的影響，并基于圖象法解釋了沖擊波在角隅的匯聚現象。C.Geretto等[5]進行了一系列實驗研究固方支板在不同空間約束程度（自由空間、半約束和完全約束）的內爆加載作用下的變形，分析了不同載荷下的變形規律和空間約束程度對變形的影響。姚術健等[6]進行了鋼箱結構在內部爆炸下的變形規律試驗研究，并擬合了箱體壁板撓厚比與爆炸當量之間的經驗公式。陳攀等[7]采用數值方法研究了艙室內爆沖擊波壁面反射特性及爆點位置對艙室內爆載荷的影響。楊亞東[8]等采用鏡像法和非線性疊加處理獲得了長方體密閉空間內爆炸沖擊波傳播和疊加分析模型。但對于內爆載荷作用的另一個重要參數，目標壁面受到的載荷總沖量尚沒有估算方法。本文采用時間及空間區域劃分的方法，建立艙內爆炸下目標壁面所受載荷的簡化模型，提出目標壁面受到的總沖量簡化工程計算方法，可為艦船艙室結構變形計算及抗內爆載荷設計提供參考。

1 艙內爆炸載荷簡化模型時間分布

對艙內爆炸載荷進行簡化處理，可以將其載荷先按作用時間順序分為爆炸沖擊波階段和準靜態壓力階段，第1階段爆炸沖擊波階段，沖擊波的多次反射作用由于都屬于作用時間短的脈沖載荷，基于沖量等效的原則，認為多次反射的波峰可以合并為一個大的尖三角脈沖載荷，且忽略其反射沖擊波強度上升的時間；第2階段準靜態壓力階段則主要受到后續準靜態壓力的作用，其特點為載荷峰值較小，但作用時間長，載荷強度的衰減速度相對較慢，故可以在簡化的瞬態尖脈沖載荷基礎上，將后續的準靜態壓力簡化為緊接著的一個衰減相對平緩、持續時間更長的長三角載荷。因此通過以上對試驗結果和有限元仿真計算結果的定性研究和分析，假設爆炸沖擊波的正壓作用時間結束時，艙內爆炸載荷再立即轉入準靜態壓力作用階段，在時間分布上將典型壁面上的內爆沖擊波載荷和準靜態壓力載荷簡化為瞬態尖脈沖載荷和長三角載荷的疊加作用，其載荷形式如圖1（b）所示，這種按時間順序的階段劃分方法與美國標準文件UFC-3-340-2[9]的中單點載荷替代整個壁面載荷的形式基本一致。其差別主要有兩點：一是沖擊波載荷和準靜態壓力是否有重疊作用的時間，UFC-3-340-2中沖擊波載荷與準靜態壓力載荷有一段同時作用的區域，而本文所用的載荷形式假設爆炸沖擊波的正壓作用時間結束時，再轉入準靜態壓力作用階段，即2種壓力載荷沒有重疊區；二是準靜態壓力載荷持續時間的不同。在UFC-3-340-2規范文件中將準靜態壓力從形成到衰減值大氣壓力的時間段內的載荷全部計入內爆載荷形式，而本文的準靜態壓力載荷只選取了準靜態壓力的前一段時間作為載荷的持續時間，因為有關研究結果表明[10]，在持續較長時間的爆炸沖擊波及后續準靜態壓力載荷作用下，金屬薄板架結構的動態響應存在一個飽和閾值。板架結構在沖擊波載荷作用下先產生大撓度變形，而后由于四周邊界條件的約束，板的中面會產生較大的拉伸膜力效應，進而使得板架結構承載能力隨其變形的增加而得到一定程度的增強，相比之下后續的載荷強度卻在逐漸衰減，這就導致結構的變形在到達一定值后不再繼續增加，其動態響應達到飽和閾值，后續的壓力載荷將不影響結構的最終變形。從結構響應動態響應飽和閾值的研究工作可見，艙內爆炸載荷的載荷形式需要考慮結構動態響應飽和閾值的影響，結構響應達到飽和閾值之后的持續作用的載荷為無效載荷，所以將結構變形達到飽和閾值的時間稱為飽和響應時間[11]，只將這段飽和響應時間內的持續載荷計入艙內爆炸作用的有效載荷且忽略該時段內準靜態壓力的衰減。因此本文最終給出的艙內爆炸載荷如圖1（b）所示。

圖1 艙內爆炸沖擊波載荷及準靜態壓力簡化模型（時間分布）Fig. 1 Simplified model of explosion shock wave load and quasi-static pressure in the cabin(time distribution)

圖1（b）中的艙內爆炸簡化載荷時間分布，為一個梯形脈沖載荷（沖擊波）與矩形脈沖載荷的聯合載荷，其載荷模型的關鍵參數有4個，分別是反射沖擊波超壓峰值p1、沖擊波超壓正壓作用時間t1、準靜態壓力值p2和準靜態壓力載荷有效作用時間t2。由于艙內爆炸的載荷特性不僅由時間分布決定，還受艙內室內爆炸載荷的空間位置的影響。因此為了得到這4個參數的簡化計算方法，需要對艙內爆炸載荷按空間進行區域劃分。

2 艙內爆炸載荷簡化模型空間分布

如圖2所示，以炸藥在壁面的垂直投影位置為原點建立壁面坐標系，壁面長邊長度為L，壁面短邊長度為H，艙室寬度為B，將艙內爆炸下壁面各點受到的爆炸沖擊波超壓峰值p1按空間分布劃分為三大類區域，第1類區域為非角隅中間區域，第2類為兩面邊界交匯的兩面角隅區，第3類為三面邊界交匯的三面角隅區。依據其各自的載荷特性，分別對三類區域的反射沖擊波超壓峰值p1建立其簡化計算方法。其中艙內爆炸載荷受角隅匯聚效應影響的區域范圍較為固定，為便于簡化計算，將邊界角隅匯聚區的統一寬度d取為壁面短邊長度H的1/10，即d=0.1H。

圖2 艙內爆炸下壁面峰值壓力空間分布模型Fig. 2 Spatial distribution model of peak pressure under the explosion in the cabin

2.1 反射沖擊波超壓峰值p1

2.1.1 非角隅中間區域

對于第1類的非角隅中間區域，將其分為8個相類似的三角區域。以圖2中區域3為例，對該區域的內爆載荷峰值p1進行簡化計算。將中間區域3內各點的爆炸沖擊波超壓峰值p1假設為一個連續的平面函數，該平面函數由三角區域的3個點p1（0，0），p1（0，x0-d）和p1（x0-d，y0-d）確定，即三點確定一個平面函數p1=f（x，y），再由所得平面函數計算出區域內各點（x，y）的超壓峰值p1（x，y），如圖3所示。為了求得該平面函數，就要先確定3個角點A點、B點和C點的超壓峰值計算方法。根據Henrych提出的經驗公式（1）可計算這3個點的爆炸沖擊波入射超壓峰值，其中比例距離x和y分為壁面坐標系下目標點的坐標，D為裝藥到壁面的垂直距離，W為裝藥質量。

圖3 中間區域3的內爆沖擊波反射超壓峰值分布Fig. 3 Peak distribution of shock wave reflection overpressure in the center area 3

為便于解釋說明，將這3個點的反射超壓峰值p1（0，0），p1（0，x0-d）和p1（x0-d，y0-d）分別記為pA，pB和pC。由于pA為炸藥正對壁面的投影位置處（0，0）的峰值超壓，因此（0，0）點處的爆炸沖擊波入射角φ為0°，沖擊波在壁面的反射為正反射，對于正反射，其峰值超壓的計算公式[12]為：

對于點（0，x0-d）和（x0-d，y0-d），要確定其反射沖擊波超壓峰值pB和pC則需要先確定其入射角φ的大小。以點（0，x0-d）為例，其入射角

當φ小于規則反射極限角φ0時，爆炸沖擊波在該點發生規則斜反射，其反射沖擊波超壓峰值仍然近似按式（2）計算，空氣的規則反射極限角近似取為40°；當φ大于等于規則反射極限角φ0時，爆炸沖擊波在該點處發生馬赫反射，其反射沖擊波峰值壓力pB按式（4）計算[12]：

其中ΔpmG為相同藥量在剛性壁面爆炸時空氣入射沖擊波超壓峰值。由于剛性壁面的限制和阻擋，爆炸沖擊波可看成是只向一半的無限空間傳播，進行鏡像對稱，則可看作是2倍的裝藥量在無限空氣域中爆炸，因此將WmG=2W代入式（1），可得：

將式（5）代入式（4）即可求得當該點發生馬赫反射時，其爆炸沖擊波反射超壓峰值p1。采用式（1）～式（5）的判定及計算方法，分別求出點（0，x0-d）和（x0-d，y0-d）處的反射沖擊波超壓峰值pB和pC。如此則分別求出了區域3的3個角點（0，0），（0，x0-d）和（x0-d，y0-d）的反射沖擊波超壓峰值pA，pB和pC。再由以上這3個點的坐標（0，0，pA），（0，x0-d，pB）和（x0-d，y0-d，pC）計算求得整個區域3內各個點的反射沖擊波超壓峰值所滿足的平面函數的表達式p1=f（x，y），其一般形式為：

其中：k1，k2，C1為根據（0，0，pA），（0，x0-d，pB）和（x0-d，y0-d，pC）三點坐標確定的常數，且對于區域3來說，C1=pA。

通過以上簡化計算方法，對于中間區域劃分成的8個小三角區域內任意一點，均可采用式（1）～式（6）的計算方法，求得其爆炸沖擊波反射超壓峰值p1。

2.1.2 兩面角隅區

對于圖2中第2類的兩面角隅區，根據其位置，也將其分成上下左右4個區域，以兩面角隅區1為例，依據仿真計算結果的角隅壓力載荷特性，建立該區域內各點的爆炸沖擊波反射超壓峰值p1的簡化分布及計算方法，如圖4所示。

圖4 兩面角隅區1的內爆沖擊波反射超壓峰值分布Fig. 4 Peak distribution of shock wave reflection overpressure in the dihedral corner region 1

對兩面角隅區的峰值及比沖量進行簡化計算時，可以認為對于兩面角隅區1，在其區域統一寬度（d）的方向上，反射沖擊波的峰值壓力和比沖量無變化，即p1（x1，y0）=p1（x1，y1）=p1（x1，y2），其中（x1，y0），（x1，y1），（x1，y2）為兩面角隅區1內x坐標相同，y坐標不同的3點。

因此在圖4中的兩面角隅區1中任意一點（x，y）的反射沖擊波超壓峰值只是與x坐標有關的函數，即p1（x，y）=f（x，y0）=f（x）；將兩面角隅區 1 的爆炸沖擊波反射超壓峰值假設為線性分段函數p1=f（x），該分段函數由圖4中D點、E點和F點這3點處的反射沖擊波超壓峰值確定，而這3個點在壁面坐標系下的坐標分別為（x0-L+d，y0），（0，y0）和（x0-d，y0）。在求出了線性分段函數后，即可采用該函數計算求得兩面角隅區1內任意點處的爆炸沖擊波反射超壓峰值。因此，需要先求解這3個點處的反射沖擊波超壓峰值。為便于解釋說明，將這3個點的反射超壓峰值p1（x0-L+d，y0）、p1（0，y0）和p1（x0-d，y0）分別記為pD、pE和pF。

以點D處反射峰值壓力pD求解為例，先根據式（1）計算D點處的爆炸沖擊波入射超壓峰值Δpm，再代入式（2）～式（5）求出點D處不考慮角隅匯聚增強效應的反射沖擊波峰值壓力pDO。點D屬于兩面角隅匯聚邊界，相鄰的另一壁面角隅點的反射沖擊波會在D點與其反射沖擊波發生疊加效應，近似相當于D點有2倍反射沖擊波超壓的作用，因此可以認為點D處反射峰值壓力pD的最終計算式為：

同樣計算步驟及方法對于點E和F有：

由獲得的D，E，F三個點的坐標（x0-L+d，y0，pD），（0，y0，pE）和（x0-d，y0，pF）確定線性分段函數p1=f（x），如下式：

其中：k3為圖4中DE直線段的斜率；k4為圖4中EF直線段的斜率；C2為（0，y0）處反射沖擊波超壓峰值。

通過式（10）可以計算出兩面角隅區1內任意一點（x，y）的反射沖擊波超壓峰值壓力p1。對于其他的兩面角隅區的任意一點，均可采用式（1）～式（5）、式（7）～式（10）的計算方法，求得其爆炸沖擊波反射超壓峰值p1。

2.1.3 三面角隅區

對三面角隅區反射超壓峰值進行簡化計算時，忽略這部分反射沖擊波壓力和比沖量差值，近似地認為區域內各點的反射沖擊波壓力和比沖量相同，然后根據沖量等效原則將三面角隅區的載荷等效成均布載荷，各個三面角隅區內均布壓力載荷的初始峰值分別由角隅位置4點（G點，M點，N點，I點）的反射超壓峰值確定，即三面角隅區內各點處的反射沖擊波峰值均與其相應的角隅點的峰值壓力相同。

以圖5所示三面角隅區2為研究對象，建立其反射沖擊波超壓峰值簡化計算方法。首先根據式（1）可計算G點處的爆炸沖擊波入射超壓峰值Δpm，再代入式（2）～式（5）求出點G處不考慮角隅匯聚增強效

圖5 三面角隅區位置分布Fig. 5 Location distribution of the trihedral corner area

應的反射沖擊波峰值壓力pGO。點G屬于三面角隅匯聚邊界，相鄰的另2個壁面角隅點的反射沖擊波會在G點與其反射沖擊波發生疊加效應，近似相當于G點有3倍的反射沖擊波超壓的作用，因此可以將點G處反射峰值壓力pG的最終計算式為：

再根據之前的均布載荷近似簡化處理，則對三面角隅區2內任意一點（x，y）的反射沖擊波超壓峰值壓力p1（x，y）均與G點處相同，即p1（x，y）為常函數，所以有：

其中對于三面角隅區2，C3=pG。對于其他的兩面角隅區的任意一點，均可采用式（1）～式（5）、式（11）～式（12）的計算方法，求得其爆炸沖擊波反射超壓峰值p1。

綜上所述，艙內爆炸載荷下目標壁面上任意一點的反射沖擊波超壓峰值p1的求解過程如圖6所示。

2.2 沖擊波超壓正壓作用時間t1

由于反射沖擊波超壓峰值p1在計算時被劃分成了3個不同的區域，沖擊波超壓正壓作用時間也需要按各自的區域分別計算。

圖6 反射沖擊波超壓峰值p1求解流程圖Fig. 6 Flow chart of reflection shock wave overpressure peak p1

2.2.1 非角隅中間區域

由于簡化模型的時間分布忽略了沖擊波的到達時間而假設各點的沖擊波同時到達，因此對于同一個區域內的各點，也假設其沖擊波超壓的持續時間相同。同樣以中間區域3為研究對象，計算其沖擊波超壓的正壓作用時間t1。

空氣中爆炸沖擊波正壓作用時間的計算公式為：

根據式（13）計算圖3中A點、B點和C點的沖擊波正壓作用時間tA+，tB+和tC+，由于之前假設同于區域內的各個點的正壓作用時間相同，因此取這3點正壓作用時間的平均值為區域內的正壓作用持續時間，即有：

由式（14）計算可得中間區域3內各點沖擊波正壓作用時間。對于中間區域的8個小區域都可采用相同的方法和式（14）計算其區域內的正壓作用時間。

2.2.2 兩面角隅區

對于兩面角隅區的沖擊波超壓正壓作用時間t1，采用和中間區域相同的簡化處理方法求得兩面角隅區域內平均化的沖擊波正壓作用時間。以兩面角隅區1為例，計算其的區域內的沖擊波超壓正壓作用時間。首先采用式（13）分別計算圖4兩面角隅區1中3點D點、E點和F點的爆炸沖擊波正壓作用時間tD+，tE+和tF+，同樣取這3點的正壓作用時間的平均值為兩面角隅區1內的平均正壓作用持續時間，即有：

由式（13）和式（15）計算可得兩面角隅區內各點沖擊波正壓作用時間。對于其余的兩面角隅區都可采用相同的方法計算其區域內的正壓作用時間t1。

2.2.3 三面角隅區

由于三面角隅區的范圍較小，在該范圍內各點的正壓作用時間隨位置的變化可近似忽略，因此對于三面角隅區內各點的沖擊波超壓正壓作用時間t1，可用角隅一點的正壓作用時間代替。以三面角隅區2為例，計算其區域內的沖擊波超壓正壓作用時間。采用式（13）計算圖5（三面角隅區2）中G點的正壓作用時間tG+，由于以G點的正壓作用時間代替整個三面角隅區2內各點的沖擊波正壓作用時間，因此該區域內各點的正壓作用時間為：

其他三面角隅區的沖擊波超壓正壓作用時間t1同樣可由式（13）和式（16）的方法計算得到。

2.3 準靜態壓力值和準靜態壓力載荷有效作用時間

對于艙內爆炸下的準靜態壓力值p2和準靜態壓力載荷有效作用時間t2，不需要采用區域劃分的方法進行簡化，而是假設壁面上任何一點的的準靜態壓力值和有效作用時間都相同。

2.3.1 準靜態壓力值

對于艙內爆炸下艙內爆炸下的準靜態壓力值p2的計算，其核心問題之一是有限密閉空間內爆炸產物形成的準靜態壓力值的求解，而研究表明有限空間內爆炸形成的準靜態壓力值是裝藥質量與有限空間容積之比（W/V）的函數，并有一系列的基于該比值的擬合經驗公式，本文則采用Carlson經驗公式計算艙內爆炸形成的準靜態壓力值：

其中：W為裝藥質量，V為艙室內部的空間體積。對于壁面上任意一點（x，y），其準靜態壓力值p2可用式（17）計算。

2.3.2 準靜態壓力載荷有效作用時間

根據第1節的時間分布載荷的簡化過程分析，準靜態壓力載荷的另一個重要參數是有效作用時間，因為結構在持續較長時間的爆炸沖擊波及后續準靜態壓力的作用下其變形存在一個飽和閾值，在變形達到飽和閾值之后，后續載荷的作用將不再引起更多的變形。因此準靜態壓力載荷的有效作用時間t2即是結構變形到達飽和閾值的響應飽和時間ts。孔祥韶等[8]針對密閉空間內爆炸載荷計算方法及結構的飽和響應時間進一步開展了系列仿真計算研究，并提出了內爆載荷下板結構的無量綱的飽和響應時間的經驗公式：

其中：ts為板結構的飽和響應時間，ms；λ為無量綱系數，取值范圍為16.0～17.5；L0為正方形板邊長，m；ρ為板材料密度，g/cm3；σ為屈服強度，MPa。從式（18）可見，板結構的飽和響應時間與結構板厚無關。但該經驗公式是基于正方形的平板結構建立的，而實際情況下艙室的壁面往往是長度和高度不相同的矩形，因此在式（18）基礎上稍作修改。取矩形壁面長邊和短邊的平均值作為長度參數，以盡量同時考慮壁面長邊和短邊對其飽和響應時間的因素，因此有：

由之前的分析和圖1（b）的時間分布模型可知，結構變形到達飽和閾值的響應飽和時間ts即是準靜態壓力載荷有效作用的時間t2。因此由式（19）即可求得內爆載荷下目標壁面上任意一點（x，y）處準靜態壓力載荷的有效作用時間t2。

2.3.3 目標壁面受到的艙內爆炸載荷的總沖量

前文建立了艙內爆炸載荷簡化時空分布模型，對于內爆載荷下目標壁面上任意一點的艙內爆炸載荷，基于經驗公式基礎上，給出了相應的求解其模型關鍵參數的簡化計算方法。通過圖1（b）的載荷時間分布形式可見，由這4個關鍵載荷參數（p1，t1，p2，t2）可計算出任意一點（x，y）處受到的內爆載荷有效比沖量I（x，y），其計算式為：

其中p1，t1，p2，t2都是與坐標位置（x，y）有關的函數。由于p1和t1的簡化計算方法都與目標壁面的區域劃分有關，因此在計算壁面受到的總沖量時，也需要先分別按各自區域求解其沖量，在進行求和得到壁面的總沖量。

圖7 艙內爆炸作用下目標壁面比沖量計算示意圖Fig. 7 Schematic diagram of calculation of target wall specific impulse under the action of explosion in cabin

對于第一類中間非角隅區內的8個區域，記中間區域i的受到的沖量為：

圖8 目標壁面總沖量計算區域劃分Fig. 8 Division of total impulse calculation area on target wall

以圖8中間區域3為例，艙內爆炸簡化載荷對整個中間區域3的面積進行二重積分，得到該區域的總沖量為：

同理可求得第1類的中間區域內4個分區域i各自的沖量I1i，則第1類中間區域受到的總沖量為：

對于第2類兩面角隅區的4個區域，記兩面角隅區域j受到的沖量為I2j，以圖8兩面角隅區1為例，艙內爆炸簡化載荷對整個兩面角隅區1的面積進行二重積分，得到該區域的總沖量為：

同理可求得第2類的兩面角隅區內4個分區域j各自的沖量I2j，則第2類兩面角隅區受到的總沖量為：

對于第3類三面角隅區的4個區域，記三面角隅區域k受到的沖量為I3k,以圖8中三面角隅區2為例，艙內爆炸簡化載荷對整個三面角隅區2的面積進行二重積分，得到該區域的總沖量為：

同理可求得第3類的三面角隅區內4個分區域k各自的沖量I3k，則第3類兩面角隅區受到的總沖量為：

采用簡化計算模型得到目標壁面受到的艙內爆炸載荷總沖量為：

3 模型結果與仿真結果對比

將所述艙內爆炸載荷簡化模型通過Matlab程序化實現，并采用文獻[13]中相同的工況和艙內爆炸數值仿真方法，建立有限元模型并計算分析，驗證所提出簡化計算模型的有效性，其計算結果對比如表1所示。從計算結果可見，對于文獻中所提供2個艙內爆炸工況的仿真計算結果，本文提出的簡化計算模型計算結果與仿真結果的誤差都在10%以內，吻合情況較好。

表1 仿真及簡化模型計算結果Tab. 1 Results of simulation and simplification

4 結 語

本文對艙內爆炸載荷相關參數計算進行了一定的簡化處理，將艙內爆炸下目標壁面受到的內爆載荷在時間分布上分為爆炸沖擊波作用和準靜態壓力作用2個階段，在空間分布上分為中間非角隅區兩面角隅區和三面角隅區，并基于此分布建立了艙內爆炸下目標壁面所受載荷的簡化計算模型，提出了目標壁面所受總沖量的簡化計算方法，其計算結果與有限元仿真結果吻合較好，對艙室結構的抗內爆載荷設計具有一定的實用價值，可以為艙內爆炸載荷強度及結構變形的快速工程估算提供參考。