基于ALE 算法的連續水射流破煤特征數值模擬研究

郝從猛 ，陳學習 ，劉清泉 ，金霏陽

（1.中國安全生產科學研究院，北京 100012；2.華北科技學院 安全監管學院，河北 三河 101601；3.中國礦業大學 安全工程學院，江蘇 徐州 221116）

龐大的煤炭產能有力支撐了我國經濟的高速發展，同時旺盛的煤炭需求使得我國礦井的開采正在以平均每年10～30 m 的速度向深部延伸，許多礦井的開采深度已超過800 m，特別是山東、安徽、黑龍江、河南等地區多對礦井已超過1 000 m[1-2]。隨著開采深度的不斷增加，煤層地應力、瓦斯壓力和含量急劇增大，而滲透率卻大幅降低，此“三高一低”的狀況使得越來越多的非突礦井正在轉變為突出、強突出礦井，煤礦企業面臨的安全生產形勢愈發嚴峻[3]。作為煤炭生產、消費大國，我國每年投入巨額經費治理煤礦瓦斯，盡管成績顯著，但瓦斯事故依然是引起煤礦重大傷亡的主要原因之一[4]。

實踐表明，瓦斯抽采是解決煤礦瓦斯事故，實現節能環保的最有效措施[5]。但普遍存在的低滲煤層使得瓦斯抽采一直面臨著抽采效率低、難度大、治理成本高等問題[6]；特別是進入深部以后，高地應力作用使得瓦斯抽采瓶頸日益顯著。為了提高煤層滲透性，廣大的科研工作者們提出了很多的增透措施和技術；其中，水力沖孔技術由于操作簡單、卸壓效果顯著而得到了廣泛的應用，相關研究成果十分豐富。劉明舉等[7]、王兆豐等[8]、陶云奇等[9]通過現場考察、實驗研究獲得了水射流破煤條件和增透效果，證明煤體卸壓增透效果和沖孔半徑成正相關；ZHANG 等[10]結合應用開發了高壓水射流鉆沖一體化造穴裝備，在新景礦試驗成功后進行了大規模推廣；葛兆龍等[11]、高亞斌等[12]通過試驗研究了水射流沖擊性能，得到水射流的破煤細觀演化機制。梳理文獻可以發現，目前水射流破煤方面的研究主要集中在破煤機理及應用優化上，對于水射流破煤特征也多是對非淹沒射流條件下進行了分析。

然而，在一些俯孔條件下進行水射流破煤時往往出現孔內積水不能及時排出，進而形成淹沒環境下射流破煤的情況[4]。由于淹沒環境下水射流能量損耗與非淹沒環境差距較大，導致兩者破煤性能具有顯著的差別，而目前的研究對兩者破煤特性認識不足，特別是對淹沒射流破煤效率及其最佳沖擊深度研究不足，這也導致在淹沒環境下開展水射流破煤造穴技術的應用存在一定的盲目性。為此，采用數值模擬的方法對淹沒射流/非淹沒射流的破煤特征進行對比分析，進一步揭示不同環境介質下水射流破煤特性的差異及其產生機制，并對淹沒環境下水射流破煤時效性及可行性展開分析。

1 基于ALE 的連續水射流破煤特性數值模型

1.1 模擬軟件及流程

多物質流固耦合法（ALE）是一種任拉格朗日-歐拉方法，可同時實現單物質和多物質的任意耦合，可以較好地實現連續射流沖擊的數值計算方法。ALE 算法兼顧了Lagrange 和Euler 方法兩者的優點，不僅具備Lagrange 方法在有效跟蹤物質在結構邊界上運動的優點，還具有Euler 方法在內部網格與物質實體之間相互獨立的強項[13]。ALE算法不僅確保了計算過程中的網格單元的完整性，同時縮短了復雜問題的求解時間并提高了精度，特別適用于沖擊作用下大變形問題的處理。

LS-Dyna 模擬軟件可實現多種算法，并且擁有豐富的材料庫和材料本構模型，可以提供50 多種材料間的接觸方式，適合分析高速碰撞、侵徹等沖擊問題。因此，采用LS-Dyna 進行水射流沖擊破煤模擬，結合Ls-Prepost 軟件的前后處理功能，進行建模、修改參數、結果顯示、結果輸出等。數值模擬流程如圖1。

1.2 幾何模型及材料參數

由于需要對非淹沒/淹沒環境下的破煤特性進行對比分析，因此幾何模型需要建立不同的射流域，即煤體分別處于空氣域和水域中。對于對稱問題，為了提高計算效率可以進行適當的精簡取其一部分進行分析。高壓水射流破煤過程就是對稱的，因此取模型的1/4 進行模擬計算，在結果分析時按對稱面呈現，非淹沒/淹沒射流破煤數值計算1/4 模型如圖2。左圖為空氣域（非淹沒）模型，右圖為水域（淹沒）模型，煤體的尺寸都為100 mm×100 mm×100 mm。

圖2 非淹沒/淹沒射流破煤數值計算1/4 模型Fig.2 Numerical calculation of 1/4 model of coal breaking by unsubmerged/submerged water jet

煤體本構模型采用高應變率、大變形下的巖石材料的Johnaon-Holmquist-Concrete（J-H-C）模型，根據煤層特點并結合相關文獻[4,14]，煤體選取J-HC 模型，研究表明該模型適合模擬水射流沖擊下煤體發生損傷的標準，即以最大剪切應變和最大拉伸應變為判斷依據的破壞條件；流體（水和空氣）選取MAT-NULL 模型。

煤體材料模型（H-J-C）參數如下：①煤體密度ρc：1.4 g/cm3；②剪切模量G：1.48 GPa；③特征化黏性強度A：0.79；④特征化壓力硬化因子B：1.6；⑤應變率影響系數C：0.007；⑥壓力硬化指數N：0.61；⑦無側限單軸抗壓強度Fc：4.8×10-3GPa；⑧抗拉強度T：4×10-4GPa；⑨累計塑性應變EFmin：0.001；⑩特征化最大強度SFmax：7；?壓潰點的壓力pc：-1.6×10-4GPa；?壓潰點的體積應變Uc：0.001；?壓實點的壓力pL：0.008 GPa；?壓實點的體積應變UL：0.1；?損傷常數D1：0.04；?損傷常數D2：1.0；?壓力常數K1：0.85 GPa；?壓力常數K2：-1.71 GPa；?壓力常數K3：2.08 GPa；?無量綱系數：3.57。

流體材料模型（NULL）空氣參數如下：①空氣密度ρa：1.9×10-3g/cm3；②空氣壓力截止pca：-1×10-11GPa；③空氣動力黏度系數Ca：0.033 1；④無量綱系數S1a：0；⑤無量綱系數S2a：0；⑥無量綱系數S3a：0；⑦狀態方程系數GAMA0a：0；⑧一階體積更正系數αa：0；⑨空氣初始內能Ea：0。

流體材料模型（NULL）水參數如下：①水密度ρw：1.0 g/cm3；②水壓力截止pcw：-1×10-5GPa；③水動力黏度系數Cw：0.148；④無量綱系數S1w：2.56；⑤無量綱系數S2w：1.986；⑥無量綱系數S3w：0.226；⑦狀態方程系數GAMA0w：0.45；⑧一階體積更正系數αw：0；⑨水初始內能Ew：0。

1.3 控制方程

模型中流體的連續性方程描述如下[15]:

式中：ρ為流體介質的密度，kg/m3；t為時間；vi為i方向的速度分量，m/s；wi為i方向的相對速度的分量，m/s；xi為i方向的位移分量，m。

射流的動量方程如下：

式中：v為射流速度，m/s；σij,j為應力張量在j方向的偏導數；bi為i方向的體積力，N；xj為j方向的位移分量，m。

射流的能量方程為：

式中：E為流體的內能密度，J/m3；σi,j為在i方向和j方向體應力矢量和；vi,j為i方向和j方向的速度矢量和；wj為j方向的相對速度的分量，m/s。

通過與狀態方程聯用而為壓力行為提供應力組件，從而為材料提供應力張量。對空氣域選取*EOD-LINEAR-POLYNOMIAL 狀態方程來提供壓力，將空氣域視為理想氣體后可得狀態方程為[14]：

式中：p為壓力；γ為空氣的比熱比，γ=CP/CV；CP、CV分別為定壓比熱容和定容比熱容；ρ0為初始密度。

對于水域，采用*EOS-GRUNEISEN 狀態方程來提供壓力p，表達式為：

式中：Cb為沖擊波波速與單元速度關系曲線的縱軸截距；S1、S2、S3為沖擊波波速與單元速度關系曲線的斜率系數；γ0為Gruneisen 系數；μ為泊松比；α為對 γ0的一階體積修正系數。

煤的本構模型選用H-J-C 模型，模型的規范化等效應力描述如下：

式中：A為無量綱黏性力；B為規范化壓縮系數；p*為無量綱壓力；為無量綱應變率；C為材料的應變系數；N為硬化系數；D為損傷因子。

式中：Δεp為等效塑性應變增量，無量綱；Δμp為等效體積應變增量，無量綱；D1（p*+T*）D2為常壓p下材料斷裂時的塑性應變，無量綱；T*為最大拉伸靜水壓力，MPa；D1、D2為損傷常數。

2 連續水射流破煤特性數值模擬結果

2.1 射流能量耗散分析

通過分析淹沒/非淹沒狀態下水射流的傳播速度可以初步獲得兩者的傳播特征和能量耗散情況。設定射流直徑為0.3 mm，射流速度為170 m/s，噴嘴與煤壁作用面的間距都為20 mm，水射流在非淹沒和淹沒條件下的傳播情況如圖3。

圖3 非淹沒/淹沒水射流傳播特征Fig.3 Propagation characteristics of unsubmerged/submerged water jet

從圖3 中可以看出：非淹沒和淹沒條件下水射流自噴嘴噴出后兩者的傳播過程差別巨大；首先，從水射流形態上看，非淹沒射流傳播過程中呈圓柱狀向前傳播，射流直徑基本上沒有變化，而淹沒射流的傳播過程中的前端由于水域的阻力作用逐漸向后分離，同時水射流內部受擠壓而產生應變，這是水射流周圍部分水速低于射流速度后導致的滯后效應，射流受阻是導致淹沒條件下水射流的傳播距離較小的根本原因；其次，射流在傳播速度上的明顯差異，非淹沒射流的傳播距離和傳播時間成正比，而相同時間淹沒射流的傳播距離大幅下降；根據射流到達煤壁的時間，淹沒條件下水射流傳播10 mm 的距離所用時間是非淹沒條件下的1.6 倍，傳播20 mm 的距離所用時間是非淹沒條件下的2.18 倍，這說明隨著傳播距離的增大淹沒條件下的水射流要不斷消耗自身的質量來抵消環境的阻力，也就是說淹沒條件下需要水射流持續不斷地噴出才能保證射流向前傳播，否則射流就會在環境水域的阻力消耗下消失。

另外，根據水射流傳播的時間變化可以推測，隨著射流噴嘴與煤體間距的增大所需射流的持續時間也會變長。

2.2 淹沒/非淹沒射流破煤特征分析

2.2.1 水射流破煤過程分析

為了進一步研究射流環境對破煤的影響，建立射流直徑為3 mm、靶距為20 mm 的破煤工況，沖擊計算時間1 000 μs，射流速度分別為170、190、210 m/s 時非淹沒/淹沒射流破煤過程如圖4。

圖4 非淹沒/淹沒射流破煤過程Fig.4 Coal breaking process of unsubmerged/submerged water jet

從圖4 中可以看出：非淹沒/淹沒條件下水射流破煤特性存在著明顯的區別；首先，3 組非淹沒射流從噴嘴噴出后到達煤壁的時間不同；另外，非淹沒射流模擬試驗的煤體全部被穿透，水射流速度為170、190、210 m/s 時煤體穿透的時間分別是970、900、830 μs，而淹沒條件下3 組試驗全部沒有穿透煤體，但隨著射流速度的增大破煤深度逐漸增加；最后，除破煤深部不同外，2 種射流的破煤基坑直徑的差別明顯，在淹沒射流模型中，隨著射流時間的延長，侵蝕坑呈漏斗狀逐漸向下延伸，但是非淹沒射流的破煤侵蝕坑始終保持圓柱狀增長，這可能是由于在水域環境下射流周圍產生強烈的湍流作用使得射流的動量不斷向周圍水分子傳遞，進而二次沖擊使得破煤范圍增大，同時，湍流作用下的返流動量會很大程度上減弱射流自身的沖擊能力，使得射流向前沖擊煤體的能力變弱，最終形成漏斗形的基坑。

2.2.2 煤體應力分布特征分析

非淹沒/淹沒射流破煤時煤體等效應力分布如圖5。

圖5 非淹沒/淹沒射流破煤時煤體等效應力分布Fig.5 Distribution of equivalent stress of coal mass during unsubmerged/submerged water jet breaking coal

從圖5 可以看出：2 種射流模式下的最大應力非常接近。這主要是因為應力的峰值主要是由于煤體自身強度決定的，與射流模式關系不大，但是在應力分布上差別明顯；非淹沒條件下的應力主要在基坑的前方，沿橫向上分布范圍較小；而淹沒條件下應力不僅在沖擊位置分布范圍較大，在煤體的底部也出現了應力的升高，這主要是淹沒射流擠壓煤體的范圍較大，而壓強達到煤體應力極限后才會出現損傷，從而導致煤體整體受力較大，而底部由于是固定邊界，因此會出現較高的應力集中現象。

2.2.3 破煤效果分析

為了獲得不同沖擊模式下水射流破煤效率，對比分析淹沒/非淹沒射流破碎煤體基坑深度隨時間的變化關系，統計了2 種射流破煤模式下每間隔20 μs 的侵蝕深度變化，非淹沒/淹沒射流破煤深度變化如圖6。

圖6 非淹沒/淹沒射流破煤深度變化Fig.6 Coal breaking depth variation of unsubmerged/submerged water jet

從圖6 可以看出：2 種沖擊模式下侵蝕深度都隨著沖擊時間的增加不斷增大，區別在于非淹沒射流的侵蝕深度基本與沖擊時間成正比，并且在沖擊計算結束之前穿透煤體。為了便于對比，按照規律將非淹沒射流的侵蝕深度延深至1 000 μs，此時射流速度為170、190、210 m/s 的侵蝕深度分別達到103、112、121 mm，對應沖擊時間時下淹沒射流的侵蝕深度分別為21、25、30 mm，分別為非淹沒射流侵蝕深度的20.4%、22.35%、24.8%。另外，淹沒條件下的侵蝕深度呈現階梯增加，并且增加的速度遠小于非淹沒條件下增加的速度。

根據淹沒條件下侵蝕深度的變化趨勢可以看出，隨著沖擊時間的延長，侵蝕深度的增加速率逐漸減小，這主要是隨著射流噴嘴距離煤體逐漸增大，射流到達煤體需要克服的行程阻力逐漸增加。可以推測隨著時間的增加射流沖擊煤體最終會停止，即一定速度下的射流通過能量不斷衰減達到煤體表面后的動能不足以破壞煤體，所形成的侵蝕坑深度就是極限破煤深度。

由于非淹沒條件下3 種射流速度的煤體均被穿透，雖然在圖6 中進行了預測補充，但卻不能對2 種射流模式下的破煤效果進行準確的分析，因此以210 m/s 穿透時的沖擊時間為上限，即對比沖擊發生后至830 μs 之間的侵蝕深度和煤體損傷單元，非淹沒/淹沒射流沖擊830 μs 時的破煤效果對比如圖7。

圖7 非淹沒/淹沒射流沖擊830 μs 時的破煤效果對比Fig.7 Comparison of coal breaking effect of unsubmerged/submerged water jet impacting 830 μs

從圖7 可以看出：在沖擊了830 μs 時，非淹沒條件下射流速度為170、190、210 m/s 時的侵蝕深度分別為83、91、100 mm，單元損傷數量分別為6 640、7 360、7 996 個；而對應速度下淹沒射流的侵蝕深度分別為16、20、25 mm，單元損傷數量分別為7 366、8 246、9 730 個，侵蝕深度分別為非淹沒射流的19.2%、21.9%、25%，單元損傷數量分別為非淹沒射流的110.9%、112%、121%。

對比2 種沖擊模式下煤體的侵蝕深度和損傷單元的數量發現，不同射流速度下的破煤深度與時間基本成正線性關系，根據2.1 節的分析可知，無論是淹沒環境還是非淹沒環境，水射流的傳播速度是確定的，但淹沒射流需要更多的能量克服介質水的阻力而使其傳播速度降低。另外，2 種沖擊模式下的損傷單元雖然隨射流速度的增加而增大，但是規律性較差，特別是淹沒環境下，損傷單元由于射流的湍流作用使得破煤面積大大增加，同時也導致了破煤深度大幅下降，但整體損傷單元的數量在增大。

2.3 淹沒/非淹沒射流破煤深度關系時變特征

為了研究淹沒射流的破煤特性，需要先獲得淹沒條件對水射流破煤的弱化作用。以非淹沒射流的侵蝕深度為基礎，求出對應沖擊時間下淹沒射流侵蝕深度和非淹沒射流侵蝕深度的關系，淹沒/非淹沒射流破煤深度時變關系如圖8。

從圖8 可以看出：在沖擊的早期淹沒射流的侵蝕坑深度與非淹沒射流的侵蝕坑深度的比值迅速增大，這是由于在淹沒射流發生沖擊破煤的時候非淹沒射流沖擊已經持續了一段時間，因此兩者的比值是從0 開始增加；之后，隨著非淹沒射流的侵蝕深度不斷增大，兩者比值也迅速增加，但3 種射流速度下的侵蝕深度變化規律一致，都是先迅速增大，之后趨于平穩。這里需要注意的是趨于平穩的是淹沒射流的侵蝕深度與非淹沒射流侵蝕深度的比值，即隨著沖擊時間的延長侵蝕深度都在增大，但兩者的比值卻趨于平穩。根據圖中擬合曲線可以推測淹沒射流的破煤深度最大約等于非淹沒射流的21%，因此根據這一結果可以結合已知非淹沒射流的破煤深度求得相應時間下的淹沒射流破煤深度。

根據研究結果可知，淹沒環境下水射流破煤效果較差，但在工程應用中通過延長沖刷時間或射流速度可以相應的提高破煤深度。另外，由于淹沒環境下水射流破煤效果難以考察，可以基于研究結果對非淹沒環境下的效果進行分析，初步評價淹沒射流的可行性。

3 結語

1）淹沒條件下水射流傳播過程中需要不斷消耗自身的能量來克服環境的阻力，因此需要水射流持續不斷地噴出才能保證射流向前傳播，否則射流會在環境水域的阻力下消耗消失。

2）相同條件下淹沒射流對煤體的損傷單元數略大于非淹沒射流，但侵蝕深度遠小于非淹沒射流，說明淹沒射流的能量在環境介質的阻礙下發散并使得破煤面積增大，同時導致破煤深度大幅下降。

3）相同沖擊時間下淹沒射流的破煤深度約等于非淹沒環境的21%，因此，在淹沒環境下開展破煤造穴時需要合理延長沖擊時間或者縮小鉆孔間距。