高温对煤炭地下气化围岩损伤的影响

王 喆，梁 杰，侯腾飞，魏永超

(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 数值模型建立

1.1 模型建立

立足河北省大城勘查区详查孔柱状图建立深部煤层地下气化几何模型，将地层结构进行整合确定模拟岩层及埋深，见表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 岩层及埋深

1.2 不同温度下岩石的热物性及力学参数

对砂质泥岩、粉砂岩、细砂岩、泥岩、中砂岩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种岩石的抗压强度和弹性模量随温度变化规律差别较大，但各自的抗压强度和弹性模量随温度变化规律相似。

这是由于岩石孔隙度、矿物组成和结晶态不同导致不同岩石高温力学性能随温度的变化规律存在差异。

2 控制方程

2.1 温度控制方程

煤层围岩的热物性及力学性质随温度发生变化，围岩层中只发生热传导，且无内热源，围岩层的温度源自燃空区与顶板交界处，本研究中假设顶板与热源接触处温度与热源相同，并以此开始进行热传导。因此，顶板的热传导方程为

(1)

式中，为岩石密度；,R为岩石比热容；为岩石温度；为时间；为岩石热导系数。

2.2 应力控制方程

在地质力学领域内，因为压缩几乎总是占主导地位，因此规定压缩应力为正。当处理用于土壤和岩石的材料模型时，为了在软件中保持一致性，还使用了“正张力”约定。

应力张量的不同不变量是建立本构模型的重要依据，也是解释应力结果的重要依据。对于任何应力张量，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)

式中，为剪切模量；为总应变的二阶形式；为拉梅常数；，，为，，方向上的应变；为线膨胀系数；为宏观体积模量；为泊松比；为弹性模量。

2.2 损伤变量

在连续损伤力学理论中，损伤变量表示裂缝扩展引起的一系列属性衰减。这个损伤变量控制了材料刚度的减弱，并且在应力和应变之间产生了非线性关系。对于线弹性材料，胡克定律将未损伤应力张量与弹性应变张量联系起来：

=+:=+:(-)

(19)

其中，为四阶弹性张量；“：”为双点张量积；为弹性应变；为总应变；为所有非弹性应变；此外，还可能存在额外的应力贡献，来自初始应力、外部应力或黏弹性应力。

对于标量损伤模型，由未损伤应力和损伤变量计算得到损伤应力张量

=(1-)

(20)

本研究中使用平滑Rankine损伤模型，使用3个未衰减主应力来定义温度作用下的等效应变。

(21)

其中，为等效应变，为弹性应变的标量度量；符号“<>”为麦考利括号。在热应力作用下，应力表达式为

(22)

式中，为，，方向上的应变之和。

定义为加载过程中的最大值，的演化遵循Kuhn-Tucker加载/卸载条件≤0，≥0，=0。

损伤模型的应变公式是基于加载函数。

=-≤0

(23)

使用线性应变软化定律，得到温度作用下的损伤变量的表达式。

(24)

()=0，<

(25)

其中，由抗拉强度和弹性模量计算，=表示损伤开始。参数由拉伸强度、特征元素尺寸、单位面积断裂能或单位体积断裂能等参数推导出。

(26)

温度作用下，抗拉强度、弹性模量均为温度的函数。

本研究中采用裂缝带正则化法，以保持其网格的客观性。裂纹带法中的裂缝带宽度是通过网格单元中的体积面积比来计算的。对于3D四面体网格，定义裂缝带宽与四面体体积的关系为

(27)

3 结果与分析

3.1 温度场分布规律

气化进行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，高温会对直接顶造成显著影响。

3.2 损伤变量分布规律

对不同燃空区长度时，围岩的损伤变量进行模拟，如图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。

3.3 温度对损伤变量的影响

在气化过程数值模拟中，围岩层的初始温度为常温，因此当时间步为0时，围岩未受高温影响，不存在温度-应力耦合，因此可以将400 d时的模拟的结果与0时的模拟结果进行对比，分析温度对损伤变量的影响。

对不同燃空区长度时，直接顶损伤变量在0和400 d的模拟结果进行对比，结果如图6所示。

图6 直接顶损伤变量分布

通过0和400 d的模拟结果对比，可以得出高温会降低直接顶在靠近热源处的损伤情况。当燃空区长度为90 m时，温度对损伤变量有着明显影响，随着燃空区长度的增长，温度对损伤变量的影响逐渐降低，燃空区长度为330 m时，温度依旧对损伤变量有一定的影响。这是因为损伤变量不同于热应力，热应力仅是温度、热膨胀系数和宏观体积模量的函数，而损伤变量与岩石的抗折强度、弹性模量和等效应变等都有关。这说明就大城勘查区的地质条件下，温度使煤层围岩力学参数改变对顶板垮落的影响要明显大于温度产生的热膨胀应力对顶板垮落的影响。

4 结 论

(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)温度对围岩的损伤变量有着显著影响，尤其在气化初期温度升高会明显降低岩石的损伤程度，使其更加难以发生垮落。