朱戰(zhàn)飛，石 全，王 琦，王廣彥，路 航

（1.解放軍63893部隊，河南 洛陽 471000；2.軍械工程學(xué)院，石家莊 050003）

小子樣條件下戰(zhàn)斗部爆炸仿真數(shù)據(jù)可信性分析

朱戰(zhàn)飛1，石 全2，王 琦2，王廣彥2，路 航2

（1.解放軍63893部隊，河南 洛陽 471000；2.軍械工程學(xué)院，石家莊 050003）

高效毀傷戰(zhàn)斗部試驗次數(shù)受到嚴(yán)格控制，如何在小子樣條件下分析仿真數(shù)據(jù)的可信性非常有必要。通過均勻試驗得到了戰(zhàn)斗部在地下目標(biāo)腔室內(nèi)部爆炸某點(diǎn)沖擊波最大超壓峰值，利用normplot及Lilliefors函數(shù)分析了仿真數(shù)據(jù)的分布特征，而后在已有檢驗方法的基礎(chǔ)上，結(jié)合順序統(tǒng)計量和均勻分布的特點(diǎn)提出了一種新的相容性驗證思路。研究表明：該方法不受樣本容量的限制，利用MATLAB程序編程可以快速判斷仿真數(shù)據(jù)的可信性。

小子樣，仿真數(shù)據(jù)，可信性，順序統(tǒng)計量，均勻分布

引言

隨著計算機(jī)技術(shù)和數(shù)值計算理論的發(fā)展，使得借助于LS-DYNA等軟件來分析戰(zhàn)斗部爆炸問題成為可能。通過數(shù)值運(yùn)算可以得到戰(zhàn)斗部爆炸過程的完整解，但仿真畢竟是對實際情況的模擬，其所得數(shù)據(jù)的可信性必須在通過真實的試驗數(shù)據(jù)驗證后，才能用于后續(xù)的毀傷分析及評估。然而在開展某些高效毀傷戰(zhàn)斗部爆炸的試驗中，受試驗周期和成本的限制，現(xiàn)場試驗次數(shù)往往受到嚴(yán)格的控制，呈現(xiàn)出明顯的小子樣特征，于是如何在小子樣條件下對仿真數(shù)據(jù)的可信性進(jìn)行驗證就顯得非常必要。

常用的可信性驗證思路是在比較相同輸入（試驗）條件下仿真與現(xiàn)場試驗結(jié)果的差異性，驗證仿真數(shù)據(jù)與實際數(shù)據(jù)的一致性程度，稱為一致性檢驗，用概率論與數(shù)理統(tǒng)計的思想可將其理解為兩組數(shù)據(jù)是否服從同一分布。通常采用的經(jīng)典顯著性檢驗方法，如柯爾莫哥洛夫-斯米爾諾夫方法、秩和檢驗法等一般要求有較為充足的試驗數(shù)據(jù)，不適合在小子樣條件下使用。本文將在小子樣條件下數(shù)據(jù)相容性研究成果基礎(chǔ)上［1-3］，進(jìn)一步改進(jìn)檢驗思路和方法，以鉆地戰(zhàn)斗部在地下目標(biāo)腔室內(nèi)部爆炸某點(diǎn)沖擊波最大超壓峰值為例，就仿真數(shù)據(jù)的可信性進(jìn)行詳細(xì)分析。

1 戰(zhàn)斗部在地下目標(biāo)腔室內(nèi)部爆炸的有限元模型構(gòu)建

地下目標(biāo)是指受到敵方威脅且需要嚴(yán)加保護(hù)的設(shè)施，從組成上來說包括防護(hù)裝置和被保護(hù)的人員及設(shè)備?？紤]到研究對象結(jié)構(gòu)復(fù)雜、工程龐大等特點(diǎn)，本文僅對戰(zhàn)斗部在某一單體目標(biāo)內(nèi)部爆炸的沖擊波傳播規(guī)律及內(nèi)部設(shè)備動態(tài)變形進(jìn)行數(shù)值計算。在利于研究又不失代表性的基礎(chǔ)上，假定單體目標(biāo)不設(shè)立通風(fēng)口和防護(hù)門，頂部防護(hù)層、四周墻壁、地面均為混凝土材料，厚度均為30 cm，內(nèi)部空間為長方形六面體結(jié)構(gòu)，長、寬、高分別為800 cm、500 cm、300 cm。腔室內(nèi)部設(shè)備簡化為材料為鋼、內(nèi)部為空心的無底長方體，固定在腔室中心位置，長、寬、高分別為400 cm、200 cm、150 cm，四周及頂部等效鋼板厚度均設(shè)定為1 cm。戰(zhàn)斗部為立方體裝藥，邊長為30 cm，上下表面中心同時起爆，裝藥質(zhì)心距腔室側(cè)面、底面的距離均設(shè)定為65 cm。根據(jù)以上分析建立的幾何模型如圖1所示。

圖1 幾何模型示意圖

TNT炸藥采用高能炸藥爆轟模型MAT_HIGH_EXPLOSIVE_BURN和EOS_JWL狀態(tài)方程來描述，各參數(shù)取值見文獻(xiàn)［4］，其中炸藥密度ρ=1.63 g/cm3；空氣采用MAT_NULL材料模型和線性多項式狀態(tài)方程EOS_LINEAR_POLYNOMIAL進(jìn)行描述，各參數(shù)取值見文獻(xiàn)［5］，其中空氣密度ρ=1.29×10-3g/cm3； 鋼 板 采 用 MAT_JOHNSON_COOK材料模型和EOS_GRUNEISEN狀態(tài)方程來描述，各參數(shù)取值見文獻(xiàn)［6］；混凝土采用MAT_JOHNSON_HOLMQUIST_CONCRETE材料模型來描述，即HJC模型，各參數(shù)取值見文獻(xiàn)［7］，其中混凝土密度ρ=2.44 g/cm3。

數(shù)值模型由炸藥、空氣、混凝土及設(shè)備4部分組成，其中炸藥、空氣及混凝土3種材料采用歐拉網(wǎng)格建模，單元使用多物質(zhì)ALE算法，使沖擊波可以通過空氣和四周混凝土介質(zhì)傳送；設(shè)備采用拉格朗日網(wǎng)格建模，建模時設(shè)備的網(wǎng)格與空氣網(wǎng)格重疊交叉，通過定義Lagrange和ALE的耦合實現(xiàn)爆炸過程中沖擊波對設(shè)備力學(xué)加載的流固耦合模擬。網(wǎng)格劃分時為保證數(shù)值模擬的精度，又不會占用大量的計算資源，試算后將炸藥、空氣及混凝土單元邊長設(shè)定為10，設(shè)備單元邊長設(shè)定為3，統(tǒng)一采用cm-g-μs-k單位制。所有模型均選用8節(jié)點(diǎn)三維實體單元SOLID164進(jìn)行映射網(wǎng)格劃分。為使沖擊波通過四周混凝土層外表面時透射出去而不產(chǎn)生反射效應(yīng)，在外表面添加無反射邊界條件。

2 采用均勻試驗設(shè)計獲取沖擊波峰值超壓仿真數(shù)據(jù)

在實際的毀傷試驗中，對于固定條件下的不同批次試驗，試驗結(jié)果是有差別的，這主要是由隨機(jī)噪聲等不可控因素造成的。下文將在合理分析影響沖擊波傳播的不可控因素基礎(chǔ)上，通過方案中因素水平的變動，采用均勻試驗設(shè)計的方法獲得一組符合實際的毀傷數(shù)據(jù)，并對所得仿真數(shù)據(jù)的特性進(jìn)行分析。

2.1 采用均勻試驗獲取沖擊波峰值超壓仿真數(shù)據(jù)

戰(zhàn)斗部在地下目標(biāo)腔室內(nèi)部爆炸具有明顯的“封閉效應(yīng)”，四周混凝土層密度的不同將直接影響沖擊波反射的強(qiáng)弱，而混凝土作為一種復(fù)合建筑材料，其密度并不是某一固定的數(shù)值，而是在某一范圍內(nèi)波動。由于系統(tǒng)誤差或人為因素的存在，實際試驗中定點(diǎn)爆炸的藥量總會存在小幅偏差，在仿真試驗中可以通過調(diào)整裝藥密度的小幅變動來控制。另外，腔室內(nèi)部空氣的溫度和壓強(qiáng)也是不可控因素，它們的改變將直接導(dǎo)致密度的變化，從而對沖擊波的傳播及衰減產(chǎn)生影響，在仿真試驗中可以直接通過調(diào)節(jié)空氣密度的大小來實現(xiàn)。

為更貼近實際，仿真試驗中將3種不可控因素均設(shè)定為上下波動，幅度為5%，每個因素劃分為15個水平。為保證合理利用各因素的水平值，取U15（158）均勻設(shè)計表來安排均勻試驗。按照均勻試驗方案，對每個試驗條件下的戰(zhàn)斗部爆炸進(jìn)行數(shù)值運(yùn)算，并提取出某點(diǎn)的最大超壓峰值，部分試驗結(jié)果如下頁表1所示，其中所得到的最大峰值超壓數(shù)據(jù)依次為：8.748 7，8.607 4，8.452 0，8.940 5，8.709 9，8.839 6，8.825 3，8.792 7，8.635 6，8.644 3，8.732 0，8.607 7，8.627 6，8.990 1，8.797 9。

2.2 仿真數(shù)據(jù)特性分析

所得峰值超壓數(shù)據(jù)呈現(xiàn)“中間多，兩頭少”這一特點(diǎn)，于是推測仿真數(shù)據(jù)服從正態(tài)分布。目前人們發(fā)現(xiàn)了一些專門的正態(tài)性檢驗方法，其檢驗效果比一般通用方法好，其中正態(tài)概率紙檢驗法是一種現(xiàn)場統(tǒng)計常用的判斷變量正態(tài)性的簡單工具，使用它可以很快地判斷變量是否服從正態(tài)分布。運(yùn)用MATLAB程序中的normplot函數(shù)繪制仿真數(shù)據(jù)的正態(tài)概率紙檢驗圖，如圖2所示。

表1 均勻試驗部分?jǐn)?shù)據(jù)

圖2 仿真數(shù)據(jù)正態(tài)概率紙檢驗圖

由圖2可以看出，仿真數(shù)據(jù)各點(diǎn)的位置近似的在一條直線附近，因而所得峰值超壓數(shù)據(jù)符合正態(tài)分布特性。但由于本文仿真數(shù)據(jù)并不是十分充足，為避免判斷失誤，進(jìn)一步采用適合小樣本條件的Lilliefors檢驗法進(jìn)行驗證［8］。

3 小子樣條件下仿真數(shù)據(jù)可信性分析

設(shè)x1t，x2t，…，xMt是通過仿真運(yùn)算獲得的數(shù)據(jù)樣本（一般M比較大），x1，x2，…，xn是現(xiàn)場試驗獲得的數(shù)據(jù)樣本（一般n比較?。?，現(xiàn)要求檢驗兩組數(shù)據(jù)是否服從同一分布。若是將仿真數(shù)據(jù)與現(xiàn)場試驗的小子樣數(shù)據(jù)簡單地混合在一起檢驗，往往檢驗結(jié)果是接受假設(shè)，體現(xiàn)不了差異，原因是大量的仿真信息淹沒了現(xiàn)場試驗信息，即使有差異也很難體現(xiàn)出來。

3.1 已有檢驗思路和方法

式中：F（x）為分布函數(shù)，f（x）為分布密度。在F（x）的分布參數(shù)已知的情況下，順序子樣相對應(yīng)的隨機(jī)變量序列記為αk：

在此隨機(jī)變量αk可被看作是從［0，1］區(qū)間內(nèi)均勻分布總體中抽取子樣的順序統(tǒng)計量，故x1，x2，…，xn與y1，y2，…，ym兩組子樣的一致性問題就轉(zhuǎn)化為檢驗隨機(jī)序列αk是否服從（0，1）均勻分布的問題。

3.2 改進(jìn)后的思路和方法

由3.1節(jié)可知文獻(xiàn)［1］的不合理處在于檢驗αk服從（0，1）均勻分布的思路及方法，本文將提出一種新的思路來檢驗αk是否服從（0，1）均勻分布。

3.3 仿真數(shù)據(jù)可信性分析

假設(shè)現(xiàn)場試驗所測得同一點(diǎn)處的最大超壓峰值數(shù)據(jù)為：8.751 5，8.832 4，8.552 9，單位為105Pa。根據(jù)3.1節(jié)與3.2節(jié)中所述的可信性分析思路：

①將仿真數(shù)據(jù)隨機(jī)分為5組，任意抽取其中的一組數(shù)據(jù)為：8.748 7，8.825 3，8.627 6，利用MATLAB程序中的sort函數(shù)將兩組子樣從小到大混合排序形成順序子樣x（k）。

②考慮到仿真信息充足以及現(xiàn)場試驗數(shù)據(jù)的小子樣特征，因而用2.2節(jié)中擬合得到的仿真子樣分布來代替總體分布F（x），即認(rèn)為總體服從N（8.731 8，0.143 92）正態(tài)分布，利用normcdf函數(shù)得到順序子樣x（k）對應(yīng)的F（x（k）），而后根據(jù)式（1）編程運(yùn)算得到對應(yīng)的αk，排序后形成順序子樣α（k），如表2所示。

④按照上述計算步驟，檢驗其余4組仿真數(shù)據(jù)的可信性。

表2 x（k）、αk、α（k）計算結(jié)果

表3 Zijk、z1、z2部分計算結(jié)果

4 結(jié)束語

本文以鉆地戰(zhàn)斗部在地下目標(biāo)腔室內(nèi)部爆炸為背景，以某點(diǎn)沖擊波最大超壓峰值為例，研究了小子樣條件下仿真數(shù)據(jù)的可信性。通過分析現(xiàn)有研究思路和方法的特點(diǎn)及不足，提出了一種新的判斷思路，該思路不受樣本容量的限制，利用MATLAB程序編程可以快速判斷仿真數(shù)據(jù)的可信性，主要包含以下幾個步驟：①將仿真數(shù)據(jù)隨機(jī)分組；②求解順序子樣α（k）；③構(gòu)造變尺度差商比序列統(tǒng)計量Zijk；④比較Zijk與z1、z2之間的大小關(guān)系，判斷該組仿真數(shù)據(jù)的可信性；⑤按照上述步驟檢驗其余組。

［1］唐雪梅.小樣本場合下相容性檢驗方法［J］.系統(tǒng)工程與電子技術(shù)，2001，23（10）:66-80.

［2］王 暉，潘高田，臧興震，等.正態(tài)分布的小樣本數(shù)據(jù)的相容性檢驗理論和方法［J］.數(shù)學(xué)的實踐與認(rèn)識，2005，35（3）:131-137.

［3］潘高田，潘 峰，王 暉，等.順序統(tǒng)計量在小樣本檢驗中的理論和方法研究［J］.工程數(shù)學(xué)學(xué)報，2005，22（2）: 301-306.

［4］鄧榮兵，金先龍，陳 峻，等.爆炸沖擊波對玻璃幕墻破壞作用的多物質(zhì)ALE有限元模擬［J］.高壓物理學(xué)報，2010，24（2）:81-87.

［5］陳 昊，陶 鋼，蒲 元.溫壓藥在有限空間內(nèi)爆炸沖擊波的實驗研究及數(shù)值模擬［J］.火炸藥學(xué)報，2009，32（5）:41-46.

［6］汪德武，高洪泉，杜海霞，等.斜侵徹彈體運(yùn)動分析與仿真［J］.彈箭與制導(dǎo)學(xué)報，2006，26（3）:121-123.

［7］Holmquist T J，Johnson G R，Cook W H.A Computational Constitutive Model for Concrete Subjected to Large Strains，High Strain Rates，and High Pressures［C］//Proceeding of the 14th International Symposium on Ballislics Quebec，1993：591-600.

［8］張德豐.MATLAB概率與數(shù)理統(tǒng)計分析［M］.北京：機(jī)械工業(yè)出版社，2010：188-191.

Compatibility for Simulation Data of Warhead Blast under Small Sample Circumstance

ZHU Zhan-fei1，SHI Quan2，WANG Qi2，WANG Guang-yan2，LU Hang2
（1.Unit 63893 of PLA，Luoyang 471000，China；
2.Ordnance Engineering College，Shijiazhuang 050003，China）

The times of experimentation for advanced destruction warhead are controlled strictly，it is very essential that how to analyze the compatibility of simulation data under small sample circumstance. The maximal pressure peak value of one point is gotten by uniform design when warhead blasted in the interior of underground target，the distribution characteristic of simulation data is studied by normplot and Lilliefors function.Then combining order statistic and uniform distribution one new provable method is advanced basing on existing methods.The research indicates that this new method couldn't be control of sample capacity and can judge the compatibility of simulation data quickly by MATLAB.

small sample，simulation data，compatibility，order statistic，uniform distribution

TB114.3

A

1002-0640（2014）10-0070-04

2013-08-05

2013-10-07

朱戰(zhàn)飛（1985- ），男，河南焦作人，碩士。研究方向：毀傷仿真、電子對抗訓(xùn)練。