， ， ，

(1.中國科學(xué)院力學(xué)研究所 流固耦合系統(tǒng)力學(xué)重點(diǎn)實(shí)驗(yàn)室 北京 100190； 2.華北水利水電大學(xué) 土木與交通學(xué)院 河南 鄭州 450045； 3.鄭州大學(xué) 水利與環(huán)境學(xué)院 河南 鄭州 450001)

0 引言

沖擊(爆炸)破巖涉及巖體中應(yīng)力波的傳播、巖體的損傷演化、巖體的破裂解體等多個(gè)過程，加載應(yīng)變率及能耗密度將直接影響巖體的動(dòng)態(tài)力學(xué)特性及破碎特性.霍普金森桿(SHPB)實(shí)驗(yàn)是研究巖石動(dòng)態(tài)力學(xué)特性及破碎特性的有效手段，國內(nèi)外的專家、學(xué)者利用SHPB對(duì)巖石的動(dòng)態(tài)力學(xué)特性進(jìn)行了大量的研究.文獻(xiàn) [1]利用SHPB對(duì)磁鐵礦的本構(gòu)關(guān)系進(jìn)行了研究，測(cè)得了磁鐵礦在不同應(yīng)變率下的應(yīng)力-應(yīng)變關(guān)系.文獻(xiàn)[2]通過動(dòng)態(tài)單軸拉伸試驗(yàn)和巴西劈裂試驗(yàn)的對(duì)比，得出巖石抗拉強(qiáng)度及彈性模量均隨應(yīng)變率的增加而增大，但彈性模量的增幅較小.文獻(xiàn)[3]研究表明，磁鐵礦、綠泥巖和混合巖的抗壓強(qiáng)度、彈性模量均隨應(yīng)變率的增加而增大，綠泥巖和混合巖的增幅大于磁鐵礦.文獻(xiàn)[4]研究表明，試件破碎能耗密度與入射能量呈線性正比關(guān)系，與試件平均應(yīng)變率呈冪函數(shù)關(guān)系；破碎能耗密度越大，試件破碎程度越劇烈，試件吸收能量主要耗散于巖石的損傷演化與變形破壞.文獻(xiàn)[5]得出磁鐵石英巖破碎能耗隨平均應(yīng)變率的增加而增大，破碎能耗密度和入射能量呈多項(xiàng)式關(guān)系，破碎后平均尺寸與能耗密度呈冪函數(shù)關(guān)系，能耗密度越大，巖石破碎程度越高.文獻(xiàn)[6]以若干硬巖試件在不同應(yīng)變率加載條件下的試驗(yàn)數(shù)據(jù)為基礎(chǔ)，分析了巖石動(dòng)態(tài)抗壓強(qiáng)度與靜態(tài)抗壓強(qiáng)度之比值與應(yīng)變率的相關(guān)性，結(jié)果表明，巖石強(qiáng)度與應(yīng)變率表現(xiàn)出較強(qiáng)或顯著相關(guān)的特性.文獻(xiàn)[7-8]對(duì)巖石在不同動(dòng)、靜組合加載下的力學(xué)特性進(jìn)行了分析，結(jié)果表明，巖石的組合加載強(qiáng)度大于其純靜載強(qiáng)度或純動(dòng)載強(qiáng)度，隨著沖擊動(dòng)載的增加，巖石組合加載強(qiáng)度增大.文獻(xiàn)[9]開展了單裂紋圓孔板砂巖試樣的沖擊壓縮實(shí)驗(yàn)，分別采用實(shí)驗(yàn)-數(shù)值法和準(zhǔn)靜態(tài)法確定了砂巖試樣動(dòng)態(tài)斷裂韌度.總體而言，研究者們利用SHPB實(shí)驗(yàn)證明了巖石動(dòng)態(tài)強(qiáng)度與應(yīng)變率的相關(guān)性，即隨著應(yīng)變率的增加，巖石的動(dòng)態(tài)強(qiáng)度將逐漸增大；此外，證明了隨著加載應(yīng)變率的增加，用于試樣破碎的能耗也逐漸增大，從而在表觀上呈現(xiàn)出隨著應(yīng)變率的增加，試樣破碎程度逐漸增高的現(xiàn)象.

然而，巖石應(yīng)變率的本質(zhì)是什么，是材料的本質(zhì)屬性，或僅僅是材料的宏觀響應(yīng)？國內(nèi)外的專家、學(xué)者并未對(duì)該問題進(jìn)行深入研究.數(shù)值計(jì)算方法能夠再現(xiàn)結(jié)構(gòu)受力及斷裂破壞的全過程[10-11].本文以SHPB實(shí)驗(yàn)的基本原理[12]為基礎(chǔ)，利用連續(xù)-非連續(xù)單元方法(CDEM)對(duì)赤鐵礦試樣的動(dòng)態(tài)力學(xué)特性進(jìn)行數(shù)值模擬研究，試圖從數(shù)值模擬的角度揭示巖石應(yīng)變率效應(yīng)的本質(zhì)，并分析能耗密度對(duì)巖石破碎程度的影響規(guī)律.

1 數(shù)值計(jì)算方法及巖體本構(gòu)

采用CDEM研究鐵礦石數(shù)值試樣在不同應(yīng)變率下的動(dòng)力響應(yīng)及破壞模式.CDEM是一種有限元與離散元耦合的顯式數(shù)值分析方法.數(shù)值模型由塊體及界面兩部分構(gòu)成.塊體由一個(gè)或多個(gè)有限元單元組成，用于表征材料的彈性、塑性、損傷等連續(xù)特征；兩個(gè)塊體間的公共邊界即為界面，用于表征材料的斷裂、滑移、碰撞等非連續(xù)特征.CDEM中的界面包含真實(shí)界面及虛擬界面兩個(gè)概念.真實(shí)界面用于表征材料的交界面、斷層、節(jié)理等真實(shí)的不連續(xù)面，其強(qiáng)度參數(shù)與真實(shí)界面的參數(shù)一致.虛擬界面主要有兩個(gè)作用：一是連接兩個(gè)塊體，用于傳遞力學(xué)信息；二是為顯式裂紋的擴(kuò)展提供潛在的通道(即裂紋可沿著任意一個(gè)虛擬界面進(jìn)行擴(kuò)展).因此，虛擬界面的強(qiáng)度參數(shù)與塊體的強(qiáng)度參數(shù)一致.

采用的巖體本構(gòu)為彈性-損傷-斷裂本構(gòu)，其中在每個(gè)有限元單元上施加線彈性本構(gòu)，輸入的參數(shù)包括密度、彈性模量及泊松比；在虛擬界面上施加損傷-斷裂本構(gòu)，輸入的參數(shù)包括法向連接剛度、切向連接剛度、黏聚力、內(nèi)摩擦角、抗拉強(qiáng)度、拉伸斷裂能及剪切斷裂能.

利用增量法表述的單元線彈性本構(gòu)為

(1)

式中：σij為應(yīng)力張量；Δσij為增量應(yīng)力張量；Δεij為增量應(yīng)變張量；Δθ為增量體應(yīng)變；K為體積模量；G為剪切模量；δij為Kronecker記號(hào)；t0表示當(dāng)前時(shí)步；t1表示下一時(shí)步.

虛擬界面上采用考慮強(qiáng)度線性軟化效應(yīng)的最大拉應(yīng)力模型及Mohr-Coulomb模型進(jìn)行損傷斷裂的計(jì)算.首先采用增量法計(jì)算虛擬界面下一時(shí)步的法向及切向試探接觸力，可以表示為

(2)

式中：Fn、Fs為法向、切向接觸力；kn、ks為單位面積上法向、切向接觸剛度；Ac為虛擬界面的面積；Δdun、Δdus為法向、切向相對(duì)位移增量.

拉伸破壞的判斷、法向接觸力及抗拉強(qiáng)度的修正可以表示為：

IF -Fn(t1)≥σt(t0)Ac,

THENFn(t1)=-σt(t0)Ac,

σt(t1)=-(σt0)2×Δun/(2Gft)+σt0,

(3)

式中：σt0、σ(t0)及σ(t1)為初始時(shí)刻、當(dāng)前時(shí)刻及下一時(shí)刻虛擬界面上的抗拉強(qiáng)度；Δun為當(dāng)前時(shí)刻虛擬界面上的法向相對(duì)位移；Gft為拉伸斷裂能.

剪切破壞的判斷、切向接觸力及黏聚力的修正可以表示為：

IFFs(t1)≥Fn(t1)×tanφ+c(t0)Ac，

THENFs(t1)=Fn(t1)×tanφ+c(t0)Ac，

c(t1)=-(c0)2×Δus/(2Gfs)+c0，

(4)

式中：φ為虛擬界面的內(nèi)摩擦角；c0、c(t0)及c(t1)為初始時(shí)刻、當(dāng)前時(shí)刻及下一時(shí)刻虛擬界面上的黏聚力；Δus為當(dāng)前時(shí)刻虛擬界面上的切向相對(duì)位移；Gfs為剪切斷裂能.

基于式(3)及式(4)，可繪制出虛擬界面上法向及切向的本構(gòu)曲線，具體如圖1所示.

2 數(shù)值模型及模擬參數(shù)

建立直徑為50 mm、厚度為25 mm的赤鐵礦數(shù)值試樣模型，如圖2所示，并采用63 900個(gè)四面體單元對(duì)該模型進(jìn)行剖分.

圖1 虛擬界面上的本構(gòu)曲線Fig.1 Constitutive curves in virtual interface

圖2 數(shù)值試樣模型Fig.2 Numerical model

數(shù)值模型執(zhí)行連續(xù)-非連續(xù)計(jì)算時(shí)，彈性單元本構(gòu)參數(shù)中，密度為3 300 kg/m3，彈性模量為60 GPa，泊松比為0.25；虛擬界面本構(gòu)參數(shù)中，法向、切向剛度為1016Pa/m，黏聚力為36 MPa，內(nèi)摩擦角為40°，抗拉強(qiáng)度為12 MPa，剪切斷裂能為2 000 Pa·m，拉伸斷裂能為200 Pa·m.

模型的底部為法向約束，模型頂部施加豎直向下的脈沖式載荷，參考SHPB的加載曲線，各應(yīng)變率下的加載持續(xù)時(shí)間均為250 μs，其中上升段及下降段的持續(xù)時(shí)間為30 μs，平臺(tái)段的持續(xù)時(shí)間為190 μs，探討了5種加載應(yīng)變率的影響，脈沖式加載時(shí)程曲線如圖3所示.

3 數(shù)值計(jì)算結(jié)果

3.1 不同應(yīng)變率下試樣動(dòng)態(tài)本構(gòu)曲線

在不同應(yīng)變率的加載作用下，對(duì)單元頂部的平均應(yīng)力進(jìn)行監(jiān)測(cè)，不同應(yīng)變率下的應(yīng)力-應(yīng)變曲線如圖4所示.由圖4可得，在線彈性階段，隨著應(yīng)變率的增大，應(yīng)力-應(yīng)變的關(guān)系基本一致；在非線性上升段及軟化段，應(yīng)變率越大，試樣的峰值強(qiáng)度越高，峰后的下降段越長.

圖3 不同應(yīng)變率下的脈沖式加載時(shí)程曲線Fig.3 Time-history curves of impulse load at different strain rates

圖4 不同應(yīng)變率下的應(yīng)力-應(yīng)變曲線Fig.4 Stress-strain curves at different strain rates

圖5 峰值強(qiáng)度與應(yīng)變率的關(guān)系曲線Fig.5 The relationship curve between peak strength and strain rate

3.2 不同應(yīng)變率下試樣峰值強(qiáng)度分析

對(duì)圖4中的峰值強(qiáng)度進(jìn)行統(tǒng)計(jì)，繪制峰值強(qiáng)度與應(yīng)變率的關(guān)系曲線，結(jié)果如圖5所示.由圖5可得，應(yīng)變率從40 s-1增加至200 s-1，峰值強(qiáng)度按冪函數(shù)形式從277 MPa增大至310 MPa(冪指數(shù)為2.33).

3.3 不同應(yīng)變率下試樣內(nèi)部能耗分析

不同應(yīng)變率下的輸入能量及應(yīng)變能時(shí)程曲線如圖6所示.可以看出，隨著加載時(shí)間的增加，外界輸入能量先緩慢增大，而后勻速增大，最后保持不變；而系統(tǒng)應(yīng)變能則呈現(xiàn)出先緩慢增加至峰值，而后緩慢下降的過程.此外，隨著應(yīng)變率的增大，輸入能量由86 J增加至340 J，而系統(tǒng)應(yīng)變能的峰值基本為40～50 J.

圖6 不同應(yīng)變率下的輸入能量及應(yīng)變能時(shí)程曲線Fig.6 Time-history curves of input energy and strain energy at different strain rates

圖7 能耗與應(yīng)變率的關(guān)系曲線Fig.7 The relationship curves between energy consumption and strain rate

300 μs加載結(jié)束后，外界輸入能量、殘存應(yīng)變能及損傷破裂能耗與應(yīng)變率的關(guān)系如圖7所示.其中，損傷破裂能耗可由輸入能量與殘存應(yīng)變能的差值獲得.可以看出，隨著應(yīng)變率的增加，外界輸入能量及損傷破裂能耗均逐漸增大，而系統(tǒng)殘存應(yīng)變能則逐漸減小.此外，外界輸入能量絕大部分用于試樣的損傷破裂，且應(yīng)變率越高，用于損傷破裂的能耗占比越大.

3.4 不同應(yīng)變率下試樣破壞程度分析

不同應(yīng)變率下的試樣破裂狀態(tài)如圖8所示.可以看出，在脈沖式加載的情況下，隨著應(yīng)變率的增加，試樣的破裂程度逐漸增大.應(yīng)變率為40 s-1時(shí)，試樣未發(fā)生破裂，應(yīng)變率為200 s-1時(shí)，試樣已經(jīng)發(fā)生了解體破壞.

定義破裂度為模型已經(jīng)破裂的虛擬界面面積與模型總虛擬界面面積的比值，則不同應(yīng)變率下的破裂度時(shí)程曲線如圖9所示.可以看出，隨著加載時(shí)間的增加，系統(tǒng)破裂度逐漸增大；當(dāng)加載時(shí)間超過250 μs后，系統(tǒng)破裂度基本保持不變.此外，由圖9還可以看出，應(yīng)變率越大，系統(tǒng)的破裂度越大.

圖8 不同應(yīng)變率下的試樣破裂狀態(tài)Fig.8 Fracture pattern of specimens at different strain rates

定義能耗密度為損傷破裂能耗與試樣體積的比值，表示為

ρw=Ef/V，

(5)

式中：ρw為能耗密度；Ef為損傷破裂能耗；V為試樣體積.

對(duì)圖9中的終態(tài)破裂度進(jìn)行統(tǒng)計(jì)，繪制破裂度與系統(tǒng)能耗密度的對(duì)應(yīng)關(guān)系曲線，對(duì)該曲線進(jìn)行線性擬合，可得破裂度與能耗密度的函數(shù)關(guān)系為

Fd=-0.21+0.068ρw，

(6)

式中：Fd為系統(tǒng)破裂度. 定義試樣的特征破碎尺寸為試樣體積與破裂總面積的比值，可以表示為

d=V/Af，

(7)

式中：d為試樣的特征破碎尺寸；Af為試樣破裂的總面積.

繪制試樣特征破碎尺寸與系統(tǒng)能耗密度的關(guān)系曲線，結(jié)果如圖10所示.可以看出，能耗密度與試樣破裂度及特征破碎尺寸密切相關(guān)，當(dāng)能耗密度增大至6.9 J/cm3時(shí)，試樣破裂度線性增大至25.5%，特征破碎尺寸則按冪函數(shù)形式迅速減小至1 mm.

圖9 不同應(yīng)變率下的破裂度時(shí)程曲線Fig.9 Time-history curves of fracture degree at different strain rates

圖10 特征破碎尺寸與系統(tǒng)能耗密度的關(guān)系曲線Fig.10 The relationship curve between characteristic size of fragments and system energy consumption density

4 結(jié)論

基于連續(xù)-非連續(xù)的數(shù)值模擬方法，以SHPB為模擬對(duì)象，研究了鐵礦石在不同應(yīng)變率下的應(yīng)力-應(yīng)變關(guān)系、峰值強(qiáng)度及破裂特征的演化規(guī)律，結(jié)果表明：在單元本構(gòu)及虛擬界面本構(gòu)中未引入應(yīng)變率效應(yīng)的情況下，數(shù)值試樣的動(dòng)態(tài)單軸抗壓強(qiáng)度表現(xiàn)出較強(qiáng)的應(yīng)變率相關(guān)性；隨著加載應(yīng)變率的增加，試樣的動(dòng)態(tài)強(qiáng)度按冪函數(shù)形式增大,數(shù)值試樣中用于損傷破裂的能耗密度也增大；能耗密度與試樣破裂度及特征破碎尺寸密切相關(guān)，隨著能耗密度的增大，特征破碎尺寸則按冪函數(shù)形式迅速減小.