梁 瑞，包 娟，周文海，2

(1.蘭州理工大學 石油化工學院，蘭州 730050； 2.蘭州大學 土木工程與力學學院，蘭州 730000)

1 研究背景

水電工程邊坡巖體爆破是水利水電工程中無法避免的關(guān)鍵問題，邊坡巖體損傷會造成巖石裂紋的激活、擴展和貫穿，導致巖體性能劣化最終破壞，繼而誘發(fā)邊坡失穩(wěn)[1-4]。爆炸過程中巖體損傷演化和裂紋分布問題復雜，對其直接理論解析難度大，前期學者常用實地測驗，但由于成本高、周期長、操作的可重復性等問題，試驗也較難開展，因此數(shù)值模擬一直是邊坡爆破損傷研究的重要方法之一[5-6]。

數(shù)值仿真可以很好地再現(xiàn)爆炸響應(yīng)的中間過程，并且能有效地描述邊坡巖土體在爆炸作用下的真實響應(yīng)情況，而損傷模型的選取直接影響爆破損傷數(shù)值模擬效果以及實測數(shù)據(jù)的驗證。劉亮等[7]結(jié)合聲波測試技術(shù)和實地試驗對白鶴灘邊坡保留巖體的損傷演化進行研究，確定了巖體處于臨界破碎狀態(tài)時的損傷變量及損傷閾值。趙小華等[8]引入RHT(Randomized Hough Transform)損傷模型，利用有限元軟件模擬沖擊載荷作用下混凝土壩的損傷情況，研究了孔口對混凝土重力壩動態(tài)響應(yīng)、抗爆性能及毀傷發(fā)展過程的影響。葉海旺等[9]選用HJC(Holmquist-Johnson-Cook)模型對分層裝藥中不同的填塞長度進行模擬，分析其對爆破效果的影響，得到實際工程最佳填塞長度。張小勇等[10]建立炸藥-空氣-混凝土的三維實體模型，以RHT模型為損傷本構(gòu)模型模擬混凝土爆破時的動態(tài)演變過程，得到襯砌結(jié)構(gòu)配筋率越高，圍巖損傷度越低。張若棋等[11]應(yīng)用有限元模擬軟件LS-DYNA和AUTODYN分析了HJC、RHT本構(gòu)模型參數(shù)，提出了確定模型失效強度參數(shù)的方法并利用實驗數(shù)據(jù)完成了參數(shù)的校驗。

總結(jié)前人研究發(fā)現(xiàn)，數(shù)值模型的選取與數(shù)值模擬結(jié)果的準確性相關(guān)聯(lián)，前人多采用傳統(tǒng)的損傷模型如HJC模型，對在HJC基礎(chǔ)上做出改善的RHT模型應(yīng)用較少。本文分別選取HJC模型和RHT模型對某水電站邊坡炮孔近區(qū)巖體建立爆破損傷模型，通過數(shù)值模擬結(jié)果和計算分析兩種模型下巖體的動力響應(yīng)，確定兩種模型下?lián)p傷范圍的分布規(guī)律、能量的衰減規(guī)律以及損傷速度閾值。研究成果可為選用不同巖石損傷模型進行爆破模擬研究提供參考。

2 數(shù)值模型的建立

2.1 損傷模型

HJC模型是一種綜合考慮應(yīng)變率效應(yīng)、損傷演化效應(yīng)、圍壓效應(yīng)和壓碎、壓實效應(yīng)影響的本構(gòu)模型。HJC模型屈服面方程有[12]

(1)

損傷因子D(0≤D≤1)為

(2)

數(shù)值模擬時HJC模型由*MAT_JOHNSON_HOLMGUIST_CONCRETE控制，詳細參數(shù)見表1。其中：Smax為名義最大強度系數(shù)；D1、D2為損傷參數(shù)；T為抗拉強度；Pc為空隙開始壓碎時的壓力；Pl為壓實時的壓力；μc為材料破壞時的體應(yīng)變；μl為固定體應(yīng)變；K1、K2、K3為壓力常數(shù)。

表1 HJC模型參數(shù)

RHT模型是在HJC模型的基礎(chǔ)上所提出的，引入了偏應(yīng)力張量第三不變量對破壞面形狀的影響。RHT模型[13]的失效面方程為

(3)

(4)

其中：

(5)

(6)

式中：*表示各物理量對于相應(yīng)材料單軸抗壓強度fc量綱歸一化后的值；Q2為拉壓子午比；rt、rc分別為拉、壓子午線處偏應(yīng)力；BQ為脆性-韌性轉(zhuǎn)化系數(shù)；Q2，0為材料常數(shù)；P為壓力。

(7)

式中B*、M為材料參數(shù)。

數(shù)值模擬時RHT模型由*MAT_RHT控制，詳細參數(shù)見表2。其中n為壓縮指數(shù)，A1、A2、A3、B1、B2、T1、T2為狀態(tài)方程參數(shù)。

表2 RHT模型參數(shù)

巖體初始密度均為2 660 kg/m3，剪切模量為17 GPa，抗壓強度為150 MPa，最小失效應(yīng)變?yōu)?0.01。

2.2 建立模型

以溪洛渡水電站巖石高邊坡的開挖為背景，其中爆破工藝參數(shù)和材料參數(shù)按照現(xiàn)場實際情況取值。圖1為邊坡以及邊坡炮孔近區(qū)巖體網(wǎng)格劃分示意圖，平臺寬度12 m，臺階坡面角65°，臺階長度 15 m，下層臺階長度8 m。為了對不同損傷模型下的巖石損傷特性進行研究，本文對其進行簡化處理，僅取邊坡炮孔近區(qū)5 m×5 m區(qū)域進行二維數(shù)值模擬，炮孔直徑為115 mm，采用PLANE162類型單元，網(wǎng)格選用四邊形劃分，炸藥和巖體網(wǎng)格單元分別為7 830和26 832。模型x和y方向施加平移約束，整體施加無反射邊界條件，以模擬無限巖體。

圖1 邊坡及邊坡近區(qū)巖體網(wǎng)格劃分Fig.1 Grid division of the slope and the rock mass in the vicinity of the slope

利用有限元軟件ANSYS/LS-DYNA和后處理軟件LS-PREPOST，選擇HJC模型和RHT模型本構(gòu)模型分別建立巖石單孔爆破損傷模型，在LS-PREPOST中利用關(guān)鍵字*ALE-MULTI-MATERIAL-GROUP設(shè)置流固耦合，將空氣、巖體、炸藥相結(jié)合。

炸藥選取2#巖石乳化炸藥，其材料類型為HIGH_EXPLOSIVE_BURN，狀態(tài)方程選取不考慮炸藥產(chǎn)物成分的JWL方程[14]，具體形式見式(8)。

(8)

式中：A、B、R1、R2、ω為試驗確定的常數(shù)；E0為初始比內(nèi)能；V為相對體積。具體參數(shù)如表3所示。

空氣選用*MAT_NULL，密度取1.29 kg/m3，狀態(tài)方程用EOS_LINEAR_POLYNOMAIAL來描述，即

P=C0+C1μ+C2μ2+C3μ3+

(C4+C5μ+C6μ2)eipvo。

(9)

式中：C0—C6為狀態(tài)方程參數(shù)；μ為比體積，取1.4；eipvo為內(nèi)能與初始體積之比。對于理想氣體，C1=C2=C3=C6=0，C4=C5=γ-1；其中γ為氣體的比熱，一般取1.4。

3 數(shù)值模擬結(jié)果與分析

水電工程邊坡的穩(wěn)定與巖體相關(guān)，本文主要對取不同損傷模型時巖體的劣化程度進行研究，故簡化模型只取邊坡炮孔近區(qū)巖體。HJC模型和RHT模型下的邊坡近區(qū)巖體在0.1、0.2、0.3、0.45、0.65、1 ms時的損傷云圖如圖2所示。圖2中巖體的宏觀損傷表明，隨著爆炸后能量在傳播過程不斷衰減，爆破由炸藥近區(qū)到中遠區(qū)的損傷程度依次降低。根據(jù)巖體損傷程度，可以簡單地將巖石的損傷區(qū)域劃分為粉碎區(qū)、裂隙區(qū)和彈性振動區(qū)。在爆炸初期，炸藥爆炸產(chǎn)生的沖擊波迅速作用于周圍巖體，并造成炸藥近區(qū)巖體粉碎，同時一部分沖擊波轉(zhuǎn)化為應(yīng)力波繼續(xù)傳播使巖體中遠區(qū)形成裂隙，并加速擴展，還有部分應(yīng)力波造成遠區(qū)振動。

圖2 不同模型下巖體損傷分布Fig.2 Distribution of rock mass damage in different models

表3 炸藥參數(shù)

HJC模型下的巖體損傷呈塊狀分布，且HJC模型下巖體的粉碎區(qū)面積較大，損傷區(qū)域分布明顯；而RHT模型下的損傷分布面積更大，近似于巖體的裂紋擴展，與Esen等[15]的實驗和模擬結(jié)果一致。

3.1 不同模型對巖體爆破開裂面積的影響

由圖2兩種模型的爆破損傷云圖可知，RHT模型下的爆破損傷分布形態(tài)更加接近巖體實際爆破時的破裂損傷形態(tài)，可以近似將其看作巖石裂紋的擴展。需要注意的是，RHT模型通過顯示材料的損傷狀態(tài)來描述材料的失效以形成視覺上的“裂紋”，但材料單元并沒有真實地開裂。同種工況下，HJC模型下巖體爆破開裂面積偏小，RHT整體爆破開裂面積更大，且分別取兩個模型不同時刻損傷分布尖端處的質(zhì)點，利用該質(zhì)點所處的位置計算損傷的面積。所選質(zhì)點至炮孔中心的距離為r，炮孔直徑為d，定義A=π(r2-d2/4)為巖體爆破開裂面積，可得如圖3所示的不同損傷模型下爆破開裂面積的時程變化曲線，以爆破作用時間為橫坐標，爆破開裂面積為縱坐標，做損傷面積和時間的擬合曲線，兩個模型的擬合度均達到0.99。

圖3 不同損傷模型下的爆破開裂面積分布Fig.3 Blasting cracked area in different damage models

對圖3中HJC模型和RHT模型的巖體爆破開裂面積進行分析，發(fā)現(xiàn)在0.4 ms前RHT模型的爆破開裂面積小于HJC損傷模型下的爆破開裂面積，而在0.4 ms后RHT爆破開裂面積迅速增加，并且最終形成的爆破開裂面積大致為HJC模型下爆破開裂面積的一倍。1 ms時RHT模型下巖體爆破開裂面積為13.84 m2，而HJC模型下巖體爆破開裂面積為7.44 m2。由此可以看出，同種工況下，巖石損傷類型決定了巖體爆破開裂面積的大小。

由上文擬合結(jié)果可建立本文HJC模型和RHT模型下巖石爆破開裂面積預(yù)測公式，其中x表示時間(ms)，y表示爆破開裂面積(m2)，即：

HJC模型：y=7.52x1.45；

RHT模型：y=0.26-3.4x+21x2-3.9x3。

3.2 不同模型對爆炸地震波能量的影響

為了進一步探究不同損傷模型下水電工程邊坡近區(qū)巖體爆破開裂面積產(chǎn)生差異的原因，從能量的角度出發(fā)，分析能量對巖體爆破損傷產(chǎn)生的影響。圖4為取兩個爆破模型0.65 ms時損傷區(qū)外沿輪廓質(zhì)點H27009的振動速度時程變化曲線。圖4中HJC模型的質(zhì)點的最大振速明顯大于RHT模型的最大振速，HJC爆破模型在0.3 ms達到速度峰值23.7 cm/s，RHT模型在0.5 ms時達到速度峰值12.2 cm/s。

圖4 不同模型質(zhì)點H27009的振動速度時程曲線Fig.4 Time-history curves of vibration velocity of particle H27009 in different models

爆破產(chǎn)生的爆炸沖擊波由于能量衰減成應(yīng)力波并在后期轉(zhuǎn)化為彈性地震波，一般認為非彈性區(qū)和彈性區(qū)交界面所包圍的空間為爆炸地震波的產(chǎn)生區(qū)。巖體損傷爆破開裂的改變勢必會對損傷區(qū)外爆炸地震波的輻射和傳播產(chǎn)生影響，爆破損傷分布形態(tài)與爆炸地震波互相作用，不同損傷模型下的巖體損傷分布不同，其能量的衰減大小也有所差異。本文的爆炸地震波能量的計算采用Sanchidrián等[16]提出的方法，即

(10)

式中：ES為監(jiān)測點處單位面積上某時刻的地震波能量；ρ為傳播介質(zhì)的密度；cL為傳播介質(zhì)的縱波波速；cT為傳播介質(zhì)的橫波波速；v1、v2、v3分別為監(jiān)測點處水平徑向、水平切向和垂直方向的振動速度。對于平面應(yīng)變模型，v3=0，故式(10)可簡化為

(11)

同樣地，分別取2個模型不同時刻損傷分布尖端處的質(zhì)點做監(jiān)測點，根據(jù)式(11)計算巖體損傷分布尖端處監(jiān)測點的地震波能量，得到不同損傷模型下?lián)p傷分布尖端處監(jiān)測點的爆炸地震波能量時程曲線，如圖5所示?？梢园l(fā)現(xiàn)，HJC模型巖體的損傷分布監(jiān)測點振動波能量峰值更大，在0.3 ms時達到峰值6.78×108J；而RHT模型巖體在0.5 ms時達到能量峰值，為3.55×108J，RHT模型下振動波能量峰值是HJC模型地震波能量峰值的52.3%。

圖5 不同損傷模型下能量衰減規(guī)律Fig.5 Energy attenuation law under different damage models

HJC模型的巖體損傷爆破開裂面積小于RHT模型下的爆破開裂分布面積，但損傷分布尖端監(jiān)測點的彈性振動波能量大于RHT模型。這是由于炸藥爆炸能向爆炸地震波能量轉(zhuǎn)化的比例取決于爆破損傷分布的范圍，損傷范圍越小，所需的炸藥沖擊波能量越小，炸藥沖擊波轉(zhuǎn)化為彈性振動波的比例越大。因此，HJC模型中，更多的炸藥爆炸能轉(zhuǎn)化為爆炸振動波能量，而RHT模型中，更多的能量被用于損傷的擴展。

3.3 不同損傷模型下PPV閾值的探究

在實際工程中廣泛利用易于監(jiān)測的PPV(質(zhì)點峰值速度)作為巖體損傷的判據(jù)，且國內(nèi)外學者大多以質(zhì)點峰值振動速度作為爆破安全控制的指標。巖體損傷模型的選取會影響巖體爆破開挖時巖體損傷面積的分布和損傷演化機理，這勢必也會造成巖體爆破損傷PPV的改變，因此本節(jié)將對不同模型下的PPV閾值進行探究。根據(jù)巖體彈性損傷理論以及損傷變量D與聲波降低率η的關(guān)系[17]可得

D=1-(1-η)2。

(12)

本文認為當η>10%時，巖體開始出現(xiàn)爆破損傷，利用式(12)計算得到其對應(yīng)的巖體損傷度閾值Dcr=0.19；當η>65%時，認為巖體已經(jīng)出現(xiàn)嚴重損傷，此時損傷度D=0.8；當損傷度D=1時，認為巖石完全損傷。根據(jù)以上關(guān)系，以炮孔為圓心，沿著同一方向，每隔0.25 m取一個質(zhì)點，共取20個，分別提取每個質(zhì)點的損傷度和質(zhì)點峰值速度，建立巖體損傷度D與質(zhì)點峰值速度PPV之間的關(guān)系，如圖6所示。經(jīng)過擬合發(fā)現(xiàn)兩種損傷模型巖體的PPV與損傷度D之間的變化規(guī)律基本一致，均符合logistics方程曲線規(guī)律。

圖6 損傷度與PPV的關(guān)系Fig.6 Relationship between damage degree and PPV

當損傷度D達到損傷度閾值Dcr時，巖體開始損傷，PPV閾值也出現(xiàn)。由圖6中巖體PPV與損傷度D的關(guān)系，找到擬合曲線上損傷度閾值0.19對應(yīng)的PPV閾值，得到HJC模型和RHT模型的巖體爆破模型出現(xiàn)損傷的PPV閾值分別為13.09 cm/s和27.29 cm/s，RHT損傷模型的損傷閾值大于HJC損傷模型的損傷閾值，損傷閾值變化率達到52%；當損傷度處于0.19～0.80，巖體處于裂隙產(chǎn)生、擴展、破碎階段，此時HJC模型的PPV的范圍為32.36～96.87 cm/s，RHT模型的PPV范圍則為73.81～116.47 cm/s。

通過圖6和表4不難發(fā)現(xiàn)，由于HJC模型損傷速度閾值小，其巖體損傷粉碎區(qū)形成更快，而由于RHT損傷模型損傷擴展持續(xù)時間更長，PPV值更大，致使RHT模型下?lián)p傷能充分擴展，最終形成更大的損傷區(qū)域。可見巖體的PPV跟巖體的材料類型息息相關(guān)，因此在實際工程中，要根據(jù)實際情況確定安全控制的指標。

表4 不同損傷模型下PPV的變化

4 結(jié) 論

本文針對水電工程邊坡炮孔近區(qū)分別建立了以HJC模型和RHT模型為損傷模型的巖石單孔爆破數(shù)值模型，探究了不同損傷模型下巖體的動力響應(yīng)和斷裂損傷特性。通過以上計算分析，得到以下結(jié)論：

(1)HJC模型損傷云圖呈塊狀分布，RHT模型損傷云圖與巖體破碎的裂紋擴展相似。同種工況下，所選的巖體損傷模型決定了損傷面積的大小，0.4 ms之前HJC模型損傷面積大于RHT模型的損傷面積；但在0.4 ms之后RHT模型損傷面積均大于HJC模型的損傷面積，最終RHT模型形成的損傷面積達到HJC模型損傷面積的一倍。

(2)HJC模型損傷分布尖端監(jiān)測點的彈性振動波能量大于RHT模型的彈性振動波能量；損傷范圍越小，炸藥沖擊波轉(zhuǎn)化為彈性振動波的比例越大。

(3)當HJC模型的PPV的范圍為32.36～96.87 cm/s，RHT模型的PPV范圍為73.81～116.47 cm/s時巖體處于裂隙產(chǎn)生、擴展、破碎階段。巖石HJC模型的損傷PPV閾值為13.09 cm/s，RHT模型的損傷PPV閾值為27.29 cm/s。

根據(jù)以上討論，發(fā)現(xiàn)RHT模型模擬巖體爆破的損傷效果更貼合實際工程爆破特點，因此水電工程邊坡爆破時周圍巖體損傷特性以及預(yù)留坡體的穩(wěn)定性研究推薦使用RHT模型。本文研究了炮孔近區(qū)均質(zhì)巖體爆破特性，未來將繼續(xù)對整個邊坡和非均質(zhì)巖體爆破損傷特性深入研究。