球形破片侵徹明膠的瞬時空腔簡化模型

田浩成， 盧芳云, 李志斌， 王 碩

(國防科技大學 文理學院，長沙 410005)

破片是造成人員傷亡的主要殺傷元，而球形破片作為一種常見的預制破片被廣泛應用于各種戰斗部中，因此研究高速球形破片對人體的殺傷效應無論是對武器戰斗部設計，還是對人體易損性分析都具有重要意義。創傷彈道學是研究諸如破片等高速投射物對生物機體創傷機理的學科，一般采用明膠等非生物靶標作為人體組織替代物來開展創傷彈道研究[1-2]。由于彈道明膠與人體肌肉組織具有相近的力學動態響應特性，本身透明，且價格相對低廉，實驗室常將其制成靶標來研究“子彈”的人體殺傷能力。破片進入彈道明膠后會沿著侵徹路徑形成空腔，該空腔的最大直徑可膨脹到破片直徑的十幾倍，而整個侵徹過程及空腔往復膨脹通常在毫秒量級，因此稱其為瞬時空腔[3-5]。醫學上對人體進行傷情評估發現，瞬時空腔是造成生物組織創傷的主要因素之一，因此瞬時空腔模型是后續開展人體易損性分析、人員毀傷效應評估等研究工作的關鍵。

隨著計算機技術的發展，各類計算軟件也得到了極大的完善，通過有限元方法以及SPH方法等無網格方法對實際問題進行仿真，為研究難以獲得試驗數據、試驗成本過高的復雜問題提供了新的解決途徑。在創傷彈道學研究領域，國內外研究工作者已經將此類技術手段廣泛應用到球形破片侵徹非生物靶標的仿真中。Dyckmans等[6]采用Autodyn對鋼球侵徹肥皂的過程進行了模擬，Frissane等[7]采用SPH方法建立了鋼球侵入彈道明膠的軸對稱模型，并且展示出在穿透深度和速度方面的顯著效果。吳建萍等[8]基于LS-DYNA軟件分析了初速700 m/s的Φ4 mm實心鋼球在侵徹明膠過程中的速度、能量和位移的變化規律。溫垚珂等[9]也基于LS-DYNA軟件分別采用Lagrange模型、SPH模型和ALE模型對鋼球侵徹明膠的同一工況進行仿真，通過與試驗結果的對比，得到對于此類高速投射物侵徹明膠的問題Lagrange算法最為合適的結論。羅少敏等[10]通過對不同材料的Φ4.8 mm的球形破片侵徹明膠的數值計算，分析得到了高速度、小質量的球形破片速度衰減快、能量傳遞大等規律，為未來戰斗部殺傷元的設計提供了參考依據。

本文基于LS-DYNA的仿真結果，分析了不同初始條件下鋼球在侵徹過程中的速度衰減、最大瞬時空腔演化規律，利用遺傳算法將仿真數據擬合得到簡化的速度衰減公式和瞬時空腔尺寸公式，通過與仿真模擬得到的瞬時空腔的比較驗證了簡化模型的有效性。

1 仿真模型

考慮到計算精度、計算時間等因素的影響，并根據溫垚珂等對不同計算模型結果的比較，本文選取Lagrange算法進行計算。為進一步減少計算時間，提高分析效率，同時考慮到原問題的軸對稱性，本文采用二維軸對稱模型對原問題進行簡化，建立的具體模型如圖1所示。

圖1中半圓形為球形破片的二維軸對稱模型，破片材料參數如表1。

表1 鋼球材料參數Tab.1 Material parameters of steel ball

球形破片在侵徹明膠過程中幾乎不發生變形，因此將其設為剛體(*MAT_RIGID)，以節省計算時間。破片直徑分別為4.8、6、7和8 mm，通過改變球形破片模型直徑的方式得到不同質量和不同等效面積的計算模型。彈丸的著靶速度分別為500、724、947和1 260 m/s，并以垂直于明膠靶面的方向入射，不考慮旋轉對彈丸的影響。結合以上設定，共得到16種計算工況分別對破片尺寸和著靶速度這兩個影響球形破片侵徹效應的主要參量進行分析。

圖1中矩形為明膠的二維軸對稱模型，其半徑15 cm，高30 cm。明膠的網格單元尺寸除中心區域外按照由著靶中心到外邊界成比例增長的方式劃分，劃分方式由粗糙到精細可得到三種網格，單元數目分別共計18 361、40 581、96 761，具體劃分如圖2。

考慮到人體厚度基本上不超過30 cm，因此在仿真中也只關心破片在侵徹30 cm以內的速度、穿透深度和最大瞬時空腔直徑等參量的規律。而在破片高速侵徹過程中，產生的高壓使明膠固體出現了流體的特性，因此可以采用LS-DYNA中自帶的多項式狀態方程(*EOS_LINEAR_POLYNOMIAL)和彈塑性材料模型(*MAT_ELASTIC_PLASTIC_HYDRO)來描述高壓狀態下明膠的動態響應，其材料及狀態方程所需參數如表2所示，其中C0、C1、C2和C3為線性多項式狀態方程中的材料常數。

表2 明膠材料參數Tab.2 Material parameters of gelatin

2 仿真結果及分析

采用粗、中、細三種網格對著靶速度為947 m/s的Φ4.8 mm鋼球侵徹明膠過程進行數值模擬，并與Wen等試驗得到的速度衰減曲線和侵徹深度曲線進行對比，結果如圖3所示。顯然，隨著網格數量的增加，數值模擬與試驗結果之間的差異也隨之減小，這說明數值計算對于網格具有收斂性，也驗證了仿真的可靠性。由于粗網格的結果與試驗結果之間存在著較大的差異，而細網格與中網格的計算結果實際上相差不算太大，中網格在侵徹深度和侵徹速度上與試驗結果分別相差不到5%和10%，卻大大減少了運算時間，提高了運算效率。因此本文選擇采用中等網格進行仿真計算。

圖4給出了Φ4.8 mm鋼球著靶速度為720 m/s時的最大瞬時空腔試驗結果與仿真結果的比較。通過對比，發現空腔的形狀和尺寸基本吻合，進一步驗證了仿真結果的準確性。

對不同工況的計算結果如圖5和圖6所示。由圖5和圖6不難看出，直徑大、入射速度高的鋼球穿透能力強，直徑小、著靶速度高的鋼球速度衰減快，并且侵徹速度的快速衰減主要在侵徹初期，鋼球的直徑和入射速度對之后的速度衰減影響不大。

圖7給出了較小尺寸(Φ6 mm)和較高初速(1 260 m/s)鋼球的速度衰減與其侵徹深度之間的關系曲線。圖7(a)表明同一直徑鋼球在穿透等厚度明膠靶標后的速度衰減程度大致相同，但最終衰減程度也會隨著入射速度的增加而小幅降低。而圖7(b)展示出同一速度下鋼球的直徑越小，其速度衰減越嚴重的現象。綜合上述兩種規律說明高速度、小尺寸鋼球穿透等厚明膠靶標的速度衰減會更加嚴重。

由表3可以發現，最大瞬時空腔直徑與鋼球直徑和侵徹速度成正比例。結合圖5和圖6得到的規律，瞬時空腔的尺寸與球形破片損失的動能也成正比例，而著靶速度高的破片速度喪失也會更大，在明膠中穿透的時間也更短，從而傳遞給明膠的動量也越大。

表3 不同工況下的最大瞬時空腔直徑Tab.3 The maximum temporary cavity diameter under different working conditions mm

侵徹過程中的瞬時空腔半徑隨速度衰減的變化規律如圖8所示。不難看出，與侵徹混凝土靶標過程類似，鋼球侵徹明膠過程也可分為“開坑段”和“穩定段”。在“開坑段”，瞬時空腔半徑差異不大，變化規律更接近二次曲線；隨著速度繼續衰減進入到“穩定段”，瞬時空腔尺寸與速度衰減程度基本呈線性關系且斜率隨入射速度增長而變大。

圖8 Φ6 mm鋼球不同著靶速度下計算結果Fig.8 Φ6 mm steel ball results at different velocities

3 遺傳算法擬合公式

仿真結果中已經得到了鋼球的直徑、質量和著靶速度對瞬時空腔的影響規律，為了便捷地對瞬時空腔輪廓進行描述，一個簡單直接的數學模型是人體易損性分析所需要的。人體損傷評估軟件HVAssess[11]中將球形破片形成的瞬時空腔近似成圓錐，這與真實情況有一定的差異。劉蘇蘇等[12]通過參數化建模提取瞬時空腔的方式進行了改進，雖保證了瞬時空腔的可靠性，但無法用于隨機破片的快速評估。因此，需要一個方便可靠的數學模型來描述損傷區域。

3.1 遺傳算法

遺傳算法(genetic algorithm,GA)是一種基于達爾文生物進化論中基因的遺傳機理進行搜索的優化算法，其通過直接對被優化的結構對象進行編碼，進而模擬人類種群中基因遺傳過程產生的選擇、交叉和變異現象[13]。在種群中，初始基因編碼經過一定代數的遺傳，子代在遺傳過程中往往會產生更優秀的基因，即種群通過自然進化的行為來適應外部環境，而我們可以用調整目標函數的方式更改其所需適應的外部環境，從而篩選出最優秀的子代，達到結構最優化的目的。

在MATLAB軟件中已經開發了GA工具包，按照上述遺傳算法的基本設置進行輸入，可以方便快捷地實現對最小值問題求解，這節省了在遺傳算法編程上面消耗的時間，僅需對適應度函數進行描述。

3.2 速度衰減公式

Dziemian定律[14]常被用來計算射彈侵徹人體組織的速度衰減，它假設作用于侵徹彈丸的阻力是彈丸和目標靶性能的函數，表示成如下形式

(1)

式中：m為彈丸質量；v為彈丸速度；t為時間；ρ為組織材料密度；A為彈丸呈現面積；μ為緊挨彈丸的液體層黏度,l為緊挨彈丸的液體層厚度，φ=μ/l;S為組織強度;gD、gV、gS為無量綱系數。

考慮如下變換

(2)

可將式(1)寫成

(3)

式中：σ為彈丸密度；a、b、c為常數；Q為等效長度，由Q=(m/σ)1/3得到。

通過求解式(3)，可得到經過侵徹深度Δx衰減后的速度公式

(4)

由于明膠與人體肌肉組織具有相似的性質，而式(4)的本質是動態阻力、黏性力和組織阻滯力對速度的函數，因此也應適用于明膠侵徹過程。基于遺傳算法對16種工況的速度衰減曲線進行參數優化，在國際單位制下得到a，b和c取值如表4所示。擬合效果如圖9所示。

表4 速度衰減系數Tab.4 Velocity attenuation coefficients

圖9 Φ7 mm鋼球速度衰減公式擬合結果Fig.9 Fitting results of velocity attenuation for Φ7 mm steel balls

3.3 瞬時空腔尺寸公式

前面分析了著靶速度和鋼球直徑對瞬時空腔尺寸的影響規律，發現瞬時空腔尺寸隨速度衰減可分為兩個階段，即侵徹初期對應的“開坑段”和侵徹進行到一定程度后對應的“穩定段”。考慮到兩個階段不同的特性，應分別選取相應的瞬時空腔公式描述侵徹過程隨侵徹速度衰減的瞬時空腔尺寸變化。

破片的呈現面積Ap是指侵徹過程中與靶標實際作用的破片截面積。呈現面積Ap的有效半徑為

(5)

式中，γ為形狀系數，球形破片的有效半徑等于其自身半徑，γ取為1.21。

對于瞬時空腔尺寸公式，考慮破片有效半徑、瞬時空腔半徑、侵徹速度、著靶速度、破片密度和明膠靶楊氏模量共6個影響因素，利用其無量綱組合建立瞬時空腔尺寸公式。

考慮到圖8中“開坑段”的曲線形狀更接近二次曲線，因此，侵徹初期瞬時空腔半徑的函數形式為

(6)

其中α、β、λ、μ和η為無量綱系數，在國際單位制下，系數的具體取值如表5。

表5 式(6)的系數取值Tab.5 Values of the coefficients of Eq.(6)

進入到“穩定段”后的瞬時空腔尺寸與速度衰減程度基本呈線性，且入射速度影響其斜率，因此，侵徹進行到一定程度后瞬時空腔半徑的函數形式為

(7)

其中ζ、ε和κ為無量綱系數，在國際單位制下，系數的具體取值如表6。

表6 式(7)的系數取值Tab.6 Values of the coefficients of Eq.(7)

通過對式(6)和式(7)分別擬合得到了兩段瞬時空腔半徑隨侵徹速度衰減的擬合曲線，將兩段曲線的交點作為兩段侵徹過程的劃分點，從而可以得到整體的擬合曲線，擬合結果如圖10所示。

由于某些工況下侵徹初期的曲線與穩定段曲線沒有交點，為了對瞬時空腔尺寸有較好的描述，需要選取適當的分割點并通過一定的過渡段將兩段曲線連接，如圖10中1 260 m/s擬合曲線中的過渡段為一條直線。

圖10 Φ6 mm鋼球對瞬時空腔尺寸公式的擬合結果Fig.10 Fitting results of temporary cavity size for Φ6 mm steel balls

3.4 瞬時空腔模型

結合3.2和3.3節得到的速度衰減公式和瞬時空腔尺寸公式，本文得到了一個計算瞬時空腔的簡易模型，可以對給定直徑、密度和著靶速度的球形破片較為快速地計算出其侵徹形成的瞬時空腔。

該簡易模型的計算流程可簡述為：第1步，將考慮的侵徹深度范圍劃分成一個個較小的間隔；第2步，根據已知破片條件，由3.2節中式(4)可計算得到破片侵徹一個間隔距離后的衰減速度；第3步，將破片的初始信息和衰減后的侵徹速度作為輸入，由3.3節中式(6)和式(7)分別得到兩個瞬時空腔半徑，通過判斷當前侵徹過程所處階段確定其中一個半徑表征該間隔上的瞬時空腔尺寸；第4步，重復第2、3步計算直至最后一個間隔；最終，經過上述四步的計算，得到了各間隔對應的瞬時空腔半徑，這些半徑實際上描述了破片侵徹明膠造成的瞬時空腔輪廓。

圖11為利用該模型計算得到的結果與LS-DYNA數值模擬結果的對比圖，可以發現本文建立的簡易模型與仿真軟件的結果吻合度較好，而在侵徹初始階段的明膠受擠壓和回彈作用在仿真中產生了鋸齒狀的空腔邊界形狀，但模型計算結果與鋸齒外沿相貼，能夠較為準確地描述空腔區域輪廓。通過16種擬合工況之外的仿真結果和計算結果的比對，也驗證了此模型的有效性。

圖11 初速800 m/s的Φ6.5 mm鋼球最大瞬時空腔模型計算結果Fig.11 Calculation results of the maximum temporary cavity model for Φ6.5 mm at a velocity of 800 m/s

4 結 論

本文基于非線性顯式計算有限元軟件LS-DYNA，對不同直徑和著靶速度的鋼球侵徹明膠過程進行了數值模擬，通過文獻中的試驗結果驗證了二維軸對稱模型的可靠性。經整理分析得到了球形破片直徑和著靶速度對侵徹速度衰減和瞬時空腔尺寸的影響規律，發現高速度、大直徑的球形破片具有更強的穿透能力和更大的空腔尺寸，而高速度、小尺寸的球形破片穿透相同厚度明膠靶標的速度衰減更嚴重。

基于仿真結果，利用遺傳算法對速度衰減公式和瞬時空腔尺寸公式進行擬合，得到了簡易瞬時空腔模型。通過與仿真結果提取的二維瞬時空腔輪廓進行對比，驗證了該模型的有效性。該模型極大縮短了計算瞬時空腔的時間，可提高相應毀傷評估的效率。