高温后无烟煤孔隙演化及多孔弹塑性变形特征研究

2022-08-16周长冰万志军

中国矿业 2022年8期

周长冰，陈鑫，张源，万志军

(1.绍兴文理学院土木工程学院，浙江绍兴 312000；2.中国矿业大学矿业工程学院，江苏徐州 221116；3.深部煤炭资源开采教育部重点实验室，江苏徐州 221116)

0 引言

煤是一种对温度、压力非常敏感的有机岩石。在煤的地下气化、煤层原位注热开采这一系列采矿工程中，煤体均处于热力耦合的环境中[1-2]。煤的热解导致煤体孔隙发育，造成煤层承载能力减弱，并在地层压力下发生了变形，而研究煤体高温下体孔隙演化及变形特征有助于解决上述工程问题。

研究煤热解导致煤体孔隙发育等过程时，获取煤体物理力学参数是一项必要的工作，对煤体力学参数的温度效应研究主要包含温度对煤体内孔隙结构的影响，王勇等[3]、周长冰等[4-5]对煤样开展热解实验，再利用CT微焦点扫描、压汞实验等手段探究其细观结构，分析热解温度对煤体的细观孔隙裂隙结构的影响；EVERSON等[6]构建了肥煤燃烧的孔隙随机扩散模型，认为初始总反应速率受颗粒内扩散控制，且随温度的升高影响增大，温度和压力对煤体变形发生的难易程度和变形量大小都有影响[7-11]；冯子军等[12]对无烟煤试样从20 ℃升至600 ℃过程中的变形规律开展研究，指出热力耦合作用和热解产气是影响煤体变形的关键因素；孟召平等[13]通过煤岩力学试验研究了煤岩物理力学性质和煤岩全应力-应变过程中的渗透规律，发现在全应力-应变过程中应变软化现象不明显或者具有应变硬化现象的煤样，最大渗透率主要发生峰值前的塑性变形阶段；赵洪宝等[14]以肥煤制作的煤试件为研究对象，对其在一定围压和瓦斯压力下的三轴压缩力学特性、不同起始应力点卸轴压时的力学特性进行了试验研究，研究表明随着起始应力点的升高，煤样变形的非线性呈增加趋势，卸轴压过程中的平均弹性模量呈下降趋势，瓦斯压力对煤样变形的影响呈增加趋势。此外，煤岩的承载能力和抗变形能力都随温度的加剧而逐步劣化[15-16]；万志军等[17]通过三轴试验对比分析了无烟煤和气煤弹性模量变化规律的异同，讨论了弹性模量突变临界温度和热解产气对煤体弹性模量的影响；FENG等[18]研究了在高温状态下煤岩体弹性模量及变形演化特征，研究发现在三轴压缩条件下，当加热到400 ℃时，弹性模量呈现非线性变化。

上述研究成果具有重要意义，揭示了煤热解后物理力学参数变化的本质。但是在表征煤体变形能力时，均以随温度升高弹性模量大幅度弱化为指标，不考虑表现为孔隙压缩的塑性阶段。针对这一问题，本文从多孔弹塑性角度对无烟煤热解后变形破坏规律开展研究。对不同温度下热解处理的无烟煤开展了压汞实验及单轴压缩实验，分析了孔隙结构的演化特征及应力-应变曲线，通过理论分析建立了多孔弹塑性模型，以期揭示煤炭热解后多孔弹塑性变形特征，研究成果对流体化采矿岩层控制具有一定的借鉴意义。

1 实验概况

1.1 实验设备

1) 绝氧加热装置。煤炭热解处理必须处于无氧环境，为此研制了绝氧加热装置，如图1所示。装置主要由三个部分组成：真空泵、氮气瓶及绝氧容器，各个设备之间通过阀门及PU软管连接，其中绝氧容器用以安置试样。

图1 试样绝氧安装装置Fig.1 Sample oxygen isolating installation device

2) 加载设备采用MTS815.03电液伺服岩石实验系统(图2)，由加载框架、增压系统、控制系统、计算机系统构成，最大轴压为4 600 kN，机架刚度为10.5×109N/m，全程计算机控制，可实现自动数据采集及处理。

图2 MTS单轴试验机Fig.2 MTS uniaxial testing machine

1.2 试样

试样采自晋城矿区的无烟煤，如图3所示。先在采煤工作面截取大块煤样，再由专用石材加工机加工成圆柱形毛坯，研磨精加工，达到试验规程要求。成品规格尺寸为Φ50 mm×100 mm，煤样整体性良好，表面平整且有轻微裂隙。

图3 无烟煤试样Fig.3 Anthracite specimens

1.3 实验方案

加热方案：马弗炉，外接控温器控制；处理温度为常温下25 ℃、100 ℃、200 ℃、300 ℃、400 ℃、500 ℃，每个温度点有3个试样。升温速率设定为15 ℃/h，达到目标温度后保温4 h以确保煤样内温度平衡及充分的热解反应，然后关闭电源，待煤样在炉腔内自然冷却至室温后取出。

加载方案：采用MTS高精度伺服试验机，加载方式为位移控制，加载速率为0.1 mm/min。

1.4 实验步骤

1) 将试样分组标号，测量并记录试样尺寸、重量并保留部分煤块用以开展常温下的压汞实验分析。

2) 将试样及煤块颗粒放入绝氧容器内，并用细砂充填试样周围空隙，避免试样直接接触容器从而导致加热不均匀。对绝氧容器进行密封处理，从预留的接口处对容器抽真空，再充入1个标准大气压的氮气，如此反复数次，最大限度的抽出容器内氧气，使其内部孔隙充满1个标准大气压的氮气，为无烟煤试样提供绝氧加热环境。

3) 将装有试样及煤块颗粒的绝氧容器分组放入马弗炉进行绝氧加热，达到目标温度后保温4 h，自然冷却至室温后打开容器取出试样，测量试样尺寸及重量。

4) 对加热过后的无烟煤颗粒开展压汞实验，分析其孔隙结构。

5) 将温度处理后的试样安装在试验机上，采用轴向位移的控制方式，以0.1 mm/min的加载速率对试样进行加载至试样破坏，得到应力-应变曲线。重复操作，测定各组煤样温后压缩实验的应力应变数据。

2 实验结果

2.1 无烟煤孔隙结构随温度变化规律

根据霍多特方案[19]将煤样的孔隙结构划分为微孔(直径d<0.01 μm)、小孔(0.01 μm100 μm)。试样在不同温度下孔隙率、热解孔隙结构及热重损失率如图4所示。

图4 不同温度下各孔型孔隙体积比、总孔隙率及热重损失率Fig.4 Pore volume ratio,total porosity and thermogravimetric loss rate ofeach pore type at different temperatures

由图4(a)可知，无烟煤总孔隙率随温度升高而升高，300 ℃之前孔隙率随温度上升平缓增加，300 ℃之后热解加剧温度升高孔隙率骤然增加。无烟煤热解过程中微孔、小孔的孔隙体积占比呈平缓下降趋势，而绝对孔隙率表现为先上升后下降的特点，到500 ℃时出现大幅度下降(表1)。

表1 各孔型孔隙绝对孔隙率Table 1 Absolute porosity of each pore type

中孔的孔隙体积占比有所增长，变化范围不大，但是绝对孔隙率增长较快；随温度的升高，大孔的体积比和孔隙率都有大幅度的升高。从热重数据可知(图4(b))，100 ℃时，煤质量损失约1%，从100 ℃升温到300 ℃时，质量仅损失约0.2%，说明在300 ℃前，主要为脱水、气体解吸等引起的物理失重。 300 ℃后，试样质量损失急剧升高，取出试样时有臭气溢出，说明试样发生了热解反应，反应阈值在300～400 ℃之间，这与冯子军等[12]得出的结论相似。对于变化原因，文献[3]和文献[16]已给出详细分析，此处不再赘述。

2.2 无烟煤热解后单轴压缩的破坏特征

对不同温度后的试样开展单轴压缩实验，每个温度点试样数量为3个，选取峰值强度为中值的应力-应变曲线，如图5所示。

图5 不同加热温度后应力-应变曲线Fig.5 Stress-strain curves after differentheating temperatures

由图5可知，压缩前期试样中原有张开性结构面或微裂隙逐渐闭合，曲线呈上凹型，且随孔隙率的升高压密段范围也有扩大的趋势，到400 ℃时，压密段应变范围占破坏前总应变的50%；同时，由室温及100 ℃的应力-应变曲线中可以看出，应力-应变曲线末端有明显的近直线段，可以认为此处为纯弹性变形。随着温度的升高，曲线末端的直线特征逐渐减弱，可以认为弹性特征逐渐弱化，直至500 ℃时，从曲线上已经很难区分出压密段和弹性段的分界点，由此可知，孔隙率升高到一定程度后，试样的压缩变形包含了孔缩压缩及弹性变形的共同作用。

从破坏的角度分析，认为应力-应变曲线可分为两种类型：第一种处于为常温至300 ℃之间，可以看到煤样具有明显的脆性特征，应力达到峰值后急剧跌落，且破坏形式为整体的大块度破坏，如图6(a)所示；第二种处于300 ℃之后，试样在达到峰值应力之后，抗压能力逐渐减弱，存在明显的塑性软化现象，韧性特征逐渐增强，且破坏后碎块多且块度较小，如图6(b)所示。

图6 破坏状态对比Fig.6 Destruction status comparison

对于多孔介质，一般认为其承载能力由骨架提供，孔隙率升高导致试样内骨架的有效面积减小，在相同的载荷下骨架内局部应力升高，即试样在压缩时表现为总体应变增加、强度降低。由此可以解释峰值强度与孔隙率高度关联的现象(图7)。

图7 孔隙率与峰值强度随温度变化规律Fig.7 Porosity and peak intensity vary with temperature

3 考虑孔隙压缩的弹塑性应变计算方法

3.1 骨架弹性模量的计算方法

以近似直线段计算不同温度下宏观弹性模量的方法可以表征煤体的宏观变形能力，具有计算简单、应用方便等优点，这种方法认为煤体在压缩过程中存在完全弹性阶段。应用此方法可以方便地计算出各温度点的宏观弹性模量(图8)。

图8 弹性模量随温度升高的变化规律Fig.8 Variation law of elastic modulus withincreasing temperature

对于多孔介质，试样的变形和承载能力主要由骨架提供。矿体在原位开采时，部分矿物热解形成孔隙，未热解部分形成骨架，其成分保持不变，骨架结构为致密状态，假定骨架的力学参数为定值，不随孔隙的发育而改变，骨架的有效作用面积计算公式见式(1)。

A′=A0(1-P0)

(1)

式中：A′为骨架的有效面积，mm2；A0为试样横截面积，mm2；P0为常温下试样面孔隙率。

因此，骨架的弹性模量可由式(2)表示。

(2)

式中：εe为轴向平均弹性应变差；E0为试样的宏观弹性模量，GPa，由应力-应变曲线的近直线段计算E0=F/(A0εe)。

进一步简化式(2)为式(3)。

(3)

由于面孔隙率难以获取，而体孔隙率可以通过压汞实验等手段测出，在计算多孔介质相关参数时通常使用的是体孔隙率代替面孔隙率[20]。需要说明的是，式(3)中以体孔隙率代替面孔隙率存在一定误差，与文献[4]中灰度处理后面孔隙率与体孔隙率相互应证(图9)，面孔隙率26.89%，压汞实验测得体孔隙率28.12%，误差不大，但是替换后极大地简化了计算过程，方便应用。

图9 灰度处理测算面孔隙率Fig.9 Measuring face gap rate by gray processing

将图8中各温度点宏观弹性模量及对应的体孔隙率代入式(3)，计算出各孔隙率下的骨架弹性模量，结果见表2。由表2可知，随着孔隙率的增加，骨架弹性模量急剧减小，这与“骨架结构为致密状态，力学参数为定值，不随孔隙的发育而改变”相互矛盾，分析认为，计算宏观弹性模量时将所有的应变都作为完全弹性应变εe，夸大的弹性应变导致了以近似直线段方法计算的弹性模量偏小，这也进一步印证了前文分析的结果，说明煤体的变形是弹性、塑性同时存在的。

表2 各孔隙率下骨架弹性模量与宏观弹性模量Table 2 Elastic modulus and macroscopic elastic modulusof skeleton under each porosity

分析可知，在孔隙发育后，试样的压缩变形包含了弹性变形及孔隙压缩变形，孔隙越发育，孔隙压缩的特征越明显。试样的变形主要为骨架变形，试样孔隙高度发育后，骨架有效面积减小，加载时骨架内部应力远大于试样宏观应力[21]，表现为局部骨架断裂导致的空隙压缩，弹性变形与骨架坍缩不可必免同时发生。

常温时，由于试样孔隙率较小，试样整体致密，从应力-应变曲线可以看出，加载初期压密段明显，并快速过渡到近直线段，可以认为此近直线段为弹性段，用以作为试样整体的宏观弹性模量计算区域。而孔隙高度发育后，试样变形曲线近直线段内的应变不是完全弹性，包含有孔隙压缩的应变，因此，不能用于计算骨架变形。由此可知，骨架弹性模量可取值为2.56 GPa。

3.2 宏观真弹性模量计算方法

由于试样中弹性、孔隙压缩变形同时存在，为了与传统方法计算的试样宏观弹性模量进行区分，将完全弹性变形部分对应的弹性模量定义为宏观真弹性模量。对于孔隙发育后的试样，已知变形包括完全弹性变形εe及孔隙压缩变形εp，对于弹性变形部分，其骨架弹性模量与宏观真实弹性模量Ee依然满足式(3)，简化后可用式(4)表示。

(4)

式中：Ee为孔隙发育后多孔介质压缩变形中的真实弹性模量；Pe为孔隙发育后对应的面孔隙率。

由于骨架弹性模量为固有属性，联合式(3)和式(4)可得式(5)。

(5)

已知常温下宏观弹性模量E0、面孔隙率P0及孔隙发育后的面孔隙率Pe，由式(5)可以计算出孔隙发育后煤体真实弹性模量。

3.3 多孔煤体的真弹性变形及孔隙压缩变形

由于压缩过程中存在孔隙压缩变形，在常温下孔隙率较小，认为孔隙压缩变形影响较小，以常温状态下宏观弹性模量E0=2.46 GPa为标准值计算，带入式(5)，得到各孔隙率下的真弹性模量，并与传统算法的宏观弹性模量对比(表3)。

表3 不同孔隙率下真弹性模量Table 3 True elastic modulus under different porosity

由表3可知，真弹性模量也随孔隙率的升高而下降，但下降幅度远小于传统方法的计算结果，因此，孔隙发育后其完全弹性变形也占总变形很小的一部分，大部分为孔隙压缩变形，只要计算出完全弹性应变就能分离出孔隙压缩应变量。

完全弹性应变量计算公式见式(6)。

(6)

孔隙压缩应变量表达式见式(7)。

εp=ε-εe

(7)

结合实验数据，分离出破坏前不同孔隙率的孔隙压缩变形，图10为不同温度的孔隙压缩应变规律，由图10可知，不同温度对应着不同的孔隙率，曲线包含压缩初期的孔隙闭合段，应力增加至3 MPa后，曲线非常陡峭，其斜率接近无穷大，说明孔隙闭合完成，可认为试样为完全弹性变形。孔隙发育后，随孔隙率的增加，曲线斜率不断减小，而压缩过程中斜率又是不断增加的，说明在孔隙压缩存在于整个压缩变形过程中，随着孔隙不断压实，压缩速率不断减小，弹性变形所占有的比重将越来越大，且存在硬化现象。破坏前的最大孔隙压缩量及其比重见表4。

图10 不同温度的孔隙压缩应变规律Fig.10 Law of pore compressive strain atdifferent temperatures

表4 孔隙压缩应变相关参数Table 4 Relevant parameters of pore compressive strain

4 塑性强化模型的建立及数值计算验证

4.1 模型建立

从分离出的孔隙压缩应变可知，在压缩过程中，破坏前其变形具有塑性硬化的特征，总体呈现指数增长。孔隙率越大，硬化能力越弱，将孔隙压缩看成塑性变形，参考Voce硬化模型并对其进行改造，建立了适合孔隙压缩的塑性模型，表达式见式(8)。

(8)

式中：σ0为初始屈服强度，MPa；k为强化系数，MPa；β为饱和指数，β=1；εp为塑性应变，εp=1；εpmax为最大孔隙压缩应变。

结合图10确定硬化模型参数，结果见表5。

表5 塑性硬化模型参数Table 5 Plastic strengthening model parameters

以孔隙发育后的300～500 ℃温度下试样参数为例建立数值模型，取对应的真弹性模量、强化系数、饱和指数参数。由于硬化模型包含压缩初期的孔隙闭合段，可认为在加载初时立即发生塑性变形，设置初始塑性强度σ0=0.01 MPa。采用轴对称模型，高度为10 cm，半径为2.5 cm，下部边界固定，上部边界施加位移载荷，位移速度为0.1 mm/s。

4.2 计算结果

图11为总应变为500 ℃试样1.95%时试样内的弹性、塑性应变分布。由图11可知，因端面效应的影响，试样两端的应变较小，中部大部分范围比较均匀，平均塑性应变为1.65%，弹性应变为0.3%，与实验数据吻合较好。图12为数值计算与实验的应力-应变曲线对比。由图12可知，总应变达到最大孔隙压缩应变后，应力不再升高，与实验结果吻合较好。