穆朝民，戎立帆

（1.安徽理工大學(xué) 能源與安全學(xué)院，安徽 淮南 232001；2.湖南科技大學(xué) 煤礦安全開采技術(shù)湖南省重點實驗室，湖南 湘潭 411201；3.安徽理工大學(xué) 煤礦安全高效開采省部共建教育部重點實驗室，安徽 淮南 232001）

1 研究背景

磨料射流是將堅硬的磨料加入高壓水射流中充分混合后形成高速磨料與水的混合物。磨料射流目前已經(jīng)廣泛地應(yīng)用于石油、煤炭和冶金等行業(yè)，對于破煤/巖、除銹/清洗、冷切割都具有很重要的意義。巖石在高壓磨料射流沖擊下的損傷破壞及裂紋演化機制是磨料射流切割巖石技術(shù)的核心問題。

光滑粒子法（smoothed particle hydrodynamics,SPH）是近20年間發(fā)展起來的一種新的計算方法，其核心思想是將物質(zhì)點代替連續(xù)的物質(zhì)體用以承載物質(zhì)體的各種物理量（質(zhì)量、動量、能量和材料自身的特性）。由于沒有網(wǎng)格的存在，克服了網(wǎng)格畸變和滑移面問題，因此，特別適合計算沖擊大變形問題。但SPH法在計算效率和處理復(fù)雜邊界的能力要低于有限單元法（FEM法）。

近年來，光滑粒子法得到了較快的發(fā)展，馬利等[1－2]采用SPH方法模擬液體射流與剛性表面的碰撞，觀察液體在碰撞后進入一個相對較短的穩(wěn)定流動過程。肖毅華等[3]運用FEM-SPH耦合法模擬了具有自由表面的不可壓流體與結(jié)構(gòu)的相互作用問題。楊秀峰等[4]利用LANS-α模型對傳統(tǒng)的光滑粒子法進行改進，并研究了物體入水問題。汪繼文等[5]運用SPH方法結(jié)合虛粒子的方式對于潰壩問題中常見的漂浮物和障礙物進行了模擬。馬理強等[6]應(yīng)用光滑粒子流體動力學(xué)方法對潰壩問題進行了數(shù)值模擬。李大鳴等[7]的研究表明，SPH法對處理以一定的初速度沖擊盛水容器內(nèi)液面等自由表面大變形問題具有一定的優(yōu)勢。張志春等[8－9]提出了一種新型SPH-FEM耦合算法，并對圓柱形鋼彈正沖擊鋼板發(fā)生沖塞破壞的過程進行了三維數(shù)值模擬。

本文運用SPH建立水體模型、FEM法建立磨料與巖石模型，磨料均勻混合在水體中。通過Fortran語言編制SPH與FEM耦合算法對不同磨料射流速度、不同磨料濃度等條件下磨料射流對巖石的損傷機制及裂紋演化形式進行了探討。

2 本構(gòu)及損傷模型

巖石在沖擊作用下會相繼出現(xiàn)彈性、塑性，然后巖石中的孔隙和裂隙在巨大沖擊荷載作用下會出現(xiàn)結(jié)構(gòu)上的損傷，此損傷不可逆轉(zhuǎn)。本文計算中引入J-H-C含損傷本構(gòu)模型。

2.1 本構(gòu)模型

巖石采用Holmquist-Johnson-Concrete本構(gòu)模型[10－13]（圖 1(a)），其計算公式為

2.2 損傷模型

此損傷模型是把等效塑性應(yīng)變和塑性體積應(yīng)變引起的損傷都累積起來，如圖1(b)所示，計算公式為

式中：Δεp為等效塑性應(yīng)變（一個計算循環(huán)內(nèi)）；Δμp為塑性體積應(yīng)變（一個計算循環(huán)內(nèi)）；為塑性應(yīng)變；為塑性體積應(yīng)變；D1、D2均為損傷因子；T*為等效最大抗拉流體靜壓，T*= Tfc；T為煤體能夠承受的最大抗拉流體靜壓。當(dāng) P*=?T*時，煤體不能經(jīng)受任何塑性應(yīng)變；P*增加時，塑性應(yīng)變也增加。

圖1 J-H-C含損傷本構(gòu)模型Fig.1 J-H-C damage constitutive model

2.3 光滑粒子與有限元交界面上的滑移計算

如圖2所示，每一個光滑粒子都由一個圓代表，當(dāng)粒子i和j侵入滑移面M時，根據(jù)動量和動量距守恒定理對該粒子（由于這種侵入是非物理的，因此，要做必要的處理）和相關(guān)有限元節(jié)點（位于滑移面上）進行速度的調(diào)整（如圖3所示）：

圖2 Lagrange單元與光滑粒子的界面模型Fig.2 SPH node sliding on Lagrange grid

圖3 交界面附近粒子的計算示意圖Fig.3 Computations of the SPH nodes near sliding interface

式中：i為光滑粒子編號，1和2為粒子i相“接觸”兩個節(jié)點 M2和M3；V1、V2和Vi分別為線段 M2M3法線方向上節(jié)點 M2、M3及粒子i調(diào)整前的速度分量；V′1、V′2和V′i分別為調(diào)整后的法向速度分量；r為矢徑。由于M2、M3和i的速度分量線性分布，單元節(jié)點和光滑粒子的位置通過速度插值方式進行調(diào)整，使沿滑移面法向方向的粒子剛好回到與滑移面相切的位置。

2.4 水、磨料的狀態(tài)方程

水射流采用Bridgman根據(jù)試驗得到水介質(zhì)的狀態(tài)方程：

式中：ρ為密度；B、P0、ρ*、k1均為常數(shù)，X=30.5 GPa，P0=103.3 kPa，k1=7.147，ρ*=1.0 g/cm3。

磨料狀態(tài)方程采用Gruneisen方程：

式中：μ為體應(yīng)變；ρ0為初始密度；E0為初始內(nèi)能；γ0為Gruneisen參數(shù)；c為聲速；a、S1、S2、S3均為材料參數(shù)。彈體的材料參數(shù)如表1所示。

表1 磨料的材料參數(shù)Table 1 Material parameters for abrasive

3 計算模型與結(jié)果分析

3.1 計算模型

運用SPH與FEM固耦合罰函數(shù)算法對高壓水射流沖擊煤體進行二維數(shù)值模擬。計算模型如圖 4所示，巖石長、寬均為20 mm的正方形，左、右和下平面為透射面，上表面為自由面；水射流長、寬分別為14、4 mm。磨料尺寸分為兩種情況（如圖5所示），一種為直徑為2 mm的圓，一種邊長為2 mm的正方形。根據(jù)磨料與水的不同比例在水射流中分別布置1～3個磨料，對應(yīng)于磨料濃度（磨料體積與磨料射流總體積比）分別為6%、12%、18%，磨料射流出口靜壓分別取15、40 MPa，對應(yīng)磨料射流速度分別為173、282 m/s（圖中網(wǎng)格線代表標(biāo)尺）。水射流采用SPH建模，巖石采用FEM 建模，采用J-H-C含損傷本構(gòu)模型具體材料參數(shù)如表2所示。水射流長度方向共有70層，每層20個粒子，共1400個粒子；巖石為40000個四邊形單元；磨料為200個四邊形單元。

圖4 磨料水射流沖擊巖石模型Fig.4 Model of Abrasive water jet impact rock

3.2 結(jié)果分析

不同濃度與形狀的磨料射流以15 MPa出口壓力沖擊巖石的破巖深度如圖6所示。巖石在15 MPa高壓水射流沖擊下的破巖深度幾乎沒有變化（圖6(a)），可見15 MPa高壓水射無法破巖。單個磨料的磨料射流破巖深度為1.0 mm（圖6(b)），2個磨料的磨料射流破巖深度為2.3 mm（圖6(c)），3個磨料的磨料射流破巖深度為4.5 mm（圖6(d)），3個方形磨料的磨料射流破巖深度為3.0 mm（圖6(e)）。由此可見，磨料射流沖擊巖石的破巖深度與磨料射流的磨料濃度和磨料形狀有關(guān)，圓形或有尖銳棱角的磨料破巖效果優(yōu)于平面磨料。當(dāng)高壓射流可以持續(xù)沖擊巖石時，由于沖擊時間較長，所以磨料濃度不需要很高也可以有較高的破巖深度，但需要很長的沖擊時間。

圖5 15 MPa磨料射流沖擊巖石模型Fig.5 Model of abrasive jet water jet impact rock with 15 MPa

表2 巖石模型的材料參數(shù)Table 2 Material parameters for rock model

圖6 15 MPa磨料射流以不同磨料濃度與形狀沖擊巖石的深度Fig.6 Depths of penetration with different shapes and concentrations of abrasive with 15 MPa

如圖7所示，圓形磨料射流以40 MPa破巖的數(shù)值模型上布置20個測點（建模是軸對稱的，所以測點也是軸對稱的）進行煤體損傷狀態(tài)的研究，測點布置在XY平面上，測點1～5、6～10、11～15、16～20 的坐標(biāo)分別為(0,0)～(0,10)、(3,0)～(3,10)、(10,0)～(10,10)、(13,0)～(13,10)且為均勻布置。

圖7 磨料射流沖擊巖石模型Fig.7 Model of abrasive jet water jet impact rock

磨料射流以40 MPa沖擊巖石，不同時刻的破巖深度如圖8所示，損傷值線灰度對應(yīng)0～1。磨料射流沖擊到巖石表面時，首先高壓水體在巖體表面形成帶有擴散邊界的剪切層。巖體在射流沖擊下處于受壓狀態(tài)，沒有受到?jīng)_擊的巖體必然會抗擊這種壓縮而產(chǎn)生剪切破壞（如圖8(a)所示）。當(dāng)磨料開始沖擊巖石時，由于磨料的動能遠大于水體，造成沿磨料方向巖石損傷迅速加劇（如圖8(b)、(c)所示）。在距離沖擊面一定距離后，由于應(yīng)力波的反射拉伸作用產(chǎn)生拉伸損傷。隨著磨料射流對巖石的沖擊，巖石在軸向與徑向損傷都迅速增大，巖石在磨料流作用下的孔洞深度與寬度都在擴大，導(dǎo)致巖石損傷最終演化為“火山口”的形狀，且磨料射流在孔洞周邊形成的損傷由里向外逐步減輕，巖石在磨料沖擊軸向方向的損傷遠大于徑向方向。隨著沖擊作用時間的增加，磨料射流對巖石的沖擊成孔的效果迅速下降。如圖8(d)所示，在23.6 μs時，巖石沿軸向孔深為4 mm；如圖8(e)在67.3 μs時，巖石沿軸向孔深為7 mm；如圖8(f)在128.0 μs時，巖石沿軸向孔深為10 mm。由此可知，磨料水射流在破巖過程中能量衰減迅速，在距離首次沖擊點較遠處，高壓水射流的沖擊能量已經(jīng)達到巖石破壞的極值，因而形成的孔洞迅速穩(wěn)定。此時巖石的損傷破壞集中在底部，這是形成“火山口”孔洞的主要原因。水射流對巖石作用由最初的沖擊作用演化為靜壓作用，但靜水壓在巖石破碎中占比例有限。根據(jù)英國里茲大學(xué)關(guān)于高壓水射流沖擊巖體的試驗結(jié)果可知，巖石在高壓水射流作用下的成孔時間為毫秒量級，而磨料射流成孔時間要比水射流小一個量級。此結(jié)果與重慶大學(xué)司鵠等[13]得到的結(jié)果基本一致。

由圖9可知，巖石在磨料水射流作用下沿徑向和軸向損傷是呈階梯式的?？锥吹膹较蚝洼S向的損傷呈現(xiàn)“快慢交替”的變化。這主要是因為以下兩點：①巖石首先受到?jīng)_擊壓縮形成壓縮應(yīng)力波，然后巖石釋放彈性能形成拉伸應(yīng)力波，這種壓縮與拉伸交替的過程需要一定的時間，這一點從圖 10(e)也可以得出，測點12所在的巖石首先受到壓縮應(yīng)力的作用，然后巖石釋放壓縮能量形成拉應(yīng)力，巖石在回彈時又受到高壓水射流的壓縮，這樣循環(huán)往復(fù)造成了巖石的損傷；②巖石在水射流的作用下的受力形式?jīng)Q定了巖石損傷累積需要一個過程。

如圖10所示，軸向測點距離沖擊表面越近，損傷曲線的梯度越大，損傷演化的階梯形式越不明顯；距離沖擊表面越遠，損傷曲線的梯度越小，損傷演化的階梯形式越明顯。這主要是由于距離沖擊表面較近處，如圖所示巖石首先受到劇烈的壓縮作用開始損傷，此時損傷值在0.4左右，由于慣性作用巖石開始反向卸載形成強拉伸應(yīng)力波，巖石在拉伸應(yīng)力波作用下?lián)p傷加速至傷值為1.0的完全損傷。距離沖擊表面較遠處，如圖所示，由于巖石首先受的沖擊壓縮并不是很劇烈，所以反向形成較弱拉伸應(yīng)力波，巖石在壓縮與拉伸應(yīng)力波反復(fù)作用下形成階梯式的損傷。

對比數(shù)值模擬結(jié)果和吳海進[14]、林柏泉[15]（如圖11所示）關(guān)于磨料射流切割巖石的試驗結(jié)果（磨料的參數(shù)如表1所列，巖石的材料參數(shù)如表2所列）基本一致，磨料射流切割巖石形成“火山口”孔洞。由于磨料切割巖石試驗時可以對巖石進行持續(xù)沖擊，時間通常以min計量，而磨料射流破巖過程通常以μs（英國里茲大學(xué)關(guān)于高壓水射流沖擊巖體的試驗結(jié)果）計量，如果要建立現(xiàn)場持續(xù)磨料射流破巖數(shù)值計算過程則計算量過大。因此，現(xiàn)場試驗巖石切割深度大于數(shù)值模擬結(jié)果但巖石成孔形態(tài)基本與現(xiàn)場試驗一致。

圖8 40 MPa磨料射流作用下煤體的損傷演化Fig.8 Damage evolution of rock under high pressure water jet impacting with 40 MPa

圖9 損傷隨時間變化Fig.9 Damage evolution with time

圖10 垂直速度與損傷隨時間的變化Fig.10 Vertical speed and damage evolution with time

圖11 磨料射流現(xiàn)場試驗Fig.11 Abrasive water jet field experiment

4 結(jié) 語

（1）介紹了一種SPH（光滑粒子法）與FEM（有限元法）相耦合的算法，水射流采用SPH建模，磨料和巖石采用FEM 建模，該算法既克服了傳統(tǒng)有限元法網(wǎng)格畸變和滑移面問題，又克服了光滑粒子法計算效率和計算精度低等缺陷，特別適合研究磨料射流大變形物體高速沖擊巖石等問題。

（2）磨料射流侵徹巖石的深度與磨料射流中磨料的濃度、磨料形狀有關(guān)，巖石在磨料水射流作用下的沿徑向和軸向損傷是呈階梯式的?？锥吹膹较蚝洼S向的損傷呈現(xiàn)“快慢交替”的變化。

（3）距離沖擊表面越近，損傷曲線的梯度越大，損傷演化的階梯形式越不明顯；距離沖擊表面越遠，損傷曲線的梯度越小，損傷演化的階梯形式越明顯。

[1]馬利,王雙連,郭乙木.金屬液體射流變形的光滑粒子流體動力學(xué)模擬[J].科學(xué)通報,2007,52(2): 134－139.MA Li,WANG Shuang-lian,GUO Yi-mu.Simulation of liquid metal jet deformation smooth particle hydrodynamics[J].Science Bulletin,2007,52(2): 134－139.

[2]馬利.無網(wǎng)格法及液體射流高速碰撞與侵徹模擬[博士學(xué)位論文D].杭州: 浙江大學(xué),2007.

[3]肖毅華,韓旭,胡德安.流體與結(jié)構(gòu)相互作用問題的FESPH耦合模擬[J].應(yīng)用力學(xué)學(xué)報,2011,28(1): 13－18.XIAO Yi-hua,HAN Xu,HU De-an.Simulation fluidstructure interaction with FE-SPH method[J].Chinese Journal of Applied Mechanics,2011,28(1): 13－18.

[4]楊秀峰,彭世鏐,劉謀斌.物體入水的光滑粒子法模擬[J].計算物理,2011,28(4): 523－528.YANG Xiu-feng,PENG Shi-liu,LIU Mou-bin.Simulation of water entry with smoothed particle hydrodynamics method[J].Chinese Journal of Computational Physics,2011,28(4): 523－528.

[5]汪繼文,林昊.幾種潰壩問題的SPH方法數(shù)值模擬[J].計算機技術(shù)與發(fā)展,2011,21(5): 128－131.WANG Ji-wen,LIN Hao.The smoothed particle hydrodynamics method for solving some dam break problems[J].Computer Techhnology and Development,2011,21(5): 128－131.

[6]馬理強,常建忠,劉謀斌,等.基于SPH方法的潰壩流動數(shù)值模擬[J].水利水運工程學(xué)報,2010,(3): 65－70.MA Li-qiang,CHANG Jian-zhong,LIU Mou-bin,et al.Numerical simulation of dam-break flows using SPH method[J].Hydro-science and Engineering,2010,(3):65－70.

[7]李大鳴,李曉瑜,林毅.液滴沖擊自由液面的SPH法數(shù)值模擬[J].中國科學(xué),2011,41(8): 1055－1062.LI Da-ming,LI Xiao-yu,LIN Yi.Numerical simulation of droplet impacting liquid surface by SPH[J].Science China(Technological Sciences),2011,41(8): 1055－1062.

[8]張志春,強洪夫,高巍然.一種新型 SPH-FEM 接觸算法及其在沖擊動力學(xué)問題中的應(yīng)用[J].爆炸與沖擊,2011,31(3): 243－249.ZHANG Zhi-chun,QIANG Hong-fu,GAO Wei-ran.A new coupled SPH-FEM algorithm and its application to impact dynamics[J].Explosion and Shock Waves,2011,31(3): 243－249.

[9]張志春,強洪夫,高巍然.SPH-FEM 接觸算法在沖擊動力學(xué)數(shù)值計算中的應(yīng)用[J].固體力學(xué)學(xué)報,2011,32(3): 319－324.ZHANG Zhi-chun,QIANG Hong-fu,GAO Wei-ran.Application of SPH-FEM contact algorithm in impact dynamics simulation[J].Chinese Journal of Solid Mechanics,2011,32(3): 319－324.

[10]穆朝民,王海露.煤體在高壓水射流作用下的損傷機制[J].巖土力學(xué),2013,34(5): 1515－1520.MU Chao-min,WANG Hai-lu.Damage mechanism of coal under high pressure water jetting[J].Rock and Soil Mechanics,2013,34(5): 1515－1520.

[11]楊軍,高文學(xué),金乾坤.巖石動態(tài)損傷特性實驗及爆破模型[J].巖石力學(xué)與工程學(xué)報,2001,20(3): 320－323.YANG Jun,GAO Wen-xue,JIN Qian-kun.Experiment on dynamic damage property of rock and a new damage model for rock fragmentation by blasting[J].Chinese Journal of Rock Mechanics and Engineering,2001,20(3): 320－323.

[12]楊軍,金乾坤.應(yīng)力波衰減基礎(chǔ)上的巖石爆破損傷模型[J].爆炸與沖擊,2000,20(3): 241－246.YANG Jun,JIN Qian-kun.A new damage model for rock fragmentation by blasting based on stress wave attenuation[J].Explosions and Shock Waves,2000,20(3): 241－246.

[13]司鵠,謝廷明,楊春和.磨料水射流作用下巖石損傷場的數(shù)值模擬[J].巖土力學(xué),2011,32(3): 935－940.SI Hu,XIE Ting-ming,YANG Chun-he.Numerical simulation of rock damage field under abrasive water jet[J].Rock and Soil Mechanics,2011,32(3): 935－940.

[14]吳海進.高瓦斯低透氣性煤層卸壓增透理論與技術(shù)研究[博士學(xué)位論文D].徐州: 中國礦業(yè)大學(xué),2009.

[15]林柏泉,呂有廠,李寶玉,等.高壓磨料射流割縫技術(shù)及其在防突工程中的應(yīng)用[J].煤炭學(xué)報,2007,32(9):959－963.LIN Bai-quan,Lü You-chang,LI Bao-yu,et al.Highpressure abrasive hydraulic cutting seam technology and its application in outbursts prevention[J].Journal of China Coal Society,2007,32(9): 959－963.