柱形裝藥空中運動爆炸沖擊波荷載計算模型

關(guān)鍵詞：柱形裝藥；空中運動爆炸；峰值超壓；最大沖量；計算模型

戰(zhàn)斗部等彈載柱形裝藥通常在高速運動狀態(tài)下爆炸（簡稱動爆），爆炸沖擊波是其主要毀傷元。受裝藥運動影響，動爆沖擊波場的分布和荷載特征與靜止爆炸（簡稱靜爆）不同。研究柱形裝藥的空中動爆沖擊波荷載對戰(zhàn)斗部爆炸威力計算、建筑目標(biāo)的毀傷評估及其防護設(shè)計等具有重要意義。

受限于試驗技術(shù)，多數(shù)研究聚焦于裝藥靜爆工況，UFC3-340-02[1]、TM5-855-1[2]和GB6722-2104[3]等規(guī)范給出了球形裝藥空中靜爆的入射沖擊波峰值超壓與比例距離Z（Z=R/W1/3，R為爆炸中心到目標(biāo)的距離，W為等效TNT裝藥當(dāng)量）的關(guān)系。Stoner等[4]、Brode[5]、Baker[6]、Henrych[7]和Mills[8]基于爆炸試驗或數(shù)值仿真提出了球形裝藥空中靜爆的入射沖擊波峰值超壓計算公式。受起爆方式（單端、中心和兩端起爆等）、長徑比L/D（L和D分別為裝藥長度和直徑）和方位角α（α為空中任一點到裝藥中心的直線與柱形裝藥軸向的夾角，規(guī)定軸向為0°方向）等因素影響，柱形裝藥的沖擊波荷載較球形裝藥復(fù)雜。Anastacio等[9]和Knock等[10-12]基于試驗提出了單/雙端起爆方式下考慮長徑比的柱形裝藥軸/徑向入射沖擊波峰值超壓和比沖量計算公式。Xiao等[13]基于數(shù)值仿真研究了起爆方式（中心、單端和雙端起爆）對柱形裝藥空中靜爆入射沖擊波場的影響，結(jié)果表明，單端起爆產(chǎn)生的入射沖擊波峰值超壓最大。Gao等[14]基于數(shù)值仿真提出了考慮長徑比（1≤L/D≤8）、比例距離（0.3m/kg1/3≤Z≤10.0m/kg1/3）和方位角（0°≤α≤360°）的入射沖擊波峰值超壓和比沖量計算公式。王明濤等[15]進一步開展了上述3種起爆方式下柱形裝藥在0.25≤L/D≤10和0.3m/kg1/3≤Z≤15.0m/kg1/3范圍內(nèi)的空中靜爆入射和反射沖擊波場數(shù)值仿真研究，提出了入射沖擊波峰值超壓和最大沖量的計算公式，以及反射沖擊波峰值超壓和等效持時的計算方法。

針對空中動爆工況，美國彈道研究實驗室（ballisticsresearchlaboratory，BRL）[16-18]開展了球形Pentolite和B炸藥的動爆試驗，發(fā)現(xiàn)動爆沖擊波入射峰值超壓隨著運動速度的增大和方位角的減小而增大，提出了動爆入射沖擊波峰值超壓的理論計算模型，并得到了試驗驗證。蔣海燕等[19]和陳龍明等[20]基于AUTODYN有限元軟件分析了球形裝藥空中動爆的沖擊波場特征，發(fā)現(xiàn)沖擊波場的時空分布受裝藥運動速度影響，入射沖擊波峰值超壓與方位角呈余弦關(guān)系，并建立了入射沖擊波峰值超壓的工程計算模型。聶源等[21]采用SPEED軟件開展數(shù)值仿真，給出了球形裝藥基于靜爆的動爆峰值超壓修正因子。Xu等[22]采用AUTODYN軟件進行數(shù)值仿真，給出了球形裝藥空中動爆的峰值超壓環(huán)/徑向分布系數(shù)，建立了球形裝藥空中動爆入射峰值超壓的工程計算模型。周至柔等[23]討論了裝藥運動速度對球形裝藥動爆沖擊波場分布的影響，發(fā)現(xiàn)沖擊波場會沿運動方向前移，移動距離與運動速度線性相關(guān)。Ma等[24]開展了柱形帶殼裝藥的空中動爆試驗和數(shù)值仿真研究，結(jié)果表明，高速運動對速度正向的入射峰值超壓起增強作用，正向的峰值超壓較負(fù)向更高，且波陣面到達時間更短。Chen等[25]基于數(shù)值仿真研究了柱形裝藥空中動爆沖擊波的傳播規(guī)律，建立了L/D=4.5的柱形裝藥空中動爆入射沖擊波峰值超壓的工程計算模型（Z≤1.5m/kg1/3）。王振寧等[26]對柱形裝藥近地動爆的沖擊波傳播規(guī)律開展了數(shù)值仿真，結(jié)果表明，裝藥運動增大了裝藥周向的沖擊波馬赫桿高度，且隨著速度的增大，運動方向的峰值超壓線性增大。相較于靜爆，動爆沖擊波場整體前移并產(chǎn)生波陣面變形，波陣面強度與運動方向呈余弦關(guān)系。

目前，針對裝藥空中動爆沖擊波荷載特征的研究存在以下不足：（1）相關(guān)研究多針對球/柱形裝藥靜爆工況[1-15]，對于動爆工況的研究多局限于球形裝藥[16-23]，由于運動方向的動爆沖擊波強度增強，導(dǎo)致基于靜爆計算公式得到的沖擊波荷載偏小，不利于工程防護設(shè)計；（2）少量針對柱形裝藥動爆的數(shù)值模擬缺乏動爆試驗的驗證[25-26]，且忽略了運動引起的擾動壓力場[24-26]，結(jié)論的可靠性有待商榷；（3）尚未建立考慮裝藥運動速度、長徑比、方位角和比例距離等影響因素的柱形裝藥空中動爆入射和反射沖擊波荷載計算模型，針對真實戰(zhàn)斗部等柱形裝藥（長徑比6≤L/D≤10、尾端起爆）的研究比較匱乏。針對上述不足，本文中，基于有限元分析程序AUTODYN提出一種有限元分析方法，并對柱形裝藥的空中靜爆試驗[15，27-28]以及球形裝藥的空中動爆試驗[16-18]進行數(shù)值模擬，通過對比試驗和數(shù)值模擬結(jié)果驗證有限元分析方法的可靠性；然后，基于驗證的有限元分析方法對柱形裝藥空中靜、動爆工況開展數(shù)值仿真計算，定量分析動爆沖擊波場的分布特征、入射和反射沖擊波荷載；最后，提出考慮裝藥運動速度、長徑比、方位角和比例距離的柱形裝藥空中動爆入射和反射沖擊波峰值超壓和最大沖量計算模型，通過簡化戰(zhàn)斗部柱形裝藥空中動爆的數(shù)值仿真驗證該計算模型的適用性。

1有限元分析方法

基于AUTODYN軟件提出柱形裝藥空中靜、動爆的有限元分析方法。鑒于柱形裝藥動爆試驗數(shù)據(jù)的缺乏，對柱形TNT和乳化炸藥空中靜爆試驗[15，27-28]以及球形Pentolite和B炸藥空中動爆試驗[16-18]開展數(shù)值模擬研究，驗證有限元分析方法的可靠性。

1.1有限元模型和分析方法

圖1給出了柱形裝藥空中靜、動爆的二維軸對稱有限元模型，對稱軸為裝藥軸。數(shù)值模擬過程中，為平衡計算效率和精度，并解決填充炸藥無法模擬裝藥高速運動及其引起的空氣擾動壓力場等問題，提出“三階段”有限元分析方法：階段Ⅰ，采用精細(xì)化網(wǎng)格建立局部空氣域模型（局部模型），包含空氣和彈體，彈體（初速度為v0）在空中運動，運動方向與裝藥軸向平行；階段Ⅱ，待彈體運動引起的擾動壓力場穩(wěn)定后（1ms）刪除彈體，通過填充命令（Fill）建立裝藥模型，并設(shè)置起爆點和初速度v0，重啟動模型進行計算并生成映射文件（.fil文件）；階段Ⅲ，采用大尺寸網(wǎng)格建立整體空氣域模型（全域模型），通過重映射策略導(dǎo)入階段Ⅱ的結(jié)果并進行后續(xù)計算，沿各方位角在不同距離處布置測點，以獲取自由場入射沖擊波的壓力時程曲線。真實戰(zhàn)斗部的長徑比取6～10，起爆方式為尾端起爆，馬赫數(shù)Ma為0～4。靜爆工況直接從階段Ⅱ開始，空氣域的邊界均設(shè)為自由流出邊界，即可以實現(xiàn)物質(zhì)和能量自由流出及壓力自由傳播的無反射邊界。

有限元模型中，階段Ⅰ為彈體運動，彈體采用AUTODYN內(nèi)置的STEEL4340材料和Lagrange網(wǎng)格，空氣采用Air材料和Euler網(wǎng)格。階段Ⅱ為裝藥爆炸，通過Fill命令填充的裝藥采用TNT材料。Air和TNT材料的狀態(tài)方程分別采用理想氣體和Jones-Wilkins-Lee狀態(tài)方程：

式中：p1和p2分別為空氣和爆轟產(chǎn)物壓力；γ為空氣絕熱指數(shù)，γ=1.4；ρ為空氣密度，ρ=1.225kg/m3；e0為空氣內(nèi)能，e0=2.068×105J/kg；A、B、R1、R2和ω分別為TNT炸藥的材料常數(shù)，其中A=373.77GPa，B=3.75GPa，R1=4.15，R2=0.9，ω=0.35；V為爆轟產(chǎn)物相對體積；E為炸藥初始體積內(nèi)能。

為確保數(shù)值仿真的可靠性，對有限元模型進行網(wǎng)格收斂性分析。以1kg球形TNT裝藥中心起爆的空中靜爆工況為例，分別對比局部模型和全域模型中不同網(wǎng)格尺寸下入射沖擊波的壓力時程曲線，確定兩類模型中的網(wǎng)格尺寸。圖2顯示了局部模型在網(wǎng)格尺寸為0.8～4mm時和全域模型在網(wǎng)格尺寸為8～25mm時的空氣壓力時程曲線。可以看出，局部和全域模型的壓力時程曲線分別在網(wǎng)格尺寸小于1和10mm時趨于收斂，因此兩類模型的網(wǎng)格尺寸分別取1和10mm。

1.2試驗驗證

王明濤等[15]開展了4kg柱形TNT裝藥單端起爆的空中靜爆試驗。圖3（a）～（b）給出了試驗設(shè)置和傳感器布置情況，其中柱形裝藥長度（L）為224mm，直徑（D）為120mm，裝藥軸平行于地面并距離地面1.0m。在距裝藥中心徑向4.85m處以及裝藥軸向4.96和7.75m處，等高布置自由場壓力傳感器Pi1、Pi2和Pi3，測量入射沖擊波的超壓時程曲線。根據(jù)試驗建立局部和全域有限元模型，如圖3（c）所示，材料模型參數(shù)和網(wǎng)格尺寸見1.1節(jié)。

圖4顯示了3個測點處入射沖擊波的超壓時程曲線，試驗[15]和數(shù)值模擬結(jié)果符合得較好。在測點Pi2和Pi3處的模擬結(jié)果較好地復(fù)現(xiàn)了裝藥形狀引起的雙波峰現(xiàn)象（柱形裝藥的端部主波和側(cè)面主波依次通過測點使超壓時程曲線產(chǎn)生2次波峰）。測點Pi3處首個入射沖擊波峰值超壓的仿真結(jié)果比試驗結(jié)果略小，這是因為測點Pi3距裝藥較遠(yuǎn)，試驗中測得的超壓經(jīng)地面反射后增強，而數(shù)值仿真僅關(guān)注自由場超壓，未考慮地面反射的影響。

Simoens等[27-28]開展了長徑比（L/D）為1和8.3時的柱形乳化炸藥空中靜爆試驗，裝藥質(zhì)量（W）為1.6kg，距離地面1.5m。如圖5（a）～（b）所示，當(dāng)L/D=1時，炸藥水平懸于空中，分別在裝藥中心和單端（規(guī)定為裝藥的0°方向）起爆，在裝藥同一水平面0°～180°范圍內(nèi)的9個方位角且距離裝藥0.8m處布置自由場超壓傳感器。圖5（b）也給出了L/D=8.3時柱形裝藥中心起爆的爆炸試驗俯視圖，裝藥懸于空中，在裝藥中心的等高水平面內(nèi)且比例距離為0.5～5m/kg1/3范圍內(nèi)布置自由場超壓傳感器。乳化炸藥關(guān)于超壓和沖量的等效TNT當(dāng)量系數(shù)分別為1.0和0.7[27-28]。試驗的有限元模型如圖5（c）所示，材料參數(shù)和網(wǎng)格尺寸見1.1節(jié)。

圖6（a）～（d）對比了L/D=1時試驗[27-28]和模擬的入射沖擊波峰值超壓pso和比沖量iso?？梢钥闯?，在中心起爆工況（見圖6（a）～（b））中，方位角α為70°和110°時，入射峰值超壓模擬值與試驗值的相對誤差δ為34.31%，其余方位角的δ均小于23.67%，而比沖量在各方位角的相對誤差均小于16.10%。在單端起爆工況（見圖6（c）～（d））中：當(dāng)方位角為70°時，峰值超壓的相對誤差最大（47.52%），其余方位角的δ均小于28.89%；方位角為160°時，比沖量的相對誤差最大（55.20%），其余方位角的δ均小于21.75%。圖6（e）對比了柱形裝藥（L/D=8.3）中心起爆時沖擊波峰值超壓-比例距離關(guān)系（pso-Z）的試驗和數(shù)值模擬結(jié)果，二者符合良好。

如圖7（a）～（b）所示，美國BRL[16-18]開展了球形Pentolite和B炸藥的空中動爆試驗，裝藥質(zhì)量（W）為0.1701kg。其中，Pentolite炸藥通過炮管以518.16～579.12m/s的速度發(fā)射并在空中預(yù)定位置引爆，起爆點位于球心。以爆心位置為原點、運動方向為0°方向，在爆心等高平面內(nèi)的15°～165°范圍內(nèi)6個方位角且距爆心0.826m處布置壓力傳感器。相應(yīng)地，B炸藥以604.11和534.31m/s的速度發(fā)射，起爆點位于球心，在預(yù)定爆心15°、45°和105°方向且距爆心0.719和0.826m處布置壓力傳感器，此外還進行了B炸藥的靜爆試驗。圖7（c）給出了“三階段”有限元模型，相關(guān)參數(shù)見1.1節(jié)。

對于Pentolite炸藥，共進行了7次動爆試驗，編號分別為No.534、No.549、No.576、No.582、No.587、No.588和No.590。圖8比較了部分動爆試驗中試驗、理論計算和數(shù)值模擬的入射沖擊波峰值超壓，可以發(fā)現(xiàn)，數(shù)值模擬與試驗結(jié)果吻合得較好。當(dāng)模擬和試驗的相對誤差δs偏大時（如No.582中44.3°、No.587中105.3°和No.590中166.3°處，δs分別為30.06%、32.80%和46.94%），模擬與理論模型的相對誤差δt均較?。é膖分別為4.15%、2.38%和1.21%）。受環(huán)境和測試設(shè)備等因素影響，部分爆炸試驗數(shù)據(jù)可能出現(xiàn)偏差，而數(shù)值模擬與理論模型的結(jié)果接近，可以認(rèn)為數(shù)值模擬結(jié)果可信。圖9對比了模擬和試驗得到的B炸藥靜、動爆試驗[17-18]入射沖擊波峰值超壓和最大沖量，二者總體吻合較好，峰值超壓和最大沖量的最大相對誤差出現(xiàn)在v0=604.11m/s和α=45°時以及靜爆時（分別為43.59%和25.03%），其余相對誤差分別小于35.74%和24.38%。

對于柱形裝藥空中靜爆和球形裝藥空中動爆工況，試驗和模擬的入射沖擊波超壓時程曲線、峰值超壓和爆炸沖量符合良好，驗證了“三階段”有限元分析方法分析柱形裝藥空中動爆沖擊波荷載的可靠性。

2運動裝藥爆炸沖擊波場分布特征

采用1.1節(jié)的有限元分析方法模擬100組1kg柱形裝藥的空中動爆工況，其中，裝藥運動的馬赫數(shù)范圍為0～4，裝藥長徑比范圍為6～10，起爆方式為尾端起爆，以揭示入射沖擊波場的分布特征，分析裝藥運動速度和長徑比對爆炸荷載的影響。

2.1運動擾動壓力場

裝藥在空中運動時會引起周圍空氣擾動，產(chǎn)生擾動壓力場，如圖10（a）所示。擾動壓力場與爆炸沖擊波場疊加形成復(fù)雜的耦合擾動沖擊波場，現(xiàn)有研究[19-26]通常忽略擾動壓力場對動爆沖擊波場的影響。為解決這一問題，模擬柱形裝藥（L/D=6）在Ma為1和4時是否考慮擾動壓力場的動爆工況，其中，考慮擾動壓力場時采用完整的“三階段”有限元模型，不考慮時則忽略階段Ⅰ（裝藥運動），有限元模型如圖10（b）所示，參數(shù)見1.1節(jié)。在裝藥的方位角為0°～180°、比例距離為0.3～0.5m/kg1/3范圍內(nèi)布置測點記錄超壓時程曲線。

圖11顯示了方位角為180°時各比例距離處擾動壓力場對超壓時程曲線的影響?？梢钥闯觯篗a為1時，Z取0.3、0.4和0.5m/kg1/3時考慮和不考慮擾動壓力場的入射沖擊波峰值超壓相對誤差δd分別為9.08%、6.46%和2.94%；Ma為4時，對應(yīng)的相對誤差分別為29.46%、14.64%和5.78%。在不同的方位角和比例距離處，考慮和不考慮擾動壓力場的入射沖擊波峰值超壓相對誤差列于表1。可以看出，擾動壓力場改變了裝藥爆炸后的沖擊波場分布，對運動方向的峰值超壓影響極大，且上述影響隨著比例距離的減小和運動速度的增大而增強。因此，在動爆工況分析中應(yīng)當(dāng)考慮擾動壓力場的影響。

2.2沖擊波場分布特征

圖12顯示了L/D=6、Ma=4、起爆方式為尾端起爆時柱形裝藥在不同時刻的沖擊波壓力云圖。引爆后0.01ms（見圖12（a）），爆轟波在裝藥內(nèi)部由爆炸區(qū)向未爆區(qū)傳播。引爆后0.1ms（見圖12（b）），爆轟波轉(zhuǎn)為沖擊波并向周圍空氣域傳播。根據(jù)沖擊波的形狀，將其分為端部主波（前、后端主波）、側(cè)面主波和橋波（前、后橋波）。前端主波和前橋波的強度高于后端主波，側(cè)面主波靠近前端的沖擊波強度強于靠近后端的沖擊波強度。引爆后1ms（見圖12（c）），沖擊波場呈橢球狀，即裝藥側(cè)面主波的范圍更大。端部主波、側(cè)面主波和橋波逐漸融合，沖擊波的波陣面更光滑，而裝藥前端的沖擊波強度仍高于后端。引爆后10ms（見圖12（d）），沖擊波場的形貌趨于球形，波陣面演化為球面，裝藥前端的強度仍略高于后端。

2.3參數(shù)影響分析

2.3.1裝藥長徑比

圖13顯示了靜爆工況下不同長徑比柱形裝藥的爆炸沖擊波壓力云圖?？梢钥闯?，隨著長徑比的增大，側(cè)面主波范圍變大，橋波和端部主波范圍減小，沖擊波波陣面整體趨于平坦，與爆轟產(chǎn)物分布一致，沖擊波場的強度分布變化較小，側(cè)面主波的峰值強度逐漸向前移動。

圖14顯示了靜爆工況下比例距離為0.4、1.0和5.0m/kg1/3時不同長徑比柱形裝藥的沖擊波峰值超壓和最大沖量?？梢钥闯觯诒ń鼌^(qū)（Z=0.4m/kg1/3）、0°～180°范圍內(nèi)的峰值超壓呈W形分布，即軸向和徑向的峰值超壓高于二者之間過渡區(qū)的峰值超壓，其中0°方向（起爆點對側(cè)）的峰值超壓高于180°方向（起爆點同側(cè)）的峰值超壓；而最大沖量則呈A形分布，即徑向的最大沖量高于軸向。在爆炸中區(qū)（Z=1.0m/kg1/3），峰值超壓和最大沖量均呈A形分布。在遠(yuǎn)區(qū)（Z=5.0m/kg1/3），峰值超壓和最大沖量在各方位角的分布基本相同，這是由于高強度的沖擊波波陣面比低強度的波陣面衰減快，超過一定距離后，各方位角的沖擊波強度趨于一致。長徑比（6～10）對峰值超壓和最大沖量的分布基本無影響，這是因為長徑比的變化較小，對較大空間范圍內(nèi)的超壓和沖量影響有限。

2.3.2裝藥運動速度

定義沖擊波場幾何中心為側(cè)面主波外凸點在裝藥軸線的投影（圖12（b））。圖15顯示了爆炸后0.1ms不同運動速度下柱形裝藥（L/D=6）的沖擊波壓力云圖?？梢钥闯?，隨著運動速度增加，沖擊波場整體前移，側(cè)面主波的前端和前橋波逐漸平滑且強度不斷增強，而后端主波、后橋波和側(cè)面主波后端的強度降低。爆轟產(chǎn)物在端部和側(cè)面呈現(xiàn)外凸?fàn)睿谶^渡區(qū)呈內(nèi)凹狀，其分布與沖擊波外輪廓基本一致，在側(cè)面前側(cè)和過渡區(qū)交界處，爆轟產(chǎn)物隨運動速度的增加呈現(xiàn)出向前的運動趨勢（溢出）。

沖擊波場的運動包含爆炸沖擊波壓力場的整體移動和波陣面的傳播，二者均受裝藥運動速度的影響。沖擊波場的移動可近似為其幾何中心的移動，以L/D=6為例，圖16顯示了Ma為0～4時沖擊波場幾何中心的位移-時間曲線?？梢钥闯?，受慣性影響，爆轟產(chǎn)物和壓縮空氣形成的沖擊波場會保持一定的速度繼續(xù)向前運動直至速度衰減為零，移動距離與裝藥初速度正相關(guān)。需要指出的是，受裝藥形狀和尾端起爆方式的影響，靜爆工況下幾何中心有較小的位移。

圖17顯示了比例距離為0.4、1.0和5.0m/kg1/3（分別對應(yīng)爆炸近區(qū)、中區(qū)和遠(yuǎn)區(qū)）時不同運動速度下裝藥的爆炸沖擊波峰值超壓和最大沖量?？梢钥闯?，存在臨界角度αc，方位角為0°～αc時峰值超壓隨運動速度的增大而增大，方位角為αc～180°時變化趨勢相反。αc隨著比例距離的增大而增大。由爆轟理論可知，爆炸沖擊波初始壓力和波速與空氣質(zhì)點速度相關(guān)：當(dāng)波陣面?zhèn)鞑シ较蚺c運動方向的夾角小于αc時，質(zhì)點速度與波陣面?zhèn)鞑シ较虺射J角，沖擊波壓力增大，且沖擊波場強度的衰減變緩；當(dāng)波陣面?zhèn)鞑シ较蚺c運動方向夾角大于αc時，質(zhì)點速度與波陣面?zhèn)鞑ニ俣确较虺赦g角，沖擊波場強度的衰減加快，沖擊波壓力減??；當(dāng)波陣面?zhèn)鞑シ较蚺c運動方向夾角等于αc（對應(yīng)側(cè)面主波的外凸點）時，對沖擊波壓力基本無影響。裝藥運動速度對最大沖量的影響與峰值超壓相似。

3運動爆炸沖擊波荷載計算模型

對典型戰(zhàn)斗部的爆炸沖擊波荷載開展數(shù)值仿真分析，提出入射和反射沖擊波載荷的計算模型。

3.1入射沖擊波荷載分布特征

圖18顯示了L/D=6時柱形裝藥在不同運動速度下的入射峰值超壓。由于L/D遠(yuǎn)大于1，柱形裝藥側(cè)表面積遠(yuǎn)大于端面，能量釋放位置集中在徑向，因此徑向的峰值超壓明顯高于軸向和過渡區(qū)的峰值超壓。在尾端起爆（α=180°）工況下，起爆端軸向的峰值超壓小于對側(cè)的峰值超壓。由Chapman-Jouguet爆轟理論[29]可知，爆轟波以起爆點為中心向外傳播，傳過的炸藥受強烈沖擊生成爆轟產(chǎn)物并釋放大量能量，起爆點對側(cè)的裝藥較多，釋放的能量也更大，因此其峰值超壓大于起爆端。隨著裝藥運動速度的增大，裝藥前端的峰值超壓增大，而后端的峰值超壓減小。這是因為裝藥運動提高了前端沖擊波波陣面的速度，由Rankine-Hugoniot方程[16]可知，波陣面上的壓力也增大，而尾端沖擊波的傳播方向與裝藥運動方向相反，沖擊波波陣面速度減小，導(dǎo)致壓力降低。

圖19顯示了L/D=6時柱形裝藥在典型方位角處的峰值超壓。峰值超壓大致與比例距離負(fù)相關(guān)，但α=0°時Ma為0和1的峰值超壓在Z=2m/kg1/3附近出現(xiàn)小幅躍升后繼續(xù)下降，α=180°時在比例距離2～3m/kg1/3處也出現(xiàn)了類似現(xiàn)象。這是因為側(cè)面主波的一部分傳播至裝藥軸向時在端部主波的后側(cè)形成第2個波陣面（見圖15（a）），第2波陣面的速度較高，追趕前一波陣面并與之合并，使波陣面壓力進一步升高。

圖20顯示了L/D=6時柱形裝藥在不同運動速度下的入射沖擊波最大沖量?？梢钥闯?，與峰值超壓類似，徑向的最大沖量高于軸向和過渡區(qū)，起爆端軸向的最大沖量高于對側(cè)，隨著裝藥運動速度的增大，運動方向的最大沖量增大，而反方向的最大沖量減小。與峰值超壓不同的是，在特定范圍（如α=0°時的0.4～1m/kg1/3）內(nèi)最大沖量隨比例距離的增大而增大，表現(xiàn)為玫瑰多葉曲線中的“穿透”現(xiàn)象，這一現(xiàn)象在軸向尤其明顯，如圖21所示。這是因為：（1）受圖15（a）中的二次沖擊波影響，最大沖量是兩個波陣面的沖量之和；（2）隨著比例距離的增大，超壓正相持時并非線性變化[1]，造成超壓關(guān)于持時的積分并非單調(diào)變化，這一現(xiàn)象在球形裝藥靜爆工況中同樣存在[1]。

3.2入射沖擊波荷載的計算模型

基于2.3節(jié)的結(jié)果，進一步考慮裝藥運動速度對沖擊波荷載的影響，修正文獻[15]中考慮裝藥長徑比、比例距離和方位角的柱形裝藥靜爆沖擊波荷載計算公式，柱形裝藥空中自由場動爆的入射峰值超壓pso可以表示為：

式中：c=340m/s為聲波在標(biāo)準(zhǔn)大氣壓中的傳播速度，Aj，1、Aj，2、Aj，3、Aj，4、Aj，5、Aj，k，1、Aj，k，2和Aj，k，3為擬合系數(shù)。基于Levenberg-Marquardt算法擬合的系數(shù)列于表2，其適用范圍為：6≤L/D≤10，0m/s≤v0≤1360m/s，0.4m/kg1/3≤Z≤15m/kg1/3，擬合優(yōu)度R2=0.975。

式中：Bj，1、Bj，2、Bj，3、Bj，4、Bj，5、Bj，k，1、Bj，k，2和Bj，k，3為擬合系數(shù)。需要指出的是，Z＜1m/kg1/3時，Iso的變化比較復(fù)雜（見圖20），為此，僅討論Z≥1m/kg1/3時的Iso。比沖量的擬合系數(shù)如表3所示，其適用的范圍為：6≤L/D≤10，0m/s≤v0≤1360m/s，1m/kg1/3≤Z≤15m/kg1/3，擬合優(yōu)度R2=0.945。

3.3反射沖擊波荷載的計算模型

基于“三階段”有限元分析方法建立100組爆炸工況的有限元模型（包括75組反射場模型和25組自由場模型），計算超壓反射系數(shù)μ（同一位置處反射和入射沖擊波峰值超壓之比），結(jié)合3.2節(jié)的入射沖擊波荷載計算模型提出反射沖擊波荷載的計算模型。

3.3.1反射場有限元模型

典型工況下戰(zhàn)斗部彈道垂直于目標(biāo)結(jié)構(gòu)表面，基于“三階段”有限元分析方法對反射場運動爆炸工況進行數(shù)值模擬。在垂直比例距離（Zv，即裝藥中心到結(jié)構(gòu)表面的垂直距離Rv與裝藥當(dāng)量的1/3次方之比，Zv=Rv/W1/3）為0.5、1.5和2.5m/kg1/3處設(shè)置剛性墻并布置測點，以測量不同入射角和比例距離處的壓力，其中入射角θ為剛性墻上任意一點到裝藥中心的連線與裝藥軸向（α=0°）的夾角，如圖22（a）所示。柱形裝藥質(zhì)量為1kg，起爆方式為尾端起爆，長徑比為6～10，Ma為0～4，剛性墻高度為15m，反射場有限元模型（以Zv=2.5m/kg1/3為例）如圖22（a）所示。為獲取剛性墻上的超壓反射系數(shù)，建立相應(yīng)長徑比和運動速度下的自由場有限元模型（見圖22（b）），并在空氣域相同位置處布置測點以測量入射沖擊波的壓力。

3.3.2超壓反射系數(shù)

以L/D=6的柱形裝藥為例，圖23顯示了不同Zv處的超壓反射系數(shù)μ。可以看出：當(dāng)Zv=0.5m/kg1/3、Ma=0時，μ隨著θ的增大呈鞍形曲線遞減，存在轉(zhuǎn)折入射角θc1和θc2；當(dāng)θ＜θc1時，μ隨θ的增大而減??；當(dāng)θc1≤θ≤θc2時，μ增大；當(dāng)θ＞θc2時，μ繼續(xù)遞減至1；θc1和θc2與v0和Zv均有關(guān)系。這是因為：（1）剛性墻上的反射沖擊波在向前傳播過程中與入射沖擊波疊加形成馬赫波，強度增大，反射系數(shù)也增大；（2）柱形裝藥在徑向釋放的能量較多，反射超壓和反射系數(shù)均增大。隨著v0的增大，θ較小時，μ與v0正相關(guān)，θ較大時μ與v0負(fù)相關(guān)。而隨著Zv的增大，μ減小，v0和θ引起的μ差異也變小，這是因為，隨著Zv的增大，入射峰值超壓差異減小，導(dǎo)致μ的差異減小。

3.3.3荷載計算模型

結(jié)合3.2節(jié)中的入射沖擊波峰值超壓pso（式（3））和3.3.2節(jié)中的超壓反射系數(shù)μ，計算尾部起爆工況下柱形裝藥作用在剛性墻（垂直于裝藥軸向）上的反射沖擊波峰值超壓pr：

式中：μ通過數(shù)值模擬插值獲得。式（6）的適用范圍為0.5m/kg1/3≤Zv≤2.5m/kg1/3。

基于文獻[15]中沖擊波超壓時程曲線為突加三角形以及同一位置處入射和反射沖擊波正相持時相等的假設(shè)，計算反射沖擊波的最大沖量Ir：

式中：td為入射和反射沖擊波超壓的正相持時。式（7）～（8）的適用范圍為：6≤L/D≤10，0m/s≤v0≤1360m/s，1.0m/kg1/3≤Zv≤2.5m/kg1/3。

4計算模型驗證

為檢驗3.2和3.3節(jié)中尾部起爆工況下柱形裝藥的動爆入射和反射沖擊波載荷計算模型的適用性，基于有限元分析方法，開展23和280kg兩種戰(zhàn)斗部的TNT柱形裝藥數(shù)值模擬，兩種裝藥的長徑比分別為8.0和6.8，裝藥運動的Ma為0～3。每種戰(zhàn)斗部開展10組沖擊波入射場和30組沖擊波反射場工況的數(shù)值模擬。

圖24顯示了23kg裝藥在不同爆炸工況下的入射峰值超壓pso、入射最大沖量Iso、反射峰值超壓pr和反射最大沖量Ir，其中V0-0°代表裝藥的馬赫數(shù)為0、方位角為0°的工況，下標(biāo)cal代表計算結(jié)果，下標(biāo)sim代表模擬結(jié)果。可以看出，數(shù)值模擬與計算結(jié)果吻合較好：對于pso和Iso，模擬和計算結(jié)果均較為接近，且其吻合程度隨裝藥運動速度的增大而降低；對于pr和Ir，75%以上數(shù)據(jù)點的相對誤差小于±25%。圖25顯示了280kg裝藥的動爆沖擊波荷載，可以看出，數(shù)值模擬和計算結(jié)果吻合較好。對于pso和Iso，α取45°和90°時的吻合程度優(yōu)于α=0°時。對于pr和Ir，超過61%pr數(shù)據(jù)點的相對誤差在±25%以內(nèi)，100%的Ir數(shù)據(jù)點位于+10%～?25%的相對誤差范圍內(nèi)。此外，彈藥的運動速度越小，計算模型預(yù)測的pr的準(zhǔn)確性越高。

數(shù)值仿真與計算模型結(jié)果吻合較好，說明提出的計算模型適用于計算柱形裝藥空中運動爆炸沖擊波荷載。

5結(jié)論

基于AUTODYN軟件，提出了考慮運動擾動壓力場的“三階段”柱形裝藥運動爆炸有限元分析方法，通過與已有靜爆和運爆試驗結(jié)果的對比，驗證了該方法的可靠性；基于驗證的有限元分析方法，對柱形裝藥空中靜、動爆工況開展數(shù)值仿真計算，定量分析了裝藥運動擾動壓力場對爆炸沖擊波場的影響，主要結(jié)論如下：

（1）裝藥運動引起的擾動壓力場與爆炸沖擊波場耦合會改變?nèi)肷錄_擊波荷載的分布，對運動方向的荷載影響最大，并且影響隨比例距離的減小和運動速度的增大而增強。

（2）入射沖擊波荷載在軸/徑向較大，過渡區(qū)較小，且隨比例距離的增大，波陣面強度逐漸趨近；受運動爆炸入射沖擊波影響，以臨界角為界的沖擊波前端波陣面強度增強而后端減弱，沖擊波場中心隨運動速度的增大逐漸前移。

（3）提出了考慮裝藥運動速度、長徑比、方位角和比例距離的柱形裝藥動爆入射和反射沖擊波荷載計算模型，并通過兩種典型戰(zhàn)斗部裝藥動爆工況的數(shù)值仿真驗證了計算模型的可靠性。

需要指出的是，本文中，僅分析了柱形裸藥的爆炸沖擊波荷載，未考慮彈藥殼體的影響，后續(xù)可針對殼體及其破片與沖擊波荷載的耦合作用開展進一步研究。