基于ALE算法的連續(xù)高壓水射流破巖數(shù)值模擬

馬文濤

(1.中煤科工開采研究院有限公司，北京 100013；2.天地科技股份有限公司 開采設(shè)計事業(yè)部，北京 100013)

近年來，高壓水射流技術(shù)得到了長足的發(fā)展[1，2]，高壓水射流是通過增壓泵將水加至高壓狀態(tài)，然后通過噴嘴流出的高能水束，其高效、節(jié)能，以及在鉆孔成型上具有的優(yōu)點使該技術(shù)深受青睞，具有廣闊的應(yīng)用前景[3，4]，在煤炭領(lǐng)域，利用高壓水射流開展煤層割縫，進(jìn)行煤層氣開采[5-7]。此外，在煤層中割縫也被作為一種解危措施用來進(jìn)行沖擊地壓防治[8，9]。但射流的理論研究仍然停留在實驗研究基礎(chǔ)之上，這是因為射流的問題復(fù)雜，是一個涉及多因素的非線性動力學(xué)問題，具有破壞瞬時性、大變形，高應(yīng)變率等特點，研究難度很大。隨著計算機(jī)技術(shù)的發(fā)展，通過數(shù)值模擬方法對相關(guān)問題進(jìn)行可重復(fù)性研究，可進(jìn)行破巖過程中任意時刻的狀態(tài)分析[10-12]。汪紅祥[13]模擬分析了射流沖擊速度、入射角、以及移動速度對巖石損傷的影響。劉佳亮[14]對比無圍壓狀態(tài)下巖石和高圍壓狀態(tài)下巖石在高壓水射流作用下破碎坑演化情況，認(rèn)為處于高圍壓狀態(tài)下的巖石損傷沿軸向的演化速率明顯低于無圍壓狀態(tài)下的巖石，沿徑向的損傷演化受圍壓影響較小。盧義玉[15]利用SPH方法模擬了射流在破巖過程中應(yīng)力波形成、傳播及衰減過程，并分析了射流速度、巖石性質(zhì)對應(yīng)力波效應(yīng)的影響。孫清德[16]利用非線性有限元法研究了高壓水射流破碎巖石規(guī)律，并提出了射流破碎巖石過程分為射流沖蝕破碎及錘沖破碎。以上成果使得射流破巖機(jī)理與特征研究進(jìn)一步深入，取得了一定的研究進(jìn)展，但上述研究通常將射流簡化為一段水柱，與真實連續(xù)射流有所差異。本文采用ALE算法以實現(xiàn)射流連續(xù)性，并建立高壓射流破巖模型，使模擬更接近實際情況，計算結(jié)果更具有參考價值，對射流破巖機(jī)理、巖石破壞特征及巖石損傷演化進(jìn)行研究分析，以期提高對高壓水射流破巖過程的認(rèn)識。

1 數(shù)值計算模型

高壓水射流破巖是一種涉及流固耦合的非線性動力學(xué)問題，在數(shù)值模擬計算過程作出如下假設(shè)：①高壓水射流可以被看作為不可壓縮流體；②射流在空氣中無擴(kuò)散效應(yīng)；③巖石簡化為各向同性。

1.1 力學(xué)模型

模型的控制方程可描述為[17]：

連續(xù)性方程：

動量方程：

能量方程：

變形率：

式中，ρ為介質(zhì)密度，kg/m3；θ為介質(zhì)速度，m/s；E為介質(zhì)內(nèi)能密度，J/m3；Bi為體力，N；ωi為介質(zhì)相對與網(wǎng)格速度，m/s；σij為kron函數(shù)；t為時間，s。

1.2 接觸準(zhǔn)則

采用ALE算法集成了拉格朗日和歐拉有限元的優(yōu)勢，一方面能夠有效捕捉介質(zhì)運(yùn)動，另一方面使得單元網(wǎng)格相互獨立于實體，因此在該算法下，當(dāng)高速射流與巖石相互作用時，由ALE算法生成的射流源為主體，由拉格朗日算法生成的巖石為從體，并采用罰函數(shù)以實現(xiàn)流固耦合。

1.3 本構(gòu)模型及材料參數(shù)

1.3.1 本構(gòu)模型

通常巖石本構(gòu)模型有TCK模型、Mazars模型和H-J-C模型，高壓水射流作用巖石時，接觸巖石表現(xiàn)為大變形以及非線性等特點，因此本文選用合理描述該類問題的H-J-C本構(gòu)模型，H-J-C本構(gòu)模型的等效應(yīng)力表達(dá)式描述如下：

σ*=[A(1-D)+BP*N](1+CInε*)

損傷因子由等效塑性應(yīng)變和塑性體積應(yīng)變累積計算，當(dāng)損傷值D為1時表示完全破壞，其表達(dá)式如下：

由于NULL模型無剪切剛度，因此NULL材料模型適用于水、空氣等流體模型的建立，以Gruneisen狀態(tài)方程表示射流狀態(tài)方程為[18]：

+(γ0+αμ)E

式中，P為壓力，Pa；ρ0為流體密度，kg/m3；C為沖擊波速度與質(zhì)點速度曲線截距，m/s；γ0為Grun系數(shù)；E為體積內(nèi)能密度，J/m3；α為一階體積修正系數(shù)；S1、S2、S3為常數(shù)。

1.3.2 材料參數(shù)

本模型中設(shè)置的材料參數(shù)見表1。

表1 材料參數(shù)

1.4 幾何模型及邊界條件

建立射流源、空氣、巖石物理模型，由于模擬連續(xù)射流，因此將射流源模型直徑設(shè)置為2mm，空氣域模型尺寸為60mm×60mm×55mm，巖石模型尺寸為50mm×50mm×50mm。確定模型尺寸后，對物理模型采用映射網(wǎng)格法進(jìn)行網(wǎng)格劃分，巖石采用拉格朗日網(wǎng)格，水和空氣采用ALE網(wǎng)格。

依據(jù)巖石真實工況設(shè)置模型邊界條件，對四周施加無反射邊界條件，對巖石底部施加全自由度約束。射流與巖石之間采用點對面的侵徹接觸算法，并施加以拉應(yīng)力和第一主應(yīng)力控制的實效準(zhǔn)則，在求解計算過程中能夠及時刪除實效的單元。巖石選用彈塑性常應(yīng)力實體單元，根據(jù)圣維南原理在射流直徑10倍以上的巖石沖擊區(qū)域受影響很小，所以在射流沖擊區(qū)域設(shè)置了較細(xì)密的網(wǎng)格。

2 射流破巖機(jī)理分析

利用上述模擬方法對高壓水射流破巖過程進(jìn)行分析，射流速度為200m/s。高壓水射流破巖機(jī)理較為復(fù)雜，相關(guān)學(xué)者進(jìn)行了一定的研究，射流破巖主要有兩種方式，其中一種作用方式為在水錘壓力以及滯止壓力下，裂隙、裂紋擴(kuò)展、貫穿，造成巖石宏觀破壞，另一種作用方式為在射流沖擊應(yīng)力波作用下，巖石發(fā)生拉伸破壞。但是這兩種作用方式哪種作用方式為主導(dǎo)因素，仍然需要深入探討。

根據(jù)實驗過程中巖石剝落情況分析破巖機(jī)理，巖石片狀剝落情況可以看出射流下方巖石發(fā)生片狀剝落，并且在射流與巖石接觸位置產(chǎn)生初始基坑。射流破巖機(jī)理[19，20]可以描述為：射流初期在沖擊壓力與水錘壓力作用下，巖石發(fā)生初始壓破壞，而在射流外圍產(chǎn)生拉、剪應(yīng)力集中，而抗拉、抗剪強(qiáng)度較低，因此在射流四周產(chǎn)生拉剪破壞，并形成拉剪裂隙，最后在水楔作用下，水射流穿透及滲入巖石，在拉、剪裂隙內(nèi)部產(chǎn)生應(yīng)力集中，達(dá)到破巖效果，使得巖石在初期發(fā)生片狀剝落，而在后期，由于深度增加，巖石破碎程度增加，并無巖石片狀剝落情況發(fā)生。初始壓破壞-拉剪破壞-水楔作用周而復(fù)始，使得破巖過程不斷進(jìn)行，巖石沖蝕深度不斷增加，高壓水射流破巖模型如圖1所示。

圖1 高壓水射流破巖模型

3 巖石破壞形態(tài)及損傷分析

3.1 巖石破壞形態(tài)分析

不同時刻巖石破壞形態(tài)如圖2所示，破碎坑最終成明顯的漏斗狀。本文分別截取了60μs、100μs、150μs、200μs下，巖石的破壞特征進(jìn)行分析。60μs時，高速高壓水柱與巖石接觸，巖石表面受到巨大的壓應(yīng)力，急劇壓縮，瞬時破壞形態(tài)呈矮圓柱形。100μs時，水柱外圍速度大大降低，但是水柱中心仍然保持的一定高的速度，繼續(xù)沖蝕巖石，瞬時破壞形態(tài)為矮“V”形。150μs時，回流水體作用加強(qiáng)，使水柱中心速度大大降低，不具有破壞作用，此刻，沖蝕深度不再增加，瞬時形態(tài)呈高“V”形。200μs時，縫槽深部寬度增加變化速度不明顯，最終形態(tài)呈漏斗狀，并與實驗結(jié)果相似，實驗破巖形態(tài)剖面如圖3所示。

圖2 不同時刻巖石破壞形態(tài)

圖3 實驗破巖形態(tài)剖面

高壓水射流破巖具有階段性，不同射流階段巖石破壞形態(tài)具有不同的特點。此外，在高壓水射流作用下，巖石在很短時間內(nèi)就發(fā)生破壞，巖石破壞一定深度、一定寬度作為最終結(jié)果，之后隨時間變化甚小。

3.2 巖石損傷分析

高壓水射流作用下巖石損傷的演化規(guī)律直觀反應(yīng)著巖石的破壞情況，本部分分別分析單元損傷演化規(guī)律以及模型損傷演化規(guī)律，進(jìn)一步了解在高壓水射流作用下巖石細(xì)觀破壞規(guī)律。

3.2.1 單元損傷分析

隨著遠(yuǎn)離射流中心，選取了三個有代表性的單元分析其損傷破壞情況，由近至遠(yuǎn)單元編號分別為H54380#單元、H54376#單元以及H54373#單元，其中H54380#單元位于射流源正下方，單元選取情況(1/4模型)如圖4所示。

圖4 單元選取情況圖(1/4模型)

圖5 距射流中心不同距離單元損傷情況

距射流中心不同距離單元損傷情況如圖5所示，由圖5可以看出，H54380#單元、H54376#單元損傷值D均達(dá)到了1，由理論部分可知H54380#單元、H54376#兩單元發(fā)生了破壞。破壞單元損傷情況具有三個階段：無損傷階段、密集損傷階段以及破壞階段，在無損傷階段，射流在空氣中未與巖石發(fā)生接觸，單元不會損傷。隨著射流接觸巖石表面，單元逐漸損傷并進(jìn)行積累。當(dāng)單元損傷積累至1時，單元發(fā)生破壞，進(jìn)入破壞階段。對比H54380#號單元可以看出，由于H54380#單元位于射流正下方，應(yīng)力集中程度較大，單元在密集損傷階段斜率較大，由初始損傷至破壞較快，2μs內(nèi)就進(jìn)入破壞階段。H54376#單元對比H54380#單元較遠(yuǎn)，在應(yīng)力波傳遞以及衰減作用下，單元在密集損傷階段斜率變小，由初始損傷至破壞時間增加。H54373#單元損傷值為達(dá)到1，單元未發(fā)生破壞，這是由于3#單元距離射流中心更遠(yuǎn)，應(yīng)力波作用下單元發(fā)生損傷而不會發(fā)生破壞。這是由于射流能量高度集中，其破巖具有明顯的局部特征。

3.2.2 模型損傷分析

由單元損傷分析可以得出，依據(jù)距離射流中心的遠(yuǎn)近，巖石損傷還具有分區(qū)特征，巖石破壞損傷分布情況如圖6所示，即由近至遠(yuǎn)巖石破壞損傷不同，可以分為巖石破壞區(qū)、損傷區(qū)以及原始狀態(tài)區(qū)，在原始狀態(tài)區(qū)，射流不會對巖石產(chǎn)生損傷作用，即損傷值為0。選取不同距離的單元得出其最大損傷值，并繪制曲線，巖石破壞損傷分布情況如圖所示，小于破巖臨界點距離為破壞區(qū)，大于損傷臨界點距離為原始狀態(tài)區(qū)，介于兩者區(qū)域為損傷區(qū)，該條件下巖石損傷半徑約為破壞半徑的2.3倍。

圖6 巖石破壞損傷分布情況

不同時刻巖石損傷演化情況如圖7所示，可以看出，在射流與巖石作用初期，巖石損傷分布較為集中，隨著射流破巖繼續(xù)進(jìn)行，應(yīng)力波在豎直及水平方向傳播，巖石損傷范圍在水平方向以及豎直方向都有不同程度的增加，在射流初期由于射流速度垂直巖石表面，因此巖石在豎直上的損傷程度大于水平方向，破碎坑豎直方向擴(kuò)展，表現(xiàn)為深度增加。在射流后期，深度增加到一定程度時，射流受到回流水體的嚴(yán)重影響，射流豎直速度降低，水平方向速度增加，因此巖石在水平方向的損傷程度大于豎直方向，破碎坑豎直方向擴(kuò)展，表現(xiàn)為直徑的增加。

圖7 不同時刻巖石損傷演化情況

4 結(jié) 論

1)射流破巖過程中巖石受到壓應(yīng)力、拉應(yīng)力以及剪應(yīng)力，巖石破碎是由沖擊壓力、水錘壓力、應(yīng)力波、水楔共同作用的結(jié)果，射流破巖具有階段循環(huán)性，在不同破巖循環(huán)時期階段各作用力具有相應(yīng)的主導(dǎo)性。

2)巖石破壞形態(tài)具有時變性，隨著時間推移，其破壞形態(tài)不斷變化，最終破壞形態(tài)呈漏斗狀或深“V”形，巖石破壞具有局部特征以及分區(qū)性，可分為破壞區(qū)、損傷區(qū)以及原始狀態(tài)區(qū)。

3)巖石損傷具有階段性，可分為無損傷階段、密集損傷階段以及破壞階段。