宮 雷， 許世蒙

(陸軍裝甲兵學(xué)院基礎(chǔ)部， 北京 100072)

彈道侵徹的數(shù)值模擬綜合利用流體力學(xué)、溫度場、電傳導(dǎo)、磁場等領(lǐng)域的相關(guān)理論,將計(jì)算機(jī)軟件、數(shù)學(xué)算法與彈道侵徹機(jī)理相結(jié)合，成功地解決了該領(lǐng)域的實(shí)際問題,具有廣泛的應(yīng)用價(jià)值，是3維編織復(fù)合材料的彈道侵徹研究的重要方法之一[1-2]。

數(shù)值模擬過程中，模型診斷是必不可少的環(huán)節(jié)[3-5]，傳統(tǒng)的模型診斷方法[6-7]基于殘差分析，通過擬合優(yōu)度檢驗(yàn)和F檢驗(yàn)，對數(shù)值模擬與試驗(yàn)結(jié)果的剩余速度在同一模型下殘差的分布進(jìn)行同一性檢驗(yàn)。其中：擬合優(yōu)度檢驗(yàn)通過決定系數(shù)判斷回歸方程對樣本數(shù)據(jù)的擬合程度；F檢驗(yàn)對模型中回歸變量與響應(yīng)變量之間的線性關(guān)系是否顯著成立進(jìn)行判斷。這些方法在一定程度上是有效的，但其條件要求較高，計(jì)算也較復(fù)雜。筆者以試驗(yàn)結(jié)果為基礎(chǔ)，在模型診斷中采用分段趨勢分析，該方法操作簡單、行之有效，較好地彌補(bǔ)了傳統(tǒng)模型診斷方法的上述不足。

1 試驗(yàn)與數(shù)值

1.1 試驗(yàn)材料及方法

靶板制備[8]時(shí)，纖維材料采用荷蘭Tarpon 1000型對位芳族聚酞胺纖維長絲束，基體使用熱固性環(huán)氧樹脂，采用四步法1×1縱橫陣列編織機(jī)進(jìn)行織造。

采用56式彈道發(fā)射槍，質(zhì)量為7.95 g的7.62 mm步槍標(biāo)準(zhǔn)彈，HG202A-Ⅱ型電子測試儀，用2組錫箔作計(jì)時(shí)觸發(fā)器測試子彈的入射速度和剩余速度。

1.2 數(shù)值模擬主要指標(biāo)

3維編織復(fù)合材料的動態(tài)侵徹?cái)?shù)值模擬采用LS-DYNA程序，采用變時(shí)步長增量解法進(jìn)行了2組有限元模擬，其中：材料模型均采用塑性動力學(xué)材料本構(gòu)模型；準(zhǔn)細(xì)觀精細(xì)模型中，彈體、樹脂和纖維束的單元類型均采用八節(jié)點(diǎn)六面體實(shí)體單元。2組有限元模擬方法的主要不同指標(biāo)如表1所示。

表1 2組有限元模擬方法的主要不同指標(biāo)

1.3 結(jié)果

3維編織復(fù)合材料彈道侵徹試驗(yàn)結(jié)果與模擬結(jié)果如表2所示。

表2 3維編織復(fù)合材料彈道侵徹試驗(yàn)結(jié)果與模擬結(jié)果

2 回歸模型診斷

在3維編織復(fù)合材料彈道侵徹試驗(yàn)和數(shù)值模擬的對比分析中，檢驗(yàn)剩余速度和入射速度的線性關(guān)系是否相同。

令x表示入射速度，y表示試驗(yàn)中剩余速度，應(yīng)用最小二乘法回歸，可以得到剩余速度與入射速度的關(guān)系：

y=1.171 3x-124.020 7。

(1)

2.1 殘差分析

將回歸模型(1)分別應(yīng)用于試驗(yàn)數(shù)據(jù)和2組有限元模擬數(shù)據(jù)中的剩余速度與入射速度，得到的殘差如圖1所示。

由圖1可以看出：模擬數(shù)據(jù)1的殘差明顯具有增加的趨勢，模擬數(shù)據(jù)2的殘差大都小于0，初步可認(rèn)為回歸模型(1)不適合2組有限元模擬數(shù)據(jù)。

2.2 擬合優(yōu)度檢驗(yàn)與F檢驗(yàn)

對試驗(yàn)數(shù)據(jù)和2組模擬數(shù)據(jù)的回歸模型進(jìn)行擬合優(yōu)度檢驗(yàn)和F檢驗(yàn)，其結(jié)果如表3所示?？梢钥闯觯?)2組模擬數(shù)據(jù)回歸模型的擬合優(yōu)度檢驗(yàn)決定系數(shù)R2>0.95，說明擬合優(yōu)度較大，響應(yīng)變量對回歸變量的解釋程度較高；2)F檢驗(yàn)中的p值均小于0.01，說明總體上響應(yīng)變量與回歸變量之間的線性關(guān)系顯著。

表3 回歸模型檢驗(yàn)結(jié)果表

上述結(jié)果表明：回歸模型(1)能夠較好地解釋試驗(yàn)和數(shù)值模擬的入射速度與剩余速度的關(guān)系，且試驗(yàn)和數(shù)值模擬結(jié)果一致。

2.3 傳統(tǒng)檢驗(yàn)方法的不足

結(jié)合圖1和表3可知：擬合優(yōu)度檢驗(yàn)與F檢驗(yàn)的檢驗(yàn)結(jié)果一致，但與殘差圖結(jié)果相矛盾。進(jìn)一步分析擬合優(yōu)度檢驗(yàn)和F檢驗(yàn)可得到以下結(jié)論：

1) 擬合優(yōu)度檢驗(yàn)中的決定系數(shù)R2僅能說明列入模型的所有解釋變量對因變量的聯(lián)合影響程度，但對存在局部系統(tǒng)誤差的模型反應(yīng)不靈敏，不能說明哪些樣本點(diǎn)不適合此模型。

2)F檢驗(yàn)的理論基礎(chǔ)是基于正態(tài)分布，若模型檢驗(yàn)中使用F檢驗(yàn)，殘差需服從正態(tài)分布。

試驗(yàn)數(shù)據(jù)、模擬數(shù)據(jù)1和模擬數(shù)據(jù)2在回歸模型(1)的殘差直方圖如圖2所示?？梢钥闯觯涸诨貧w模型(1)下，試驗(yàn)數(shù)據(jù)的殘差具有正態(tài)分布特征，而模擬數(shù)據(jù)1的殘差則不具有正態(tài)分布特征，模擬數(shù)據(jù)2的殘差具有厚尾特征。采用Kolmogorov檢驗(yàn)對試驗(yàn)數(shù)據(jù)、模擬數(shù)據(jù)1和模擬數(shù)據(jù)2進(jìn)行正態(tài)檢驗(yàn)，得到的p值分別為0.025、0.575、0.817 5，說明試驗(yàn)數(shù)據(jù)的殘差服從正態(tài)分布，可進(jìn)行F檢驗(yàn)；模擬數(shù)據(jù)1和模擬數(shù)據(jù)2的殘差不服從正態(tài)分布，不可進(jìn)行F檢驗(yàn)。

3 分段趨勢診斷方法

趨勢性分析直觀上較容易理解，并不需要太多的數(shù)學(xué)和統(tǒng)計(jì)理論，將其引入到模型檢驗(yàn)中，可有效檢驗(yàn)?zāi)Ｐ偷挠行浴?/p>

3.1 分段趨勢統(tǒng)計(jì)量的構(gòu)造

檢驗(yàn)假設(shè)：

H0：任意2組殘差無上升或下降趨勢。

H1：存在2組殘差有上升或下降趨勢。

令{ε1,ε2,…,εn}表示殘差數(shù)據(jù)集，n0表示每組相互比較的數(shù)據(jù)集中的元素個(gè)數(shù)，取數(shù)對集 {(εi,εi+c),(εi+1,εi+1+c), …,(εi+n0-1,εi+n0+c-1)}，其中i=1,2,…,n-c，0

構(gòu)造分段趨勢統(tǒng)計(jì)量：

(2)

(3)

pst=P{T(1)≤t(1)}。

(4)

式中：t(1)為T(1)的觀察值。Tm能夠檢驗(yàn)任意2組殘差是否存在上升或下降趨勢，T(1)與pst均能檢驗(yàn)?zāi)Ｐ褪欠翊嬖诰植肯到y(tǒng)誤差：當(dāng)模型局部系統(tǒng)誤差不存在時(shí)，T(1)的值接近組中元素個(gè)數(shù)的1/2，此時(shí)模型擬合程度越好；pst為模型不存在系統(tǒng)誤差的概率。

3.2 分段趨勢統(tǒng)計(jì)量的分布律

定理1：連續(xù)型隨機(jī)變量X1、X2相互獨(dú)立且具有相同的分布，其密度函數(shù)均為f(x)，則P{X1≥X2}=P{X1≤X2}=1/2。

證明：由于X1、X2相互獨(dú)立，由卷積公式可得X1-X2的密度函數(shù)

(*)

令x+z=u，則有x=u-z，這時(shí)式(*)可改寫為

fX1-X2(z),

說明X1-X2的密度函數(shù)關(guān)于z軸對稱。有

P{X1≥X2}=P{X1≤X2}=1/2。

證畢。

(5)

3.3 n0對檢驗(yàn)的影響

表4為根據(jù)式(5)計(jì)算的不同元素個(gè)數(shù)n0下pst的值?？梢钥闯觯喝★@著性水平α=0.05，當(dāng)t=0，n0=6時(shí)，pst=0.016<α，說明當(dāng)n0過小時(shí)，趨勢檢驗(yàn)不夠靈敏，實(shí)際應(yīng)用中應(yīng)使n0>6。

分段趨勢診斷過程中，需依據(jù)試驗(yàn)條件和數(shù)據(jù)量選擇n0，然后依次計(jì)算Tm、T(1)與pst的值。當(dāng)顯著性水平為α?xí)r，若pst≤α，說明殘差具有上升或下降趨勢，需重新調(diào)整模型；若pst>α，說明殘差無上升或下降趨勢，回歸模型能夠解釋回歸變量與響應(yīng)變量的關(guān)系。

表4 不同元素個(gè)數(shù)n0下pst的值

4 分段趨勢診斷

為保證分段趨勢檢驗(yàn)的靈敏度，令組中元素個(gè)數(shù)n0=7，依據(jù)分段趨勢統(tǒng)計(jì)量的構(gòu)造方法，比較每組7對殘差的大小關(guān)系并計(jì)算T(1)和pst的值，如表5所示。

表5 分段趨勢診斷表

由表5可以看出：取顯著性水平α=0.05，模擬數(shù)據(jù)1中分段趨勢統(tǒng)計(jì)量T(1)=0，pst=0<α，試驗(yàn)數(shù)據(jù)和模擬數(shù)據(jù)2中分段趨勢統(tǒng)計(jì)量T(1)=1，pst=0.062 5>α，說明模擬數(shù)據(jù)1沒有通過分段趨勢檢驗(yàn)，殘差局部存在增加或減少的趨勢，而試驗(yàn)數(shù)據(jù)和模擬數(shù)據(jù)2通過了檢驗(yàn)。

5 結(jié)論

分段趨勢診斷方法能診斷出模型中的局部系統(tǒng)誤差，可更好地區(qū)分模型對數(shù)據(jù)的不合理解釋；分段趨勢診斷方法本質(zhì)上是非參數(shù)方法，沒有參數(shù)方法中殘差的正態(tài)假設(shè)或漸進(jìn)正態(tài)假設(shè)，更具普適性；在3維編織復(fù)合材料彈道侵徹試驗(yàn)與數(shù)值模擬結(jié)果的對比分析中引入分段趨勢分析方法，能夠更有效地診斷試驗(yàn)與數(shù)值模擬結(jié)果對回歸模型的同一性。