張 婧，袁 海，王春雨，閆 巖

(江蘇科技大學(xué) 船舶與海洋工程學(xué)院，江蘇 鎮(zhèn)江 212003)

水下爆炸載荷的統(tǒng)計(jì)特性

張 婧，袁 海，王春雨，閆 巖

(江蘇科技大學(xué) 船舶與海洋工程學(xué)院，江蘇 鎮(zhèn)江 212003)

為研究水下爆炸載荷的概率密度，選取炸藥密度、水密度作為隨機(jī)變量，基于乘同余組合發(fā)生器，利用Fortran語言編寫隨機(jī)數(shù)生成程序，產(chǎn)生100組隨機(jī)變量的樣本；采用有限元程序LS-DYNA對水下爆炸沖擊波進(jìn)行仿真計(jì)算，得到?jīng)_擊波的峰值壓力；并檢驗(yàn)其是否服從正態(tài)分布，驗(yàn)證最大熵法程序的正確性，采用最大熵法擬合其概率密度函數(shù)，研究爆炸載荷的統(tǒng)計(jì)規(guī)律，得到不同爆距，不同時(shí)刻的沖擊波峰值壓力的概率密度。該方法可較好地解決爆炸載荷的隨機(jī)數(shù)理統(tǒng)計(jì)模型，為進(jìn)一步進(jìn)行結(jié)構(gòu)的可靠性分析提供理論依據(jù)。

水下爆炸載荷；峰值壓力；最大熵法；概率密度

0 引 言

目前國內(nèi)外對水下爆炸沖擊波載荷的研究主要有實(shí)驗(yàn)法、數(shù)值仿真法。由于實(shí)驗(yàn)研究中測量數(shù)據(jù)困難、危險(xiǎn)大、實(shí)驗(yàn)費(fèi)用高昂等原因，所以大部分研究采用的是數(shù)值仿真方法。水下爆炸過程非常不穩(wěn)定，很難重復(fù)進(jìn)行，很大程度上體現(xiàn)了爆炸載荷的隨機(jī)性。雖然結(jié)構(gòu)在爆炸載荷作用下的響應(yīng)分析有許多學(xué)者做過相關(guān)的研究[1 – 4]，并取得豐富成果，但對爆炸載荷的統(tǒng)計(jì)模型只有美國的USACE進(jìn)行了43組爆炸試驗(yàn)，選擇不同的爆距、裝藥量得到爆炸波壓力時(shí)程曲線。而從國家安全性考慮，詳細(xì)的載荷數(shù)據(jù)，尤其是裝藥量和爆距等參數(shù)都不對外公開[5]。

考慮到缺少關(guān)于爆炸載荷統(tǒng)計(jì)規(guī)律的信息，假設(shè)參數(shù)服從正態(tài)分布來實(shí)現(xiàn)對爆炸載荷的統(tǒng)計(jì)研究，目前水下爆炸載荷的計(jì)算公式僅限于Cole提出的經(jīng)驗(yàn)公式[6]。本文采用蒙特卡羅方法，隨機(jī)生成參數(shù)數(shù)組，通過有限元分析，得到?jīng)_擊波的峰值壓力，并檢驗(yàn)其是否服從正態(tài)分布，利用最大熵法擬合其概率密度函數(shù)，得到不同爆距、不同時(shí)刻的沖擊波峰值壓力的概率密度函數(shù)，分析爆炸載荷參數(shù)對統(tǒng)計(jì)規(guī)律的影響。

1 最大熵法

對于一個(gè)隨機(jī)變量，熵定義為[7]

而對于離散隨機(jī)變量，熵為

式中：f（X）為隨機(jī)變量X的概率密度函數(shù)；f（xi）為離散點(diǎn)的概率密度函數(shù)。

最大熵法的基本方程為

式中，mi為X的第i階原點(diǎn)矩。

由駐值條件

可以得到

由式（8）可知，要使等式成立，被積函數(shù)必定為0，所以有

由式（9）可以得到最大熵概率密度函數(shù)為

下面確定拉格朗日乘子λ0和λi的值，將式（9）代入式（4）中，得到

將式（12）對λi求導(dǎo)，可得

將式（13）對λi求導(dǎo)，可得

由式（14）和式（15）可得

根據(jù)式（16）可求得λ1，λ2，…，λi，再由式（14）可得到λ0，為了方便求解，將式（16）作如下變形：

然后利用最小二乘法原理，使余量的平方和為最小的條件來確定λi，即

求解上式較好的方法是采用非線性優(yōu)化方法。

采用Fortran語言編制最大熵法擬合概率密度函數(shù)的程序，分別用標(biāo)準(zhǔn)正態(tài)分布、對數(shù)正態(tài)分布、Weibull分布，極值Ⅰ型分布[8]來驗(yàn)證自編程序的正確性。

標(biāo)準(zhǔn)正態(tài)分布密度函數(shù)為

對數(shù)正態(tài)分布密度函數(shù)為

Weibull分布密度函數(shù)為

極值Ⅰ型密度函數(shù)為

最大熵法與解析解的比較如圖1所示，A實(shí)線曲線為采用最大熵法擬合得到的密度函數(shù)，B虛線曲線為概率密度的真實(shí)分布。從圖1可看出，采用最大熵法擬合得到的密度函數(shù)與真實(shí)的密度函數(shù)分布比較，兩者吻合較好，說明自編程序最大熵法能夠很好地?cái)M合密度函數(shù)。

2 水下爆炸載荷的統(tǒng)計(jì)特性

2.1 有限元模型的建立

水域的網(wǎng)格尺寸大小為50 mm，建模時(shí)避免三角形單元，采用Solid單元，水域添加無反射邊界條件。有限元模型如圖2所示。

2.2 隨機(jī)變量的確定

2.2.1 確定隨機(jī)變量的均值和變異系數(shù)

選用TNT炸藥，假設(shè)其密度ρ服從正態(tài)分布，水的密度服從正態(tài)分布。各隨機(jī)變量的均值和變異系數(shù)如表1所示。

表 1 各隨機(jī)變量的均值μ和變異系數(shù)VTab. 1 Mean and coefficient of variation of the random variables

2.2.2 隨機(jī)數(shù)的生成

利用乘同余法產(chǎn)生均勻分布的隨機(jī)數(shù)[9]：

式中：α和M均為預(yù)先選定的常數(shù)；ri為[0，1]之間均勻分布的隨機(jī)數(shù)。

采用隨機(jī)數(shù)產(chǎn)生程序得到100個(gè)均勻分布隨機(jī)數(shù)。

2.2.3 數(shù)值仿真計(jì)算結(jié)果

沖擊波在不同時(shí)刻的傳播過程如圖3所示，沖擊波峰值壓力呈指數(shù)衰減趨勢。

水中沖擊波峰值壓力的平均值與經(jīng)驗(yàn)公式對比，如圖4 所示。通過對峰值壓力的統(tǒng)計(jì)分析，可看出：

1）隨著爆距的增加，水下爆炸沖擊波的峰值壓力衰減較快，隨著沖擊波的傳播，當(dāng)爆距大于0.5 m時(shí)，沖擊波的峰值壓力的衰減速率相對減??；

2）水下爆炸沖擊波的峰值壓力與裝藥量呈線性相關(guān)，隨著裝藥量的增加，峰值壓力相應(yīng)增大，同時(shí)傳播介質(zhì)（水）也影響著沖擊波峰值壓力大小，裝藥量相同的情況下隨著水的密度增加沖擊波峰值壓力也相應(yīng)增大，這表明了傳播介質(zhì)也影響了沖擊波的壓力峰值。

2.2.4 結(jié)果統(tǒng)計(jì)分析

2.2.4.1W檢驗(yàn)

W檢驗(yàn)，又稱Shapiro-Wilk檢驗(yàn)[9]，是一種基于相關(guān)性的算法。通過計(jì)算可得到一個(gè)相關(guān)系數(shù)，如果數(shù)值越接近1，就表明數(shù)據(jù)和正態(tài)分布擬合得越好。

計(jì)算式為：

其檢驗(yàn)步驟如下：

②計(jì)算上式分母；

③計(jì)算α值；

④計(jì)算檢驗(yàn)統(tǒng)計(jì)量W；

⑤若W值小于判斷界限值Wα，按顯著性水平α舍棄正態(tài)性假設(shè)；若W＞W(wǎng)α，接受正態(tài)性假設(shè)。

利用origin8.0程序直接進(jìn)行W檢驗(yàn)計(jì)算，若W＞W(wǎng)α，則服從正態(tài)分布。

2.2.4.2 峰值壓力概率密度統(tǒng)計(jì)

對100組數(shù)據(jù)統(tǒng)計(jì)分析，選取爆距為0.1 m處的壓力峰值進(jìn)行頻數(shù)統(tǒng)計(jì)，如圖5所示。由于頻數(shù)直方圖的形狀類似正態(tài)分布，假設(shè)0.1 m處沖擊波峰值壓力服從正態(tài)分布，則均值和標(biāo)準(zhǔn)差的估計(jì)值為：

擬合概率密度函數(shù)為：

利用origin8.0程序直接進(jìn)行W檢驗(yàn)計(jì)算可知，W=0.982 0。對于n=50，顯著性水平α=0.05時(shí)，W的Wα臨界值=0.947。W＞W(wǎng)α，服從正態(tài)分布。概率密度曲線如圖6所示。

2.2.4.3 不同時(shí)刻概率密度統(tǒng)計(jì)

對0.002 48 ms沖擊波壓力峰值進(jìn)行頻數(shù)統(tǒng)計(jì)，如圖7所示。由于頻數(shù)直方圖的形狀類似正態(tài)分布，假設(shè)0.024 8 ms沖擊波峰值壓力服從正態(tài)分布，則均值和標(biāo)準(zhǔn)差的估計(jì)值為：

擬合概率密度函數(shù)為：

利用origin8.0程序直接進(jìn)行W檢驗(yàn)計(jì)算可知，W=0.972 5。對于n=50，顯著性水平α=0.05時(shí)，W的Wα臨界值=0.947。W＞W(wǎng)α，服從正態(tài)分布。沖擊波峰值的概率密度曲線如圖8所示。

對0.007 47 ms沖擊波壓力峰值進(jìn)行頻數(shù)統(tǒng)計(jì)，如圖9所示。利用origin8.0程序直接進(jìn)行W檢驗(yàn)計(jì)算可知，W=0.972 5。對于n=50，顯著性水平α=0.05時(shí)，W的Wα臨界值=0.947。W＜Wα，不服從正態(tài)分布。利用最大熵法求解沖擊波峰值的概率密度曲線如圖10所示。

選取0.125 ms的壓力峰值進(jìn)行頻數(shù)統(tǒng)計(jì)，如圖11所示。利用origin8.0程序直接進(jìn)行W檢驗(yàn)計(jì)算可知，W=0.912 0。對于n=50，顯著性水平α=0.05時(shí)，W的Wα臨界值=0.947。W＜Wα，不服從正態(tài)分布。利用最大熵法求解沖擊波峰值的概率密度曲線如圖12所示。

3 結(jié) 語

采用蒙特卡羅和有限元相結(jié)合的方法，對水下爆炸載荷的統(tǒng)計(jì)規(guī)律進(jìn)行研究，得出如下結(jié)論：

1）水（傳播介質(zhì)）的密度影響了沖擊波壓力峰值，水（傳播介質(zhì)）的密度越大沖擊波峰值壓力也相應(yīng)增大。

2）建立基于較少數(shù)量的樣本數(shù)據(jù)得到水下爆炸載荷的隨機(jī)信息的樣本擬合法。給出了隨機(jī)變量正態(tài)分布類型的檢驗(yàn)方法以及利用最大熵法擬合非正態(tài)分布隨機(jī)變量概率密度函數(shù)的方法。最大熵法可較好地?cái)M合爆炸載荷的隨機(jī)數(shù)理統(tǒng)計(jì)模型。

3）本文方法能夠?qū)^少數(shù)量的樣本進(jìn)行擬合的基礎(chǔ)上，充分利用隨機(jī)變量的具體信息，采用最大熵法擬合水下爆炸載荷的概率密度，方法簡便，避免了成本較高的水下爆炸試驗(yàn)。

[1]WANG Z Y, LIANG X, FALLAH A S. A novel efficient method to evaluate the dynamic response of laminated plates subjected to underwater shock[J]. Journal of Sound and Vibration, 2013, 332(21): 5618–5634.

[2]姚熊亮, 郭君, 曹宇. 在水下爆炸沖擊波作用下的新型沖擊因子[J]. 中國造船, 2008, 49(2): 52–60.YAO Xiong-liang, GUO Jun, CAO Yu. A new impulsive factors on the underwater shock load[J]. Shipbuilding of China, 2008,49(2): 52–60.

[3]ZONG Z, ZHAO Y J, LI H T. A numerical study of whole ship structural damage resulting from close-in underwater explosion shock[J]. Marine Structures, 2013, 31: 24–43.

[4]CHEN Y, CHEN F, DU Z P, et al. Protective effect of polymer coating on the circular steel plate response to near-field underwater explosions[J]. Marine Structures, 2015, 40:247–266.

[5]EAMON C, BAYLOT J T, O'DANIEL J L. Modeling of concrete masonry walls subjected to explosive loads[J]. Journal of Engineering Mechanics, 2004, 130(9): 1098–1106.

[6]庫爾. 水下爆炸[M]. 北京: 國防工業(yè)出版社. 1960.

[7]陳虬, 劉先斌. 隨機(jī)有限元法及其工程應(yīng)用[M]. 成都: 西南交通大學(xué)出版社, 1993.

[8]何水清, 王善. 結(jié)構(gòu)可靠性分析與設(shè)計(jì)(第一版)[M]. 北京: 國防工業(yè)出版社, 1993.

[9]肖剛, 李天柁. 系統(tǒng)可靠性分析中的蒙特卡羅方法(第二版)[M]. 北京: 科學(xué)出版社, 2003.

Statistic characteristic of underwater explosion load

ZHANG Jing, YUAN Hai, WANG Chun-yu, YAN Yan
(School of Naval Architecture and Ocean Engineering, Jiangsu University of Science and Technology, Zhenjiang 212003, China)

In order to study probability density of underwater explosion load, density of explosive and density of water are considered as basic random variables. 100 samples of random variables are obtained using the Fortran program by which random numbers could be generated. The simulation of underwater explosion load is carried out by LS-DYNA, peak pressure of shock wave is obtained and verify whether it obey the normal distribution. Verify the correctness of the maximum entropy method, maximum entropy method is used to fit the probability density function of the peak pressure in different blast distance and different time. The method can be used to solve the random statistics model of explosion load, and it can provide the theoretical basis for the reliability analysis of the structure.

underwater explosion load；peak pressure；maximum entropy method；probability density

O38

A

1672 – 7649(2017)09 – 0012 – 05

10.3404/j.issn.1672 – 7649.2017.09.003

2016 – 09 – 05；

2016 – 11 – 01

國家自然科學(xué)基金資助項(xiàng)目(51509113)；江蘇省高校自然基金資助項(xiàng)目(16KJA580003，14KJB580005)；江蘇省重點(diǎn)實(shí)驗(yàn)室開放基金資助項(xiàng)目(CJ1403)

張婧(1983 – )，女，博士，副教授，研究方向?yàn)楸W(xué)。