JH-2模型參數確定及花崗巖重復侵徹數值分析

王志亮，李允忠，黃佑鵬

(合肥工業(yè)大學 土木與水利工程學院,合肥 230009)

彈體對地質材料的侵徹屬于沖擊動力學研究的一個主要內容,其實驗和理論研究已有三百多年的歷史.因其在軍事方面有著廣泛的應用，目前仍是一個非常熱門的研究領域.當今，隨著侵徹深度和爆炸規(guī)模大的精確制導武器的迅速發(fā)展，無論是防護結構計算理論與防護技術，還是制導武器戰(zhàn)斗部研制均亟需解決巖石中侵徹局部破壞效應問題.

隨著現(xiàn)代科學技術的飛速發(fā)展，現(xiàn)代常規(guī)武器的性能日臻完善.對于深埋地下的重要防護工程，實施單次打擊并不能將其一次性摧毀，故考慮對同一目標實施多彈重復打擊.而且，信息與制導技術的廣泛應用使常規(guī)武器的精確多彈重復打擊成為可能.迄今，國內外學者對單次打擊巖石及混凝土的相關問題開展大量研究.楊剛等[1]基于三維FEM-SPH自適應耦合算法模擬了子彈侵徹混凝土靶跳飛問題，指出作用于子彈的不對稱阻力所產生的力矩是子彈跳飛的根本原因；Wan等[2]通過鋼管約束混凝土(STCC)目標的典型損傷模式，分析了沖擊位置、約束鋼管直徑和厚度對STCC靶材抗穿透性能的影響；劉崢等[3]基于AUTODYN開展了鎢桿彈對成層式防護結構的超高速打擊數值計算分析，得到了結構的破壞特征和能量分配情況；王海兵等[4-6]通過彈丸侵徹花崗巖和石灰?guī)r靶板實驗，分析了靶板裂紋長度、侵徹深度及靶板的破壞類型；宋春明等[7-8]對超高速彈體侵徹巖石靶板進行研究，認為彈體對巖石靶板的侵徹機制隨撞擊速度的增加發(fā)生“剛體侵徹→半流體侵徹→流體侵徹”的轉變；王明洋等[9-10]通過分析侵爆近區(qū)巖石介質的動態(tài)可壓縮性行為，創(chuàng)造性地構建了流體彈塑性內摩擦侵徹理論模型，并基于界定超高速侵徹速度內涵范圍及鉆地彈固體侵徹、擬流體侵徹和流體侵徹的最小動能閾值，系統(tǒng)地提出了超高速動能彈打擊侵深、成坑及地沖安全厚度的計算方法.Gomez等[11-13]對多彈重復打擊破壞效應實驗及數值模擬分別進行了一些研究，但多以侵爆為主.考慮到侵爆耦合及各個過程的復雜性，用控制變量法分析單一因素對重復打擊破壞效應的影響則較為困難.目前，國內外在多彈重復打擊破壞效應方面的公開研究成果很少，可能歸結于這類實驗研究危險較大、成本高，且往往只能獲得被打擊目標破壞后的最終響應結果及部分測點信息.

基于上述研究現(xiàn)狀，本文將多彈重復侵徹過程及其動力響應作為研究重點.根據已有實驗數據和理論推導，初步確定出花崗巖的Johnson-Homquist(JH-2)材料模型與拉伸斷裂軟化模型耦合模型的相關參數.通過模擬花崗巖單次侵徹試驗驗證并優(yōu)化已確定的參數.最后，將優(yōu)化后的模型參數用于描述花崗巖的動力響應，運用非線性動力分析軟件AUTODYN[14]并基于FEM-SPH耦合法對花崗巖多彈重復侵徹過程進行模擬，力求得出對防護工程設計與建設有參考價值的結論.

1 花崗巖損傷本構模型

非線性動力分析軟件AUTODYN是由美國世紀動力公司研發(fā)的，已成為國際上爆炸力學、高速沖擊碰撞領域研究最著名數值模擬軟件之一[14]，為模擬高速碰撞問題提供了豐富的材料模型，如巖石、混凝土、金屬和陶瓷等.材料模型通常由兩部分組成[15]：1)強度模型，可利用強度模型通過塑性理論確定材料的畸變響應；2)狀態(tài)方程 (EOS)，可通過EOS描述材料的體積響應.本文引入JH系列模型中JH-2模型[16]，并將其與拉伸斷裂軟化模型相耦合[18]，用以模擬花崗巖靶板內高應力區(qū)塑性壓縮和剪切破壞效應以及低應力區(qū)靶板在主拉伸應力作用下產生的損傷和裂紋擴展.其耦合機理詳見文獻[4,17]，此處不再贅述.

1.1 模型介紹

1.1.1 強度模型

JH-2模型的應力和損傷是壓力和其他變量的分析函數，允許以更系統(tǒng)的方式對常數進行參數變化；當損傷開始累積時，JH-2模型材料開始軟化，允許材料在塑性應變增加的情況下逐漸軟化(JH-1模型在材料完全破壞前不會軟化，而在完全破壞時軟化立即發(fā)生).JH-2材料強度模型描述如圖1(a)所示.可以看出，JH-2材料模型強度是完整強度、斷裂強度、應變速率和損傷的平滑變化函數.強度的歸一化等效應力可表述為[16]

(1)

完整材料和完全破壞材料的歸一化等效應力表述為[16]

(2)

(3)

損傷因子D是巖石材料從完整到破裂強度的過渡，損傷演化表達式為[16]

(4)

(5)

1.1.2 狀態(tài)方程

為了表征巖石類材料在高壓沖擊下的非線性變形特征，材料發(fā)生損傷前的靜水壓力與體積應變關系用三次多項式狀態(tài)方程描述為[16]

P=K1μ+K2μ2+K3μ3.

(6)

式中：K1為體積模量，K2、K3為壓力常數；體積應變μ=ρ/ρ0-1，ρ0為初始密度；當μ<0，巖石受拉伸應力作用，此時P=K1μ.當材料損傷開始累積(D>0)時，體積膨脹，靜水壓力附加一個壓力增量ΔP，即[16]

P=K1μ+K2μ2+K3μ3+ΔP.

(7)

式中：0≤ΔP≤ΔPMax，ΔP在D=0和D=1時分別取到兩端點值，如圖1(c)所示.

圖1 JH-2材料模型

1.2 參數確定

1.2.1 狀態(tài)方程參數確定

花崗巖基本物理力學參數[18]：初始密度ρ0=2 657 kg/m3，泊松比ν=0.29，彈性模量E=80 GPa.則剪切模量G=E/2(1+ν)=30.01 GPa，體積模量K1=E/3(1-2ν)=63.49 GPa.確定巖石P-μ關系所需的高壓數據可通過平板試驗獲得，如圖2所示.利用式(7)擬合花崗巖不同體積應變率下靜水壓力[19-20]及花崗巖Hugoniot實驗數據[21]可得K1=63.44 GPa，K2=448 GPa，K3=-2 556 GPa.

圖2 P-μ關系

1.2.2 強度參數確定

Hugoniot彈性極限是用于確定材料完整強度的重要參數，包含了壓力和偏應力分量，其定義為

(8)

(9)

1)D=0時，完整巖石材料強度在不考慮應變率影響的情況下，由三軸壓縮實驗結果[19]經式(2)擬合得常數A=1.075，N=0.721，如圖3所示.由擬合方程解得曲線與橫坐標交點P*=-T*=-0.016 7，則最大靜拉伸應力T=-T*·PHEL=-49.43 MPa.

圖3 完整材料等效應力

圖4 應變率系數及其測定[24]

1.2.3 損傷參數確定

巖石材料從未損傷強度到完全破壞強度的過渡由損傷因子D描述.完全破壞脆性材料的塑性應變量非常小且不可能直接測量，導致無法求解損傷參數D1、D2.相反，可通過數值模擬試算調整以達到可接受的破壞效果獲得可靠的D1、D2的值.模擬中暫取初始損傷參數D1=0.005，D2=0.7[4].

1.2.4 拉伸斷裂軟化模型參數確定

在拉伸斷裂軟化模型中，主拉伸失效應力Tf為35 MPa[4]；而在AUTODYN中拉伸斷裂軟化模型只涉及一個參數，即斷裂能Gf，默認值Gf=70 J/m2[4].

2 損傷常數優(yōu)化及模型參數驗證

2.1 SPH算法原理

在爆炸與沖擊問題的數值模擬中，材料的高應變率、大變形及失效破壞均要求數值算法能合理描述邊界畸變、裂紋擴展等，同時必然對網格精度及計算成本提出較高要求.對于Lagrange算法，大變形的網格畸變會使得時間步長過小而導致計算困難甚至終止；對于Euler算法，其很難精確捕捉固體材料的變形響應；而ALE(Arbitrary Lagrange-Euler)算法也存在著因網格不協(xié)調導致的Euler物質非物理穿透問題[25].為克服上述各種缺陷，一種無網格算法——光滑粒子流體動力學方法(SPH)被提出[26]，用于模擬天體物理、連續(xù)體結構的解體、碎裂及固體的層裂、脆性斷裂等問題，并逐步解決了可變光滑長度、人工黏性、零能模式、邊界處理與多物質界面等相關技術問題.

SPH法通過核函數W(x-x′,h)逼近支持域Ω內任意點的場變量函數f(x)，實現(xiàn)連續(xù)體離散成一系列具有質量、坐標、速度與內能等的粒子[25].任意場變量函數f(x)近似值表達式為

(10)

式中：h為光滑長度，x為質點位置.

在由Monaghan[27]提出的B-Spline光滑函數或插值核心函數鄰域內，如圖5所示，粒子i的某個變量可以通過粒子集j的加權求和得到

(11)

式中：Ψj=mj/ρj為與粒子j相關的體積，mj、ρj分別為粒子集j的質量與密度(j=1,2,3,…,L).

圖5 二維粒子鄰域與核函數

2.2 數值模擬驗證分析

為得到完整可靠的JH-2與拉伸斷裂相耦合模

型參數，采用上述初步確定的狀態(tài)方程和部分強度模型參數，利用AUTODYN模擬卵形彈丸沖擊花崗巖靶板實驗[29].首先，將數值計算的損傷效果與實驗結果進行對比；然后，通過數值調整法分析花崗巖在彈丸沖擊后的損傷狀態(tài)，從而得到可接受的損傷效果及損傷參數D1、D2.為了增加可比性，二維數值模型中卵形彈丸和靶板尺寸與實驗相同，如圖6所示.

卵形彈丸射彈頭型系數CRH=3.0，材料為4340鋼，采用Johnson-Cook材料模型描述，相關參數均采用AUTODYN材料庫中參數；花崗巖靶板材料采用JH-2模型與拉伸斷裂軟化模型描述，并在靶板兩端施加固定約束.卵形彈丸采用 Lagrange算法，花崗巖靶板采用SPH算法，通過AUTODYN中的FEM-SPH耦合法實現(xiàn)靶板沖擊過程.

圖6 靶板與彈丸尺寸示意

一個完整沖擊過程中彈丸速度與加速度變化規(guī)律如圖7所示.彈丸初始彈速為279 m/s，在彈丸侵徹開坑階段，侵徹深度在2倍彈徑范圍內(0

數值模擬所得花崗巖靶板在卵形彈丸沖擊下的最終損傷形態(tài)與實驗結果對比如圖8所示，可見彈丸沖擊初始階段，在迎沖面著彈點附近花崗巖的損傷由剪切應力引起的剪切損傷為主導；隨著沖擊波向前傳播，損傷的產生機理為主拉應力，裂紋通過拉伸應變萌生并向外擴展，最終在花崗巖體中形成徑向拉伸裂紋.將由數值模擬所得的靶板迎沖面彈坑張角α和直徑Df及迎沖背面錐形沖擊塞高度H和最大直徑Dr的值與相應的實驗值列于表1.可以看出，相關變量的模擬值與實驗值相比，其誤差絕對值均在10%以內.同時，通過圖8中花崗巖靶板裂紋區(qū)分布情況對比可知，模擬所得結果與實驗結果基本一致，說明本文模擬方法能夠較好地重現(xiàn)實驗中花崗巖真實破壞狀態(tài)，亦證明了采用JH-2模型與拉伸斷裂軟化模型相結合模擬花崗巖靶板侵徹的合理性及模型參數的可靠性.基于上述分析，JH-2材料模型與拉伸斷裂軟化模型參數列于表2.

圖7 彈丸速度與減速度曲線

圖8 花崗巖靶板沖擊后裂紋分布情況[29]

表1 損傷靶板相關變量的模擬值與實驗值誤差分析

表2 JH-2模型與拉伸斷裂軟化模型參數

3 彈丸重復沖擊數值分析

為得出花崗巖在彈丸重復沖擊下的相關規(guī)律，用彈徑為30 mm，CHR=3.0，長徑比為4.7的彈體沖擊直徑為1 500 mm，高為800 mm的花崗巖靶板，待前次彈體速度出現(xiàn)負值，或彈體有回彈趨勢，則停止計算，并重新賦予彈體相同彈速，實現(xiàn)重復侵徹.彈體材料用4340鋼，參數與2.2節(jié)所用材料相同；花崗巖用JH-2模型與拉伸斷裂軟化模型耦合模型，并采用表1中參數.考慮模型尺寸較大，為節(jié)約計算時間并提高計算效率，靶板徑向0～150 mm范圍采用SPH光滑粒子算法，150～750 mm范圍用FEM算法，并將SPH粒子與FEM單元用共節(jié)點連接，如圖9所示，與節(jié)點到中面節(jié)點以及任意節(jié)點到面節(jié)點兩種耦合方式相比[26]，F(xiàn)EM-SPH采用共節(jié)點耦合，其界面處相關變量的連續(xù)性更好.為消除反射波影響，靶板模型左側面及底面設為無反射邊界(即透射邊界)，用以模擬半無限靶板侵徹.

圖9 數值模型及FEM-SPH耦合方式

以不同彈速Vs(220,270,320,370及420 m/s)重復侵徹花崗巖靶板，每個彈速重復侵徹7次.以彈速Vs=220 m/s為例，數值模擬所得部分重復侵徹過程如圖10所示.彈體第一次沖擊花崗巖靶板的破壞情況如圖10(a)所示，迎沖面整體呈“錐形開坑區(qū)+隧道區(qū)”組合的破壞形態(tài).介于初次形成的彈道及其對后續(xù)侵徹彈體的橫向約束作用，后續(xù)各侵徹過程相對穩(wěn)定，即隨沖擊次數的增加，雖然錐形彈坑直徑有一定增加，但變化范圍不大；由圖10(b)～10(d)可知，與初次侵徹相比，后續(xù)侵徹各錐形彈坑與彈道交界點，即彈坑底部到自由面距離的增幅基本一致，對迎沖面的影響不大.同時，由彈速Vs=220 m/s重復侵徹7次條件下連續(xù)侵徹深度時程曲線可知，如圖11所示，初次侵徹為關鍵的開坑階段，其侵徹深度及歷時最大；隨侵徹次數n的增加，相對侵徹深度及相應的歷時逐漸減小，并漸趨于穩(wěn)定值，亦說明后續(xù)侵徹過程處在一個相對穩(wěn)定侵徹階段.

圖10 彈速Vs=220 m/s重復侵徹過程

圖11 彈速Vs=220 m/s重復侵徹深度時程曲線

圖12為相對侵徹深度(relative depth of penetration，R-DOP)隨沖擊次數變化的規(guī)律曲線圖.可以看出，不同彈速下相對侵徹深度隨沖擊次數的增加逐漸減小，且n=1時的相對侵徹深度與后續(xù)相比相差較大，從n=2時起，這種變化關系較為平緩，彈體呈現(xiàn)為較穩(wěn)定的侵徹.整體上相同沖擊次數下(除n=2外)彈速越高相對侵徹深度越大，僅當n=2時，相對侵徹深度變化不大，穩(wěn)定在80 mm左右.相對侵徹深度隨沖擊次數逐漸減小的原因在于所形成的隧道區(qū)對后續(xù)各彈體的橫向約束作用；加之隧道區(qū)內壁密度增大及彈體與隧道區(qū)內壁存在摩擦作用使得彈體動能耗散，且沖擊次數越大，彈體所經歷的隧道區(qū)越長，以上作用越大.

圖12 相對侵徹深度與沖擊次數關系

從侵徹深度與沖擊次數變化規(guī)律上看，不同彈速下侵徹深度整體上隨沖擊次數的增加逐漸增大，如圖13所示.彈速較低時，如220 m/s≤Vs≤320 m/s，在重復侵徹次數n<4時，侵徹深度變化率變化較大，表現(xiàn)出較為明顯的非線性關系；而在重復侵徹次數n≥4時，侵徹深度變化率變化不大，線性關系明顯；彈速較高時，如370 m/s≤Vs≤420 m/s，侵徹深度與沖擊次數整體呈較強的線性相關關系.由上述可知，在所討論速度范圍內，當彈體重復沖擊次數越大或重復沖擊的速度越高時，侵徹深度與沖擊次數線性關系越明顯，相對侵徹深度值越穩(wěn)定，即彈體侵徹花崗巖靶板的過程越穩(wěn)定.

圖13 侵徹深度與沖擊次數關系

Forrestal等[30]給出經典的卵形彈丸單次沖擊靶板的侵徹深度公式：

(12)

(13)

式中：m′、a為射彈的質量(kg)和半徑(m)；ρ為靶板密度(kg/m3)；N為CRH(曲徑比Ψ)的函數，N*=(8Ψ-1)/(24Ψ2)；S為材料強度修正因子，S=82.6(fc×10-6)-0.544，fc為靶板材料單軸抗壓強度，Pa；z為單次侵徹中錐形彈坑底部與所形成的彈道過渡點深度(m)，取z=4a.

Gomez等[11]認為：彈體以同等彈速第2次射擊靶板的同一目標點，此時靶板初始條件與第1次沖擊相比已發(fā)生變化，其通過重復侵徹實驗反演出不同彈速相同沖擊次數下材料強度修正因子S給出重復侵徹深度計算方法.花崗巖靶板無側限抗壓強度fc被沖擊后在整個響應區(qū)內并非均勻性變化，如在迎沖面臨近空腔[7]的粉碎區(qū)內花崗巖無側限抗壓強度變化最大，幾乎已無抗沖擊能力，而在彈性區(qū)及未擾動區(qū)內無側限抗壓強度基本沒變.由此可見，Gomez通過反演S將材料強度整體修正的做法有不足之處.本文以經典的Forrestal公式作為不同彈速下重復侵徹的初始侵徹深度值，記為P0，并通過定義相對侵徹系數函數φ(n)給出重復侵徹深度的表達式：

P(n,Vs)=φ(n)P0,n=1,2,3,…,7,

(14)

(15)

式中：DOPn為第n次侵徹深度，取n=1時其為第1次侵徹深度.

式(14)顯示P(n,Vs)僅在因式P0項中涉及因子S，且為材料初始因子S0，避免了將靶板材料損傷后因子S均質化問題.考慮到花崗巖材料的非均勻性使得初始因子S0存在一定差異，由不同彈速下初始侵徹深度，通過式(12)、(13)解得均值S0=3.226，使P0僅為速度Vs函數，記為P0(Vs),則式(14)變?yōu)?/p>

P(n,Vs)=φ(n)P0(Vs).

(16)

將不同彈速下隨沖擊次數變化的靶板侵徹深度(DOPn)按式(15)進行歸一化處理，如圖14所示.相對侵徹系數φ(n)與侵徹次數的關系采用如下表達式

φ(n)=A*+B*·ek*·n.

(17)

式中：A*、B*和k*為相關常數.通過擬合數據得A*=4.696，B*=-4.192，k*=-0.126.擬合公式的可決系數R2=0.986，說明擬合效果很好.

圖14 相對侵徹系數

由上表述，聯(lián)立式(12)、(13),(16)、(17)得關于速度Vs與沖擊次數n的重復侵徹深度表達式:

P(n,Vs)=(4.696-4.192e-0.126n)P0(Vs),
n=1,2,3,…,7.

(18)

式中：各參數含義與前述相同；P0(Vs)由式(12)、(13)計算可得，并在式(13)中取S=S0.

將按式(18)計算值作為理論計算值，以曲面形式繪于圖15中，并與數值模擬獲取的數據(圖中實心圓點)形成對比.為合理顯示理論計算值與數值模擬結果之間偏離度的大小，其偏離度用模擬數據點到理論計算曲面的垂線段表示，垂線段越長，則模擬數據點到目標值的偏離度越大；反之，偏離度越小，以此來驗證式(18)的合理性.由圖15可知,該式能較準確地預測花崗巖重復侵徹深度.

圖15 侵徹深度隨沖擊次數和彈速變化關系

4 結 論

基于FEM-SPH耦合算法，在220 m/s≤Vs≤420 m/s彈速下，對花崗巖多彈重復侵徹問題進行了數值研究，得如下主要結論：

1)拉伸斷裂軟化模型與JH-2模型相耦合，可實現(xiàn)JH-2模型模擬靶板內低應力區(qū)在主拉伸應力作用下徑向裂紋的產生與發(fā)展，且確定的材料模型參數能較準確模擬巖石靶板破壞形態(tài).

2)多彈重復侵徹巖石靶板，靶板自由邊界對初始侵徹破壞形態(tài)的影響較大，在迎沖面呈“錐形開坑區(qū)+隧道區(qū)”破壞形態(tài)，隧道區(qū)對彈體的橫向約束使后續(xù)各侵徹過程處在一個相對穩(wěn)定侵徹階段.

3)不同彈速下侵徹深度隨沖擊次數的增大而增大，而相對侵徹深度則減小；當彈體重復沖擊次數越多或侵徹速度越高，彈體侵徹花崗巖靶板越穩(wěn)定.

4)基于所提出的相對侵徹系數函數給出重復侵徹深度計算式P(n,Vs)，能較準確描述侵徹深度隨沖擊次數與沖擊彈速變化的規(guī)律.

彈體侵徹花崗巖是個復雜的過程，本文探討了剛性侵徹條件下重復侵徹深度的變化規(guī)律，對于半流體侵徹及流體侵徹條件下的超高速重復侵徹深度的變化規(guī)律尚需進一步研究.