王 喆,梁 杰,侯騰飛,魏永超
(1.太原理工大學 安全與應急管理工程學院,山西 太原 030024;2.中國礦業大學(北京)化學與環境工程學院,北京 100083;3.河北省煤田地質局物測地質隊,河北 邢臺 054000)
煤炭地下氣化(UCG)就是將煤炭在原位進行有控制的燃燒,通過煤的熱解以及煤與氧氣、水蒸氣、二氧化碳發生一系列化學反應,產生一氧化碳、氫氣和甲烷等可燃氣體的原位流態化開采技術。該技術可以回收老礦井廢棄煤炭資源,對深部、急傾斜及傳統采煤技術難以開采的煤炭資源進行原位清潔轉化。被列為國家《能源技術革命創新行動計劃(2016—2030)》中煤炭無害化開采技術創新戰略方向。
相較于傳統井工開采,煤炭地下氣化最主要的特點就是氣化爐的高溫環境。1 000~1 200 ℃的高溫會使巖石產生熱膨脹進而產生熱應力,更為重要的是巖石的物理性質也會隨溫度發生改變。燃空區形成帶來的結構應力和高溫造成的熱應力共同作用對巖石造成損傷,損傷到達一定程度時,會使巖石中的原始裂縫擴展,并誘發新的裂縫生成,最終導致巖石破裂。進而造成頂板冒落,情況嚴重可能會導通含水層,造成氣化爐涌水,影響氣化過程。
近年來國內外學者對溫度-應力耦合過程和地下氣化過程中的圍巖穩定性進行了大量研究。席建奮等使用COMSOL軟件對煤炭地下氣化過程中頂板應力場變化過程及頂板穩定性進行實驗研究,表明熱應力的最大值可達1.5 MPa。李懷展等采用室內實驗、理論分析和數值模擬相結合的方法,研究了地下氣化過程中的溫度影響范圍。OTTO等建立了一種熱-力學耦合模型,以評估煤炭地下氣化氣化爐附近由于氣化和熱-力學效應造成的滲透性變化。辛林等研究煤炭地下氣化覆巖在熱固耦合條件下溫度、應力以及塑性區分布演化規律。JIAO等基于非連續性變形分析法,提出了一種模擬熱破裂過程的溫度-應力耦合模型,使用Mohr-Coulomb失穩準則來判定相鄰物體間是否有裂隙產生。但未對地下氣化過程中的巖石損傷特性進行研究。
徐小麗等根據基于不可逆熱力學理論和損傷理論,提出了熱-力耦合黏彈性損傷余能釋放率的理論表達式,建立了巖石熱-力耦合損傷破壞的能量準則。陸銀龍等根據彈性損傷理論建立了溫度-應力耦合作用下的巖石損傷本構方程,并對上覆巖層的拉伸損傷和裂隙發育進行了模擬。筆者先測得巖石不同溫度下的熱物性及力學基本參數,提出基于平滑Rankine損傷模型的高溫巖石損傷方程,使用COMSOL Multiphysics多物理場耦合軟件對深部地下氣化過程圍巖溫度、主應力、損傷變量進行模擬研究。
立足河北省大城勘查區詳查孔柱狀圖建立深部煤層地下氣化幾何模型,將地層結構進行整合確定模擬巖層及埋深,見表1。擬氣化的10號煤層為近水平煤層,厚度8 m左右,煤質為氣肥煤。使用COMSOL Multiphysics建立了尺寸為長()×寬()×高():490 m×80 m×282 m的模型,如圖1所示。從距左邊界90 m處向右進行后退式氣化,設有2條注氣管路,采寬為80 m。在燃空區工作面處設置2個8 m×8 m×8 m的熱源模擬氣化過程注氣管口的高溫。兩熱源的初始溫度為1 573.15 K,其余域初始溫度為273.15 K,所有外邊界設置為熱絕緣。采用施加載荷的方式補足模型上方未建模的1 100 m地層質量。地層的平均密度按2 600 kg/m計算,模型頂部應施加單位面積力為28 MN,因此模型頂部的邊界條件為邊界載荷,考慮模型自重,模型底部為固定約束。

圖1 模型試驗幾何示意

表1 巖層及埋深
對砂質泥巖、粉砂巖、細砂巖、泥巖、中砂巖5種典型的巖樣在不同溫度下的比熱容、導熱系數、彈性模量、單軸抗壓強度等基本物理力學參數進行測試。在COMSOL Multiphysics軟件中將材料的上述參數設置成隨溫度變化的函數。
使用耐馳LFA427型激光導熱儀對5種巖石不同溫度下的比熱容和導熱系數進行測試。使用加熱爐以100 ℃/h的速度加熱至預定溫度(常溫,200,400,600,800,1 000 ℃),加熱完成后保溫1 h,打開激光燈源進行測試。將測試結果進行擬合,擬合方程見表2。由表2可知:隨溫度升高,5種巖石導熱系數大體呈先下降后上升的趨勢,比熱容整體呈上升趨勢。以往研究也得出類似結果。在升溫初期,巖石中自由水、結合水的揮發是導致巖石導熱系數下降的主要原因;溫度繼續升高會使巖石發生熱膨脹,導致孔隙率增加,進而降低了導熱系數;溫度到達800 ℃之后,巖石中礦物質的結晶態發生變化,導致導熱系數出現升高。

表2 巖石比熱容、導熱系數隨溫度變化擬合方程
使用高溫電液伺服巖石力學實驗臺對不同溫度環境下巖石力學性能進行測試,先使用與伺服機配套的高溫加熱爐將樣品環境溫度升至預定值,在保持預定溫度20 min后,使用三軸伺服機進行加載實驗。將測試結果進行擬合,擬合方程見表3。

表3 巖石抗壓強度、彈性模量隨溫度變化擬合方程
從結果可知,5種巖石的抗壓強度和彈性模量隨溫度變化規律差別較大,但各自的抗壓強度和彈性模量隨溫度變化規律相似。
這是由于巖石孔隙度、礦物組成和結晶態不同導致不同巖石高溫力學性能隨溫度的變化規律存在差異。
煤層圍巖的熱物性及力學性質隨溫度發生變化,圍巖層中只發生熱傳導,且無內熱源,圍巖層的溫度源自燃空區與頂板交界處,本研究中假設頂板與熱源接觸處溫度與熱源相同,并以此開始進行熱傳導。因此,頂板的熱傳導方程為

(1)
式中,為巖石密度;,R為巖石比熱容;為巖石溫度;為時間;為巖石熱導系數。
在地質力學領域內,因為壓縮幾乎總是占主導地位,因此規定壓縮應力為正。當處理用于土壤和巖石的材料模型時,為了在軟件中保持一致性,還使用了“正張力”約定。
應力張量的不同不變量是建立本構模型的重要依據,也是解釋應力結果的重要依據。對于任何應力張量,3個基本不變量,,為
()=trace()
(2)

(3)
()=det()
(4)
其中,為應力。偏應力張量的不變量,,為
()=trace(dev())=0
(5)

(6)
()=det(dev())=

(7)

(8)
式中,為應力張量;為平均應力的影響;為剪應力的大小(≥0);包含剪應力的方向。
主應力是應力張量的本征值,由本征值方程計算得到。
-=0
(9)
3個主應力的順序為
≥≥
根據主應力,應力不變量為
()=++
(10)
()=++
(11)
()=
(12)
主應力是特征方程(Cayley-Hamilton定理)的根。

(13)
當對材料力學行為進行數學表述時,還需著重考慮材料的加載方式和環境條件。即使加載時間不變,黏彈性材料也具有時間依賴性。在本研究中,假設材料變形的黏性部分是不可壓縮的,因此體積變形純粹是彈性的。
假定頂板相似材料是均質各向同性的彈性體,依據彈性力學理論,巖體應力初始平衡方程為
+=0
(14)
其中,為應力二階張量;為方向上的受力。考慮溫度影響產生的熱膨脹和熱應力,此時用應變表示的應力方程為
=2+(++)+3
(15)

(16)

(17)

(18)
式中,為剪切模量;為總應變的二階形式;為拉梅常數;,,為,,方向上的應變;為線膨脹系數;為宏觀體積模量;為泊松比;為彈性模量。
在連續損傷力學理論中,損傷變量表示裂縫擴展引起的一系列屬性衰減。這個損傷變量控制了材料剛度的減弱,并且在應力和應變之間產生了非線性關系。對于線彈性材料,胡克定律將未損傷應力張量與彈性應變張量聯系起來:
=+:=+:(-)
(19)
其中,為四階彈性張量;“:”為雙點張量積;為彈性應變;為總應變;為所有非彈性應變;此外,還可能存在額外的應力貢獻,來自初始應力、外部應力或黏彈性應力。
對于標量損傷模型,由未損傷應力和損傷變量計算得到損傷應力張量
=(1-)
(20)
本研究中使用平滑Rankine損傷模型,使用3個未衰減主應力來定義溫度作用下的等效應變。

(21)
其中,為等效應變,為彈性應變的標量度量;符號“<>”為麥考利括號。在熱應力作用下,應力表達式為

(22)
式中,為,,方向上的應變之和。
定義為加載過程中的最大值,的演化遵循Kuhn-Tucker加載/卸載條件≤0,≥0,=0。
損傷模型的應變公式是基于加載函數。
=-≤0
(23)
使用線性應變軟化定律,得到溫度作用下的損傷變量的表達式。

(24)
()=0,<
(25)
其中,由抗拉強度和彈性模量計算,=表示損傷開始。參數由拉伸強度、特征元素尺寸、單位面積斷裂能或單位體積斷裂能等參數推導出。

(26)
溫度作用下,抗拉強度、彈性模量均為溫度的函數。
本研究中采用裂縫帶正則化法,以保持其網格的客觀性。裂紋帶法中的裂縫帶寬度是通過網格單元中的體積面積比來計算的。對于3D四面體網格,定義裂縫帶寬與四面體體積的關系為

(27)
氣化進行10,50,400 d時的溫度場分布如圖2所示。

圖2 不同氣化時間溫度場分布
氣化10 d時,溫度場的擴展范圍很小,影響的頂板范圍僅為3.27 m。氣化50 d時,溫度在頂板中的擴展范圍依舊較小,達到5.73 m。最終氣化400 d時,溫度在頂板中的影響范圍達到18.20 m。不同氣化時間熱源正上方覆巖溫度隨高度變化曲線如圖3所示,隨氣化進行,圍巖中的溫度下降梯度逐漸降低,從481.48 ℃/m降至71.43℃/m。溫度影響范圍隨氣化時間變化規律如圖4所示。對曲線進行擬合,擬合結果為指數函數,擬合效果良好。擬合結果表明,隨氣化進行溫度場在圍巖中的擴展速率在逐漸降低。

圖3 不同氣化時間圍巖溫度隨高度變化曲線

圖4 溫度影響范圍隨氣化時間變化擬合曲線
在實際地下氣化過程中,普遍采用控制注氣點后退氣化法,因此高溫區會隨連續油管不斷移動,巖石處于高溫區的時間在40 d左右。根據溫度場模擬結果可知,40 d內溫度場對圍巖的影響范圍約為4.7 m,高溫會對直接頂造成顯著影響。
對不同燃空區長度時,圍巖的損傷變量進行模擬,如圖5所示。損傷變量為“0”表示未出現損傷,損傷變量降低表明損傷加劇,損傷變量為“-1”則表示完全損傷。

圖5 不同燃空區長度時的損傷變量分布
燃空區上方及兩端均出現損傷區,且相互連接形成呈“凹”字型的大面積損傷區。燃空區長度90 m時,損傷區高度為72.2 m,損傷區高度與燃空區長度比為0.802。燃空區長度170 m時,損傷區高度為114.1 m,損傷區高度與燃空區長度比為0.673。燃空區長度250 m時,損傷區高度為148.8 m,損傷區高度與燃空區長度比為0.595。燃空區長度330 m時,損傷區高度為162.6 m,損傷區高度與燃空區長度比為0.493。
在氣化過程數值模擬中,圍巖層的初始溫度為常溫,因此當時間步為0時,圍巖未受高溫影響,不存在溫度-應力耦合,因此可以將400 d時的模擬的結果與0時的模擬結果進行對比,分析溫度對損傷變量的影響。
對不同燃空區長度時,直接頂損傷變量在0和400 d的模擬結果進行對比,結果如圖6所示。

圖6 直接頂損傷變量分布
通過0和400 d的模擬結果對比,可以得出高溫會降低直接頂在靠近熱源處的損傷情況。當燃空區長度為90 m時,溫度對損傷變量有著明顯影響,隨著燃空區長度的增長,溫度對損傷變量的影響逐漸降低,燃空區長度為330 m時,溫度依舊對損傷變量有一定的影響。這是因為損傷變量不同于熱應力,熱應力僅是溫度、熱膨脹系數和宏觀體積模量的函數,而損傷變量與巖石的抗折強度、彈性模量和等效應變等都有關。這說明就大城勘查區的地質條件下,溫度使煤層圍巖力學參數改變對頂板垮落的影響要明顯大于溫度產生的熱膨脹應力對頂板垮落的影響。
(1)隨溫度升高,砂質泥巖、粉砂巖、細砂巖、泥巖、中砂巖的比熱容整體呈上升趨勢,導熱系數整體呈下降趨勢。5種巖石的抗壓強度和彈性模量隨溫度變化規律差別較大,但各自的抗壓強度和彈性模量隨溫度變化規律相似。
(2)根據對圍巖溫度場的數值模擬,得到溫度影響范圍隨氣化時間呈指數變化。氣化10 d時,溫度影響范圍僅為3.27 m;氣化50 d時,溫度影響范圍達到5.73 m;氣化100 d時,溫度影響范圍為8.21 m;氣化400 d時,溫度影響范圍達到18.20 m。
(3)燃空區上方及兩端均出現損傷區,且相互連接形成呈“凹”字型的大面積損傷區。燃空區長度90 m時,損傷區高度為72.2 m;燃空區長度170 m時,損傷區高度為114.1 m;燃空區長度250 m時,損傷區高度為148.8 m;燃空區長度330 m時,損傷區高度為162.6 m。
(4)溫度對圍巖的損傷變量有著顯著影響,尤其在氣化初期溫度升高會明顯降低巖石的損傷程度,使其更加難以發生垮落。