劉俊新，楊春和，劉育田

(1.西南科技大學(xué) 土木工程與建筑學(xué)院,四川 綿陽 621010；2. 中國(guó)科學(xué)院 武漢巖土力學(xué)研究所 巖土力學(xué)與工程國(guó)家重點(diǎn)實(shí)驗(yàn)室,武漢 430071)

鉆地武器指攜帶鉆地彈頭鉆入地下目標(biāo)引爆的精確制導(dǎo)武器。該種武器鉆入地下目標(biāo)過程及目標(biāo)內(nèi)部后的延時(shí)爆炸均涉及對(duì)目標(biāo)組成部分在沖擊荷載作用下動(dòng)力響應(yīng)。目標(biāo)主要由巖土介質(zhì)與混凝土組成。對(duì)壓實(shí)黏土受沖擊荷載作用的動(dòng)力響應(yīng)在防護(hù)工程及武器設(shè)計(jì)中有重要意義。

在沖擊荷載作用下土體處于復(fù)雜應(yīng)力狀態(tài)，土體自身特征導(dǎo)致其與金屬材料不同特性[1-3]，如：① 材料剪切破壞前會(huì)萌生大量微孔隙，引起體積變形增大而產(chǎn)生剪脹特性；② 材料應(yīng)力-應(yīng)變曲線及屈服強(qiáng)度表現(xiàn)出的明顯靜水壓力與應(yīng)變率相關(guān)性；③ 強(qiáng)度不服從各向同性假定，抗拉強(qiáng)度遠(yuǎn)低于抗壓強(qiáng)度；④ 土體內(nèi)部的含大量水或不含水孔穴的力學(xué)行為十分復(fù)雜。為此提出各種理論模型描述土的力學(xué)性能。

為分析爆炸與沖擊荷載下土與結(jié)構(gòu)物的相互作用，Schwer等[4-10]采用失效土壤與可壓扁泡沫本構(gòu)模型、土/混凝土本構(gòu)模型、地質(zhì)帽子模型(*Mat _Geologic _Cap_Model)及連續(xù)光滑帽子模型等四種本構(gòu)模型進(jìn)行數(shù)值模擬。該四種模型均為與壓力相關(guān)的本構(gòu)模型，只在描述剪切破壞時(shí)應(yīng)力偏量第二不變量J2與壓力p 采用不同擬合關(guān)系。第一、二種本構(gòu)模型并未考慮剪切面與硬化面銜接，用體積應(yīng)變與壓力關(guān)系表述，且在較低靜水壓力狀態(tài)時(shí)數(shù)值迭代不穩(wěn)定，其應(yīng)用有一定局限性；第三種本構(gòu)模型則克服第一、二種的不足，在p-q表示的子午平面上增加一個(gè)帽子面，考慮剪切面與硬化面銜接，但二者相接時(shí)存在一階不連續(xù)性，造成數(shù)值計(jì)算時(shí)收斂困難等；第四種本構(gòu)模型基于Schwer 等[4]的光滑“cap”模型，克服第三種本構(gòu)模型的不足，并引入Rubin縮比函數(shù)，考慮拉壓強(qiáng)度差異及基于duvaut-lions類粘性應(yīng)變率效應(yīng)及能量理論的軟化損傷效應(yīng)。

國(guó)內(nèi)對(duì)土體在沖擊荷載作用下動(dòng)力學(xué)響應(yīng)研究成果較少，所用本構(gòu)方程大多為第一、二種，如劉彥等[11-12]的研究。理論上光滑帽子模型能較好描述土體動(dòng)態(tài)力學(xué)性能，但實(shí)際應(yīng)用如對(duì)侵徹深度影響驗(yàn)證及不同本構(gòu)模型、數(shù)值算法對(duì)侵徹深度影響尚未見報(bào)導(dǎo)。本文基于此進(jìn)行研究。

1 光滑帽子模型

光滑帽子模型由地質(zhì)帽子模型[4]發(fā)展而來。即將剪切屈服面與硬化帽子面用光滑曲線相連，克服因剪切面與硬化面相接時(shí)存在一階不連續(xù)造成數(shù)值計(jì)算收斂困難等問題，并考慮應(yīng)力張量第三不變量、應(yīng)變率效應(yīng)、損傷軟化影響，見圖1。

圖1 光滑帽子模型二維示意圖[10]

1.1 應(yīng)力不變量

屈服面由應(yīng)力張量第一不變量I1、偏應(yīng)力張量第二不變量J2及第三不變量J3表示。

1.2 塑性面

屈服面函數(shù)由三個(gè)應(yīng)力不變量及帽子面硬化參數(shù)表示[4,10]為

(1)

式中：Ff為剪切屈服面；Fc為硬化帽子；R為Rubin應(yīng)力縮放函數(shù)；κ帽子硬化參數(shù)，即帽子面與剪切面交點(diǎn)對(duì)應(yīng)的I1值。

1.3 剪切破壞面

低圍壓下材料剪切破壞[4,10]可表示為

Ff(I1)=α-λexp-βI1+θI1

(2)

式中：α，λ，β，θ為由常規(guī)三軸壓縮試驗(yàn)結(jié)果擬合的參數(shù)。

1.4 硬化面

由于地質(zhì)材料應(yīng)力-應(yīng)變曲線在低圍壓下表現(xiàn)出軟化特征，高圍壓下表現(xiàn)出硬化特征，當(dāng)材料處于由低圍壓至圍壓區(qū)域時(shí)，材料強(qiáng)度需結(jié)合帽子面及剪切面進(jìn)行模擬。帽子模型由兩分段函數(shù)表示。處于低圍壓狀態(tài)時(shí)帽子函數(shù)為“1”；由低圍壓至高圍壓區(qū)域時(shí)帽子函數(shù)為橢圓，屈服函數(shù)完全依靠帽子面及剪切屈服面函數(shù)[4,10]，可表示為

Fc(I1,κ)=

(3)

帽子面硬化規(guī)律[4,10]表示為

D2(X(κ)-X0)2]}

(4)

1.5 Rubin縮尺函數(shù)

Q1=α1-λ1exp-β1I1+θ1I1

(5)

Q2=α2-λ2exp-β2I1+θ2I1

(6)

式中：α1,λ1,β1,θ1,α2,λ2,β2,θ2均為擬合參數(shù)。

1.6 損傷軟化

光滑帽子模型中含兩類損傷[14-15]，即應(yīng)力張量第一應(yīng)力不變量為壓應(yīng)力時(shí)引入延性損傷使應(yīng)力折減；應(yīng)力張量第一應(yīng)力不變量為拉應(yīng)力時(shí)引入脆性損傷使應(yīng)力折減。損傷應(yīng)力定義為

(7)

損傷累積受拉時(shí)用脆性損傷方程表示；受壓時(shí)用延性損傷方程表示。

1.7 單元侵蝕

損傷參數(shù)d趨近于1時(shí)，單元會(huì)失去強(qiáng)度及剛度。為防止因剛度減小導(dǎo)致的計(jì)算困難，單元可據(jù)用戶需要進(jìn)行刪除或保留。d>0.99及最大主應(yīng)變超出用戶設(shè)定值時(shí)，單元會(huì)被刪除。

1.8 粘塑性率效應(yīng)

(8)

1.9 參數(shù)確定

光滑帽子模型參數(shù)主要通過靜水壓力試驗(yàn)、常規(guī)三軸壓縮試驗(yàn)、三軸拉伸試驗(yàn)、單軸應(yīng)變?cè)囼?yàn)、真三軸試驗(yàn)及SHPB試驗(yàn)確定。壓實(shí)度為95%、含水量12.54%(最大干密度γd,max=2.02 g/cm3，最佳含水量ωopt=12.87%)西南紅層泥巖粉碎土(土粒粒徑小于2 mm)的光滑帽子模型參數(shù)[6,16]見表1。

表1 光滑帽子模型參數(shù)

2 動(dòng)態(tài)本構(gòu)模型合理性分析

為論證提出的動(dòng)態(tài)本構(gòu)模型能正確反映土體本身力學(xué)性能進(jìn)行彈丸侵徹試驗(yàn)。采用數(shù)值模擬對(duì)彈丸侵徹試驗(yàn)重演，并分析本構(gòu)模型合理性。

2.1 彈丸侵徹試驗(yàn)

試驗(yàn)在空氣炮裝置上進(jìn)行，土體壓實(shí)度為95%，含水率12.54%，外側(cè)用鋼筒約束，內(nèi)徑尺寸為150 mm× 1 000 mm，壁厚5 mm，物理模型制備采用分層靜壓法。鋼質(zhì)彈丸直徑9.52 mm，侵徹速度v=200～300 m/s。試驗(yàn)試驗(yàn)裝置包括靜力及空氣炮加載裝置等。彈丸侵徹及開坑統(tǒng)計(jì)結(jié)果見表2。由表2看出，侵徹深度隨侵徹速度增加而增加，無論速度大小均存在開坑現(xiàn)象，且開坑大小與速度相關(guān)。侵徹速度為202 m/s時(shí)其平均侵徹深度為37.06 mm、平均開坑直徑19.71 mm、平均開坑深度5.00 mm；侵徹速度為310 m/s時(shí)其平均值分別為58.92 mm、41.98 mm、10.38 mm。

2.2 數(shù)值模型

圖2 數(shù)值計(jì)算模型

計(jì)算模型取物理模型的1/4(無特殊說明，計(jì)算模型下同)，見圖2，土體高度200 mm，半徑75 mm，彈丸與約束鋼筒為剛性體，土體、彈丸及約束鋼筒均為L(zhǎng)agrange網(wǎng)格(土體完全損傷后自動(dòng)刪除)，底部為透射邊界，彈丸與土體為侵蝕接觸，土體與約束鋼筒為面-面接觸，罰因子取2.0，土體壓實(shí)度95%，含水率12.54%，撞擊速度分202 m/s，310 m/s兩種，土體計(jì)算模型參數(shù)見表1。

表2 侵徹情況統(tǒng)計(jì)表

整理后彈丸侵徹速度與深度關(guān)系曲線見圖3。由圖3看出，撞擊速度310 m/s、202m/s彈丸剩余速度25.79 m/s、40.24m/s時(shí)，其對(duì)應(yīng)的侵徹深度為63.69 mm、44.74 mm；而模型試驗(yàn)確定的侵徹深度為58.92 mm、37.06 mm，計(jì)算結(jié)果基本合理。

圖3 不同撞擊速度下彈丸侵徹?cái)?shù)值計(jì)算結(jié)果

3 不同本構(gòu)模型對(duì)計(jì)算結(jié)果影響

通常采用Drucker-Prager模型、土壤與泡沫模型(Soil_and_Foam)、光滑帽子模型(Schwer_Murry_Cap_ Model)分析土壤在沖擊荷載下的動(dòng)力響應(yīng)，其剪切屈服面采用不同擬合關(guān)系。

Drucker-Prager屈服面[17-18]為

(9)

Drucker-Prager模型的破壞曲面為圓錐體，其大小可通過α，k調(diào)整。據(jù)錐面與摩爾-庫倫受拉子午線相互關(guān)系可確定α，k[17-18]為

(10)

式中：φ為摩擦角；c為粘聚力。

Soil_and_Foam屈服面[17-18]為

(11)

式中：a0，a1，a2為擬合參數(shù)；p=I1/3。

由以上分析知，光滑帽子模型考慮材料應(yīng)變率效應(yīng)、高圍壓下材料硬化及低圍壓下材料損傷軟化、應(yīng)力狀態(tài)對(duì)屈服強(qiáng)度影響及壓力導(dǎo)致的塑性體積應(yīng)變。應(yīng)變率效應(yīng)導(dǎo)致材料抵抗變形能力增大，而因材料硬化引入的帽子面、屈服面、材料損傷軟化及塑性體積應(yīng)變導(dǎo)致材料抵抗變形能力減弱；土壤與泡沫模型及Drucker-Prager模型均為理想的彈塑性本構(gòu)模型，其屈服面在π平面的形狀為圓形，均考慮因剪切破壞引起的剪脹效應(yīng)而未考慮應(yīng)力狀態(tài)對(duì)屈服面影響及應(yīng)變率效應(yīng)；兩模型不同處在于土壤與泡沫模型考慮壓力導(dǎo)致的塑性體積應(yīng)變。因此為研究不同本構(gòu)模型對(duì)侵徹過程影響，對(duì)彈丸侵徹試驗(yàn)進(jìn)行數(shù)值模擬，彈丸侵徹速度為310 m/s時(shí)計(jì)算模型同前，計(jì)算結(jié)果見圖6。由圖6看出，侵徹至相同深度時(shí)光滑帽子模型消耗能量最小，土壤與泡沫模型次之，Drucker-Prager模型最大，且土壤與泡沫模型與光滑帽子模型結(jié)果較接近，說明因材料硬化引入的帽子面、屈服面、材料損傷軟化及塑性體積應(yīng)變導(dǎo)致材料抵抗變形能力減弱效應(yīng)大于因應(yīng)變率引起材料抵抗變形增強(qiáng)效應(yīng)及塑性體積應(yīng)變對(duì)侵徹深度影響。由于Drucker-Prager模型及土壤與泡沫模型對(duì)剪切面采用線性擬合，所致剪脹效應(yīng)基本相同，但Drucker-Prager模型未考慮因壓力導(dǎo)致的土體積減小，因此Soil_and_Foam模型較Drucker-Prager模型引起的剪脹小(因屈服準(zhǔn)則參數(shù)不同)，故侵徹阻力亦小。Drucker-Prager模型中當(dāng)侵徹速度下降至54.8 m/s時(shí)，數(shù)值迭代出現(xiàn)不穩(wěn)定。

圖4 土壤與泡沫模型剪切面擬合參數(shù)

表3 Drucker-Prage、土壤及泡沫、 光滑帽子模型剪切面參數(shù)統(tǒng)計(jì)表

4 不同數(shù)值算法對(duì)計(jì)算結(jié)果影響

由于土體破壞為延性破壞，與混凝土與巖石等脆性材料相比，殘余強(qiáng)度與峰值強(qiáng)度相差較小，土體在沖擊荷載作用下易發(fā)生大變形，破壞后土體對(duì)侵徹過程影響不可忽略。Lagrange算法中以物質(zhì)坐標(biāo)為基礎(chǔ)，網(wǎng)格單元依附于物質(zhì)網(wǎng)絡(luò)，物質(zhì)不會(huì)在單元與單元間發(fā)生流動(dòng)。該方法優(yōu)點(diǎn)為能精確描述結(jié)構(gòu)邊界運(yùn)動(dòng)，但處理大變形問題時(shí)會(huì)出現(xiàn)嚴(yán)重的網(wǎng)格畸變現(xiàn)象。為使計(jì)算進(jìn)行須采用網(wǎng)格刪除技術(shù)；Euler方法以空間坐標(biāo)為基礎(chǔ)，劃分的網(wǎng)格及分析的物質(zhì)結(jié)構(gòu)相互獨(dú)立，網(wǎng)格在整個(gè)分析過程中始終保持最初的空間位置不動(dòng)，在整個(gè)數(shù)值模擬中各迭代過程計(jì)算數(shù)值精度不變；但該方法在物質(zhì)邊界的捕捉較困難, 用該方法時(shí)網(wǎng)格與網(wǎng)格間物質(zhì)可流動(dòng)；ALE兼具Lagrange方法與Euler方法特長(zhǎng)，即在結(jié)構(gòu)邊界運(yùn)動(dòng)的處理上引進(jìn)Lagrange方法特點(diǎn)，因此能有效跟蹤物質(zhì)結(jié)構(gòu)邊界運(yùn)動(dòng)；在內(nèi)部網(wǎng)格劃分上吸收Euler長(zhǎng)處，即使內(nèi)部網(wǎng)格單元獨(dú)立于物質(zhì)實(shí)體而存在，但其又不完全與Euler網(wǎng)格相同，網(wǎng)格可據(jù)定義的參數(shù)在求解過程中適當(dāng)調(diào)整位置，使網(wǎng)格不致出現(xiàn)嚴(yán)重畸變。此有利于分析大變形問題，用此方法時(shí)網(wǎng)格與網(wǎng)格間物質(zhì)也可流動(dòng)；SPH 算法優(yōu)點(diǎn)為無需網(wǎng)格、無單元網(wǎng)格畸變問題，能模擬爆炸引起的大變形，在模擬炸藥爆炸物質(zhì)飛濺方面更形象逼真；但SPH 算法精度不夠高、界面處理不成熟，對(duì)鄰粒子搜索需占用較多計(jì)算資源。

為研究不同數(shù)值算法對(duì)侵徹過程影響，對(duì)彈丸侵徹試驗(yàn)進(jìn)行數(shù)值模擬。彈丸侵徹速度310 m/s，計(jì)算模型同上，土體分別采用SPH、ALE及Lagrange網(wǎng)格，土體用光滑帽子模型為本構(gòu)模型，參數(shù)見表1，計(jì)算結(jié)果見圖7。由圖7看出，在MM_ALE算法中侵徹曲線斜率逐漸變小又逐漸增大，隨侵徹速度減小，土體應(yīng)變率效應(yīng)逐漸減弱，說明土體抵抗變形能力變?。欢?dāng)侵入到一深度時(shí)鋼筒約束作用逐漸增強(qiáng)，此時(shí)土體抵抗變形能力開始增大；而SPH，MM-ALE算法侵徹曲線斜率逐漸變緩，說明侵徹后期侵入深度相同，所需能量減小，與事實(shí)不符，其原因?yàn)長(zhǎng)agrange算法采用網(wǎng)格刪除技術(shù)，使其阻力減??；SPH算法雖未刪除粒子，但其計(jì)算精度較差。而采用Lagrange算法時(shí)適當(dāng)提高罰因子可有效彌補(bǔ)因網(wǎng)格刪除導(dǎo)致的阻力減小，使侵徹深度達(dá)到合理范圍，但并未改變侵徹速度、深度總體變化趨勢(shì)。

圖7 不同數(shù)值算法對(duì)侵徹過程影響

5 結(jié) 論

(1) 通過模型試驗(yàn)與數(shù)值方法對(duì)彈丸侵徹試驗(yàn)分析結(jié)果表明，光滑帽子模型能較好模擬在一定范圍內(nèi)沖擊荷載作用下土體的動(dòng)態(tài)力學(xué)響應(yīng)。

(2) 侵徹到相同深度時(shí)光滑帽子模型消耗能量最小，土壤與泡沫模型次之，Drucker-Prager模型最大，土壤與泡沫模型與光滑帽子模型結(jié)果較接近。說明因材料硬化引入的帽子面、屈服面、材料損傷軟化及塑性體積應(yīng)變導(dǎo)致材料抵抗變形能力減弱效應(yīng)大于因應(yīng)變率引起材料抵抗變形增強(qiáng)效應(yīng)及塑性體積應(yīng)變對(duì)侵徹深度影響。

(3) Drucker-Prager模型及土壤與泡沫模型對(duì)剪切面用線性擬合引起的剪脹效應(yīng)基本相同，但Drucker-Prager模型未考慮因壓力導(dǎo)致土體積減小，Soil_and_ Foam模型較Drucker-Prager模型所致剪脹小(因屈服準(zhǔn)則參數(shù)不同)，故侵徹阻力亦小。

(4) 在Drucker-Prager模型中侵徹速度下降到一定程度時(shí)，數(shù)值迭代出現(xiàn)不穩(wěn)定。

(5) 采用MM-ALE算法，侵徹速度-侵徹深度曲線斜率基本開始逐漸變小，后又逐漸增大，而Lagrange算法及SPH算法曲線斜率逐漸變緩；SPH算法計(jì)算精度較差。

(6) 用Lagrange算法時(shí)適當(dāng)提高罰因子可有效彌補(bǔ)因網(wǎng)格刪除導(dǎo)致的阻力減小，使侵徹深度達(dá)到合理范圍，但侵徹速度、深度總體變化趨勢(shì)未改變。

[1] 劉文韜.巖石含損傷本構(gòu)模型和地下爆炸效應(yīng)研究[D].合肥:中國(guó)科技大學(xué),2002.

[2] 王志亮,鄭田中,李永池.巖石時(shí)效損傷模型及其在工程爆破中應(yīng)用[J].巖土力學(xué),2007, 28(8): 1615-1620.

WANG Zhi-lang, ZHENG Tian-zhong, LI Yong-chi. Time effect damage model of rock mass and its application to blasting engineering[J]. Rock and Soil Mechanics, 2007, 28(8):1615-1620.

[3] 勞俊.地下空腔解耦爆炸的數(shù)值模擬[D].合肥:中國(guó)科技大學(xué),2008.

[4] Schwer L E,Murray Y D. A three-invariant smooth cap model with mixed hardening[J].International Journal for Numerical and Analytical Methods in Geomechanics,1994(18): 657-688.

[5] Fossum A F,Fredrich J T. Cap plasticity models and compactive and dilatant pre-failure deformation[J].Proceedings of the Fourth North American Rock Mechanics Symposium, 2000(7): 1169-1176.

[6] Livermore software technology corporation (lstc). Ls-dyna? keyword user’s manua[M]. Livermore: Livermore Software Technology Corporation, 2009.

[7] Katona M G,Mulert M A. A viscoplastic cap model for soils and rocks[C].//Proceedings of the International Conference on the Constitutive Laws for Engineering Materials. New York, NY: John Wiley & Sons Ltd., Mechanics of Engineering Materials, 1984: 335-350.

[8] Murray Y D. Numerical simulation of damage in concrete[R]. Alexandria: APTEK, Inc, Technical Report DNA-TR-94-190, Prepared for the Defense Nuclear Agency under Contract DNA 001-91-C-0075, 1985.

[9] Sandler I, dimaggio F L,Barron M L. An extension of the cap model: inclusion of pore pressure effects and kinematic hardening to represent an anisotropic wet clay[C].//Proceedings of the International Conference on the Constitutive Laws for Engineering Materials. New York, NY: John Wiley & Sons Ltd., Mechanics of Engineering Materials, 1984: 547-560.

[10] Schwer L E. Continuous surface cap model for geomaterial modeling: a new ls-dyna material type[C].// Proceedings of the 7th International LS-DYNA Users Conference, Dearborn, Michigan: LSTC & ETA, 2002: 16-50.

[11] 劉彥,黃風(fēng)雷,張振宇. 土壤覆層對(duì)混凝土中爆炸毀傷破壞的影響[J]. 北京理工大學(xué)學(xué)報(bào), 2007,27(9):765-768.

LIU Yan, HUANG Feng-lei, ZHANG Zhen-yu. Influence of soil cladding on explosion in concrete[J]. Transactions of Beijing Institute of Technology, 2007, 27 (9):765-768.

[12] 李小雷. 土壤混凝土分層介質(zhì)中爆炸的數(shù)值模擬[D].長(zhǎng)沙:國(guó)防科技大學(xué),2006.

[13] Simple R M.Convenient isotropic failure surface[J].Journal of the Engineering Mechanics Division, American Society of Civil Engineers, 1991(117): 348-369.

[14] Simo J C, Ju J W. Strain and stress based continuum damage models-Ⅰ formulation[J]. International Journal of Solids and Structures, 1987, 23(7):821-840.

[15] Simo J C, Ju J W. Strain and stress based continuum damage models-Ⅱ computational aspects[J].International Journal of Solids and Structures, 1987, 23(7): 841-869.

[16] 劉俊新.壓實(shí)粘性土動(dòng)態(tài)力學(xué)性研究[R].四川:中國(guó)工程物理研究院總體工程研究所,2011.

[17] 鄭穎人,沈珠江,龔曉南.巖土塑性力學(xué)原理[M].北京:中國(guó)建筑工業(yè)出版社,2002.

[18] 高紅,鄭穎人,馮夏庭. 巖土材料能量屈服準(zhǔn)則研究[J].巖石力學(xué)與工程學(xué)報(bào), 2007, 26(12): 2437-2443.

GAO Hong, ZHENG Ying-ren, FENG Xia-ting. Study on energy yield criterion of geomaterials[J].Chinese Journal of Rock Mechanics and Engineering, 2007, 26(12): 2437-2443.