網格劃分對沖擊波波形的影響分析

索 強，徐 鵬，尤文斌

(中北大學 a.理學院；b.儀器與電子學院，太原 030051)

沖擊波超壓是武器彈藥最主要的毀傷元，對絕大部分目標起主要毀傷作用，沖擊波效應一般以沖擊波超壓峰值、正壓持續時間和比沖量進行表征，通常關注的是沖擊波各項參數中的某一項，而沖擊波的整體波形往往被忽略。但事實上，沖擊波波形對毀傷過程的研究具有重要作用，其直接與毀傷威力大小緊密相關。并且在對沖擊波仿真計算分析時，為簡化計算，往往將沖擊波脈沖等效為三角波或矩形波，沖擊波波形也是壓力簡化參數設置的重要依據[1]。對爆炸沖擊波進行數值模擬時，有限元網格尺寸對數值仿真的精度具有很大的影響[2-5]。目前，已有部分學者提出或研究過爆炸波模擬中的網格尺寸效應。崔瑩等通過開展鋼管混泥土柱爆炸試驗，分析不同空氣和炸藥網格尺寸對沖擊波傳播及數值分析結果的影響，最終確定網格尺寸為20 mm符合折合距離不超過1.1 m·kg-1/3試驗條件下的爆炸數值模擬要求[6]；Luccioni等[7]研究了利用流體力學軟件模擬預測爆炸荷載時的網格尺寸效應，認為100 mm的網格尺寸就可以較為精確的模擬爆炸荷載的傳播規律，而較粗的網格尺寸則僅僅可以用來定性的模擬爆炸荷載在城市復雜環境中的傳播規律；石磊等[8]研究了與炸藥當量相關的網格劃分方法，并應用有限元軟件LS-DYNA進行了計算對比，結果表明按照3/80倍炸藥體邊長進行網格劃分在計算精度與計算效率上取得了很好的平衡；可見，不同研究者在研究中使用的網格尺寸有較大差別且只關注沖擊波的個別參數信息，對沖擊波的波形關注不夠，模擬結果的精度不夠。已有文獻對網格尺寸效應的研究往往針對某一特定情況，故提出的網格尺寸的使用有較大的局限性。而且，由于爆炸問題的復雜性，多數爆炸沖擊波的網格效應的研究均集中在理論分析和數值模擬的對比上，缺乏有效的實驗驗證。

本研究運用ANSYS/LS-DYNA軟件對炸藥在空氣中的爆炸進行模擬，對比分析不同網格尺寸對沖擊波波形的影響，得出與經驗計算誤差小于5%，并且計算時間合理的最佳網格劃分尺寸。并進行爆炸實驗，用本文得出的網格尺寸對實驗進行模擬，與實驗結果進行比較，模擬結果與實驗結果的誤差小于5%，證明本文結果的正確性與適用性。

1 爆炸沖擊波及其經驗計算方法

自由空氣中的理想沖擊波超壓曲線，如圖1所示。在沖擊波到達之前，該處的壓力等于大氣壓力P0，壓力經過時間Tc由大氣壓力突躍至最大值。壓力最大值與P0的差值，通常稱為入射超壓峰值P∞。波陣面通過后壓力即迅速下降，經過時間Td，壓力經指數衰減到大氣壓力并繼續下降，直至出現負超壓峰值，在一定時間內又逐漸地回升到大氣壓力[9]。

圖1 理想沖擊波超壓曲線

目前常用的方法，是用比例距離表達沖擊波的各種參數。比例距離為Z=R/W1/3。其中，R為測點與爆心之間的距離，W為等效TNT藥量[10-13]。參照文獻[8]和文獻[9]，選取以下經驗公式。

Henrych根據實驗提出沖擊波的超壓峰值(MPa)表達式為：

(1)

Henrych根據實驗提出沖擊波超壓峰值(MPa)的表達式為:

0.05≤Z≤3

(2)

Chengqing Wu和Hong Hao給出了沖擊波超壓上升段(從大氣壓上升到峰值)持續時間的表達式為：

TC=0.001 9Z0.13C

(3)

薩多夫斯基建議爆炸沖擊波超壓持續時間為：

(4)

2 單元網格尺寸對爆炸沖擊波波形的影響

分別取0.84 kg，13.12 kg，105 kg炸藥作為計算模型，利用計算模型的對稱性取1 /8模型進行計算，對稱面采用對稱邊界，非對稱面采用透射邊界，炸藥和空氣均采用LS-DYNA Explicit 3D SOLID164模型，炸藥邊長為d，d分別取80 mm，200 mm，400 mm，模型如圖2。當比例距離小于1.3時,網格尺寸對沖擊波波形計算結果影響較大，這里計算并比較了比例距離小于1.3時的沖擊波結果。

圖2 計算模型

2.1 材料的狀態方程和參數

空氣采用*MAT_NULL材料模型，壓力P用線性多項式狀態方程*EOS_LINEAR_POLYNOMIAL來描述[14]。

(6)

式中：E為單位初始體積內能；ρ為空氣質量密度；ρ0為參考質量密度；線性多項式狀態方程描述空氣時遵守γ定律；V為初始相對體積；C0、C1、C2、C3、C4、C5、C6為實常數。材料參數如表1所示。

炸藥采用ANSYS/LS-DYNA中的*MAT_HIGH_EXPLOSIVE_BURN材料模型，定義壓力P為相對體積初始能量E的函數[15]。

(7)

式中:V為相對體積；E為單位體積內能；A、B、R1、R2、ω為表征炸藥材料特性的常數。JWL狀態方程參數及炸藥參數如表2所示。

通過對爆炸沖擊波已有研究成果進行對比分析，參照文獻[6-8]和文獻[11-13]，發現當沖擊波超壓峰值的誤差較小時，網格劃分在炸藥體邊長的1/8至炸藥體邊長的3/80之間。所以，本文模擬時的單元劃分尺寸w分別采用1/4、1/8、1/16、1/32倍炸藥體邊長(d)。

表1 空氣參數

表2 炸藥參數

2.2 模擬結果及分析

圖3(a)，圖3(b)，圖3(c)分別是d=80 mm，200 mm，400 mm時，分別按照1/4、1/8、1/16、1/32倍炸藥體邊長d進行單元劃分的超壓峰值曲線，圖4、圖5分別為圖3相對應的正壓時間曲線和超壓上升時間曲線?？梢钥闯觯孩?隨著比例距離的增大，超壓峰值的模擬結果與式(1)之間的相對誤差越來越小；② 相同網格劃分時，炸藥的當量越大，在相同的比例距離，模擬結果的超壓峰值相同，而超壓上升時間和正壓時間與式(3)和式(4)的誤差越大；③ 爆炸沖擊波的超壓峰值，正壓時間和超壓上升時間隨著單元尺寸的減小而越來越接近經驗公式值,沖擊波的從平滑曲線變為上升和下降明顯的三角形曲線；④ 在比例距離為0.1～0.2時，網格劃分為1/16d的沖擊波各項參數的誤差都大于5%，發現當進一步劃分網格至1/32d時，0.1-0.2比例距離的沖擊波各項參數很接近經驗公式，誤差小于5%。導致這種情況的主要原因可能是當網格劃分的單元比較大時，提取點距爆心的比例距離小于0.1，這導致了計算數值提取點與相應比例距離的偏差，而在0.1～0.2比例距離內沖擊波高頻分量作用顯著，受爆心距影響非常明顯，這種微小的偏差足以導致計算結果的巨大差異。當網格劃分的單元比較小時，這種偏差也就減小了，計算結果將更為精確。

進一步探究網格劃分對沖擊波結果的影響，對比分析不同網格劃分對沖擊波模擬結果和經驗公式的誤差，確定準確的單元網格劃分方法。圖6為比例距離為0.2時，炸藥體邊長d分別為80 mm，200 mm，400 mm時，不同網格劃分的沖擊波超壓峰值和式(1)的誤差曲線。在網格劃分小于1/27d以后，誤差均小于5%，認為是符合模擬的精度要求的。圖7為比例距離為0.3時，d分別為80 mm，200 mm，400 mm時，不同網格劃分的正壓時間模擬結果和式(4)的誤差曲線。在網格劃分小于1/29d以后，誤差均小于5%符合模擬精度要求。

圖3 不同當量炸藥不同網格尺寸的超壓峰值曲線

圖4 不同當量炸藥不同網格尺寸的正壓時間曲線

圖5 不同當量炸藥不同網格尺寸的超壓上升時間曲線

圖6 不同網格劃分的超壓峰值和經驗公式的誤差曲線

圖7 不同網格劃分的正壓時間和經驗公式的誤差曲線

隨著炸藥當量的增大，相同網格劃分的正壓時間模擬結果和經驗公式的誤差逐漸減小。導致這種情況的主要原因可能是當炸藥當量比較大時，相同比例距離的位置距離炸藥的實際距離較遠，當單元比較小時，計算數值提取點與相應比例距離的偏差減小了，計算結果將更為精確。

圖8為比例距離為0.2時，d分別為80 mm，200 mm，400 mm時，不同網格劃分的超壓上升時間曲線和式(3)的誤差曲線。在網格劃分小于1/30d以后，誤差均小于5%，符合模擬的精度要求。但隨著炸藥當量的增大，相同網格劃分的超壓上升時間的模擬結果和經驗公式的誤差逐漸增大。于是進一步研究網格劃分為炸藥體邊長的1/32時，不同炸藥當量對模擬結果的影響，沖擊波超壓上升時間的模擬結果和式(3)的誤差曲線如圖9所示。可以得出：炸藥質量在1 000 kg以內，網格劃分為炸藥體邊長的1/32時，模擬得到的沖擊波各項參數都符合精度要求，沖擊波的波形也符合精度要求。

圖8 不同網格劃分的超壓上升時間和經驗公式的誤差曲線

圖9 不同炸藥當量對模擬結果和經驗公式的誤差曲線

3 爆炸實驗及對實驗的模擬

3.1 實驗樣品和靶場布置

實驗采用彈體進行試驗，由于沖擊波正壓時間一般在從幾毫秒到數十毫秒，上升前沿僅為數微秒，測試裝置采用高頻響測試系統，主要由壓力傳感器和瞬態數據采集儀組成。傳感器選用美國PCB公司的壓電式傳感器113 A。實驗樣品為彈藥，等效TNT藥量m為515 kg，裝藥密度為1.64 g·cm-3，裝藥長徑比為1∶1。傳爆藥為JH-14，質量為30 g，為減小沖擊波能量對地作用損耗，采用8號銅雷管從裝藥底端起爆。

實驗時，將彈藥放在距地面3 m的支架上，以裝藥在地面上的垂直投影點為爆心，在爆心區周圍200 m內無建筑物、較開闊的地面上布置地面傳感器，如圖10所示。壓力傳感器安裝在Φ200 mm的鋼質基礎上，傳感器敏感面與基準地面平齊安裝并用土夯實進行防護，防止熱及沖擊振動干擾。實驗現場布置如圖11所示。

圖10 傳感器

圖11 實驗現場布置示意圖

3.2 對實驗的模擬

土選取SOIL_AND_FOAM_FAILURE材料模型，該模型在某些方面具有流體性質，其應用于土或泡沫被限制在結構中或有幾何邊界存在的情況下。壓力是正壓縮，在負壓縮情況下，體積應變是相對體積的自然對數，相對體積是計算開始時當前體積與初始體積之比。其塑性屈服極限函數φ根據應力偏量第二不變量J2描述[16]。

φ=J2-(a0+a1p+a2p)

(8)

式中：J2=SijSij/2，a0，a1，a2為常數；p為壓力。主要材料參數參見表3。

表3 主要材料參數

取515 kg炸藥作為計算模型，利用計算模型的對稱性取1/8模型進行計算，對稱面采用對稱邊界，非對稱面采用透射邊界。炸藥和空氣均采用LS-DYNA Explicit 3DSOLID164模型，炸藥邊長為680 mm，炸藥離地3 m，利用模型的對稱性取1/8模型進行計算，網格劃分為炸藥體邊長d的1/32，模型如圖12。

圖12 實驗的模擬模型

3.3 模擬結果分析

爆炸實驗的測點1處的實驗結果如圖13所示，圖14為對應的模擬結果，表4列出了本次實驗3個測點處沖擊波的實驗值，模擬值和計算值。

圖13 實驗結果

圖14 模擬結果

表4 沖擊波的實驗值，模擬值和計算值的偏差分析

從表4可以看出，網格劃分為炸藥體邊長d的1/32時，模擬結果的各項參數對實驗值和計算值的誤差都小于5%，實驗結果和模擬結果的沖擊波波形基本一致，能夠滿足模擬的精度要求。

4 結論

不同網格尺寸對沖擊波波形有影響。網格尺寸過大時，壓力小于真實值，隨著網格的不斷減小，超壓峰值逐漸增大，接近真實值。在網格尺寸小于炸藥體邊長的1/27后，超壓峰值符合精度要求。網格尺寸越小，正壓時間越小，越符合真實值。在網格尺寸小于炸藥體邊長的1/29后，正壓時間符合精度要求。網格尺寸越小，超壓上升時間越小，越符合真實值。在網格尺寸小于炸藥體邊長的1/32后，超壓上升時間符合精度要求。數值模擬與實驗值和理論值的誤差隨著網格的減小而不斷減小，沖擊波波形從平緩曲線變為上升段和下降段明顯的三角形曲線。實驗結果表明，在炸藥小于1 000 kg，單元網格劃分小于炸藥體邊長的1/32時，模擬結果與計算結果，模擬結果和實驗結果的各項參數誤差都小于5%。