CN102567649B - 一种煤层地下燃空区的体积建模方法 - Google Patents

Abstract

本发明公开一种煤层地下燃空区的体积建模方法，从热传导基本理论出发，基于燃空区的高温与常规地下温度场的相对恒温，通过分析渗流对温度场的影响机理，推导出渗流作用下的温度场控制方程；引入贝塞尔函数，推导出燃烧通道周围温度场分布方程，并给出有限元数值求解，借助现场温度探测资料，用有限元数值软件反演渗流作用下的温度场分布，并确定燃烧通道的高温影响边界；通过调整燃空区的截面尺寸，获得与现场探测拟合性最佳的煤层温度分布，反演出煤炭地下燃烧燃空区形状及尺寸，进而计算出燃空区的近似体积。此方法能够较为准确地判定煤炭地下燃烧通道的形状及尺寸，比较其他计算手段，所提出的燃烧通道的燃空区热源模型具有一定的实用价值。

Description

一种煤层地下燃空区的体积建模方法

技术领域

本发明属于地质工程领域，特别涉及一种确定煤层燃烧后地层空洞区域的形状及尺寸，进而确定燃空区体积的方法。

背景技术

煤层地下燃烧是世界煤炭开采的研究方向之一，是从根本上解决传统开采方法存在的一系列技术和环境问题的重要途径。在实际运行中，由于燃空区埋藏较深以及地质条件的复杂性，仅依靠现有的地质测量手段很难对煤层地下燃空区的形状与大小做出准确计算。因此，目前还没有一种有效的方法来计算地下燃空区的体积，以此来解决在燃烧采煤过程中诸如如何进一步提高燃烧热值、有效控制燃烧、防止地面塌陷(即燃空区的稳定性)以及解决塌陷导致的大量地下水流失等技术问题。

另外，渗流的热学效应分析是水工、采矿等工程中不可忽略的一部分。渗流的裂隙岩体地区，运动的地下水与其所接触的岩体之间会发生热量交换(即表面对流换热)，使得岩体温度降低，温度场分布形式发生变化。另外，渗透系数以及渗流速度的变化，引起了岩体与水流交换的热量变化，同样影响了岩体温度场的分布。在燃烧采煤区，由于受高温热源的影响，煤层内的温度分布出现异常；地下水的渗流必将带走一部分热量，使得燃空区附近热值降低，同时改变了煤层内的温度分布状态。热源体积大小、强度的不同导致煤层内温度分布状态不同。

假设岩土体均质且各向同性，基于深度地层的渗流效应影响较小，因此忽略燃烧通道周围渗流场对温度场的影响。因此，通过对地下温度场的研究进而来反推燃烧区体积、评价燃烧区的安全稳定性等重大问题对燃烧采煤技术的长远发展至关重要，本发明即是基于此种考虑，提出了一种地下燃空区的建模方法，本案由此产生。

发明内容

本发明的目的，在于针对现有的地质测量手段很难对煤层地下燃空区的形状与大小做出准确计算，提供一种煤层地下燃空区的体积建模方法，其能够较为准确地判定煤炭地下燃烧通道的形状及尺寸，比较其他计算手段，所提出的燃烧通道的燃空区热源模型具有一定的实用价值。

为了达成上述目的，本发明的解决方案是：

一种煤层地下燃空区的体积建模方法，包括如下步骤：

(1)根据岩体内部单元体热量平衡原理，在岩层平面内建立渗流作用下的一维非稳定温度场方程：

$\mathrm{c\ρ}\frac{\∂T}{\∂t}=-c_w{\mathrm{\ρ}}_w\frac{\∂\left(\mathrm{vT}\right)}{\∂x}+\frac{\∂}{\∂x}\left(\mathrm{\λ}\frac{\∂T}{\∂x}\right)$

假设岩体内部有一椭圆形热源，且该热源温度处处相等，渗流影响下的二维非稳定温度场控制方程如下式：

$\mathrm{c\ρ}\frac{\∂T}{\∂t}=q_0-c_w{\mathrm{\ρ}}_w[\frac{\∂\left[v_x\left(H\right)T\right]}{\∂x}+\frac{\∂\left[v_y\left(H\right)T\right]}{\∂x}]+\mathrm{\λ}(\frac{{\∂}^{2}T}{{\∂x}^2}+\frac{{\∂}^{2}T}{{\∂y}^2});$

其中：T-温度，℃；q₀-恒定热源单位面积发热功率，w/m²；c_w-水的比热，J/(kg·℃)；ρ_w-水的密度，kg/m³；c-岩石的比热，J/(kg·℃)；ρ-岩石密度，kg/m³；λ-岩石导热系数，J/(m·s·℃)；v-地下水在岩体中渗透流速；v_x、v_y-地下水在岩体中分别沿x、y方向的渗透流速，m/s；K-渗透系数，m/s；H-水头差；

对于稳定流来讲满足：

$\frac{{\∂}^{2}H}{{\∂x}^2}+\frac{{\∂}^{2}H}{{\∂y}^2}=0$

将 代入到上式得：

$\mathrm{\λ}(\frac{{\∂}^{2}T}{{\∂x}^2}+\frac{{\∂}^{2}T}{{\∂y}^2})-c_w{\mathrm{\ρ}}_w[v_x{T}_x+v_y{T}_y]+q_0=\mathrm{c\ρ}\frac{\∂T}{\∂t};$

(2)以形函数N_k为权函数，其中k＝1，2，L，将恒定热源量放于定解条件中，在二维平面内对渗流影响下的温度场方程用Galerkin加权余量法得到：

$\underset{s}{\∫\∫}{N}_k\{\mathrm{\λ}(\frac{{\∂}^{2}T}{{\∂x}^2}+\frac{{\∂}^{2}T}{{\∂y}^2})-c_w{\mathrm{\ρ}}_w(v_x{T}_x+v_y{T}_y)-\mathrm{c\ρ}\frac{\∂T}{\∂t}\}\mathrm{dxdy}+\underset{\mathrm{\Γ}}{\∫}{N}_k\frac{\∂T}{\∂n}\mathrm{d\Γ}=0$

上式中的第一项为在求解域内的面积分，第二项是满足边界条件的线积分。对于无需满足第三类边界条件的情况，此时认为材料表面放热系数为零；

(3)在时间域内，用向后差分法进行二维离散，得到渗流场影响下的温度场有限元求解方程：

$(L+P){\left\{T\right\}}_t+\frac{\left[G\right]}{\mathrm{\Δt}}({\left\{T\right\}}_t+{\left\{T\right\}}_{t+\mathrm{\Δt}})=0$

其中： $L_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{l}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\left\{\underset{s}{\∫\∫}\right[\frac{{\∂H}_i}{\∂x}\·\frac{{\∂H}_j}{\∂x}\·\frac{{\∂H}_i}{\∂y}\·\frac{{\∂H}_j}{\∂y}\left]\right|J\left|\mathrm{d\ξd\η}\right\};$

$p_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{p}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\frac{{\mathrm{Kc}}_w{\mathrm{\ρ}}_w}{\mathrm{\λ}}\underset{s}{\∫\∫}{N}_i[\frac{{\∂H}_i}{\∂x}\·\frac{{\∂H}_j}{\∂x}\·\frac{{\∂H}_i}{\∂y}\·\frac{{\∂H}_j}{\∂y}]\left|J\right|\mathrm{dxdy};$

$G_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{g}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\frac{\mathrm{c\ρ}}{\mathrm{\λ}}\underset{s}{\∫\∫}{N}_i{N}_j\left|J\right|\mathrm{dxdy};$

$\text{|J|=}\left[\begin{array}{cc}\frac{\∂x}{\∂\mathrm{\ξ}}& \frac{\∂x}{\∂\mathrm{\η}}\\ \frac{\∂y}{\∂\mathrm{\ξ}}& \frac{\∂y}{\∂\mathrm{\η}}\end{array}\right]=\left[\begin{array}{c}\mathrm{\Σ}\frac{{\∂N}_i}{\∂\mathrm{\ξ}}{x}_i\mathrm{\Σ}\frac{{\∂N}_i}{\∂\mathrm{\ξ}}{y}_i\\ \mathrm{\Σ}\frac{\∂N_i}{\∂\mathrm{\η}}{x}_i\mathrm{\Σ}\frac{\∂N_i}{\∂\mathrm{\η}}{y}_i\end{array}\right],$ 是雅克比矩阵的行列式。

采用上述方案后，本发明具有以下有益效果：

(1)本发明对燃烧采煤过程中煤炭地下燃烧产生的燃空区形状及尺寸探测具有较高的实用价值；

(2)本发明可以实现煤炭地下燃烧中燃空区体积的计算，有利于解决燃空区上部岩层的稳定性问题，并对如何避免煤层顶板大区域脱落导致的地面塌陷问题具有一定的理论指导意义；

(3)本发明不仅能计算出燃烧后燃空区的体积，也能计算煤炭在燃烧过程中孔洞生长的体积。

附图说明

图1是本发明所采用的热源模型简化示意图；

图2是剖面内温度场的数值模拟示意图；

图3是本发明中的计算曲线和探测曲线；

图4是短轴b＝3m时煤层在渗流影响下的温度场分布图；

图5是短轴b＝4m时煤层在渗流影响下的温度场分布图；

图6是短轴b＝5m时煤层在渗流影响下的温度场分布图；

图7是本发明的流程图。

具体实施方式

以下将结合附图，对本发明的技术方案进行详细说明。

如图7所示，本发明提供一种煤层地下燃空区的体积建模方法，包括如下步骤：

(1)根据岩体内部单元体热量平衡原理，在岩层平面内建立渗流作用下的一维非稳定温度场方程：

$\mathrm{c\ρ}\frac{\∂T}{\∂t}=-c_w{\mathrm{\ρ}}_w\frac{\∂\left(v\right)T}{\∂x}+\frac{\∂}{\∂x}\left(\mathrm{\λ}\frac{\∂T}{\∂x}\right)$

假设岩体内部有一椭圆形热源，且该热源温度处处相等，其示意图如图1。渗流影响下的二维非稳定温度场控制方程如下式：

$\mathrm{c\ρ}\frac{\∂T}{\∂t}=q_0-c_w{\mathrm{\ρ}}_w[\frac{\∂\left[v_x\left(H\right)T\right]}{\∂x}+\frac{\∂\left[v_y\left(H\right)T\right]}{\∂x}]+\mathrm{\λ}(\frac{{\∂}^{2}T}{{\∂x}^2}+\frac{{\∂}^{2}T}{{\∂y}^2})$

其中：T-温度，℃；q₀-恒定热源单位面积发热功率，w/m²；c_w-水的比热，J/(kg·℃)；ρ_w-水的密度，kg/m³；c-岩石的比热，J/(kg·℃)；ρ-岩石密度，kg/m³；λ-岩石导热系数，J/(m·s·℃)；v-地下水在岩体中渗透流速；v_x、v_y-地下水在岩体中分别沿x、y方向的渗透流速，m/s；K-渗透系数，m/s；H-水头差。

对于稳定流来讲满足：

$\frac{{\∂}^{2}H}{{\∂x}^2}+\frac{{\∂}^{2}H}{{\∂y}^2}=0$

将 代入到上式得：

$\mathrm{\λ}(\frac{{\∂}^{2}T}{{\∂x}^2}+\frac{{\∂}^{2}T}{{\∂y}^2})-c_w{\mathrm{\ρ}}_w[v_x{T}_x+v_y{T}_y]+q_0=\mathrm{c\ρ}\frac{\∂T}{\∂t}$

(2)以形函数Nk为权函数，其中k＝1，2，L，将恒定热源量放于定解条件中，在二维平面内对渗流影响下的温度场方程用Galerkin加权余量法得到：

$\underset{s}{\∫\∫}{N}_k\{\mathrm{\λ}(\frac{{\∂}^{2}T}{{\∂x}^2}+\frac{{\∂}^{2}T}{{\∂y}^2})-c_w{\mathrm{\ρ}}_w(v_x{T}_x+v_y{T}_y)-\mathrm{c\ρ}\frac{\∂T}{\∂t}\}\mathrm{dxdy}-\underset{\mathrm{\Γ}}{\∫}{N}_k\frac{\∂T}{\∂n}\mathrm{d\Γ}=0$

上式中的第一项为在求解域内的面积分，第二项是满足边界条件的线积分。对于无需满足第三类边界条件的情况，此时可认为材料表面放热系数为零。

(3)在时间域内，用向后差分法进行二维离散，可得到渗流场影响下的温度场有限元求解方程：

$(L+P){\left\{T\right\}}_t+\frac{\left[G\right]}{\mathrm{\Δt}}({\left\{T\right\}}_t+{\left\{T\right\}}_{t+\mathrm{\Δt}})=0$

其中： $L_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{l}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\left\{\underset{s}{\∫\∫}\right[\frac{{\∂H}_i}{\∂x}\·\frac{{\∂H}_j}{\∂x}\·\frac{{\∂H}_i}{\∂y}\·\frac{{\∂H}_j}{\∂y}\left]\right|J\left|\mathrm{d\ξd\η}\right\};$

$p_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{p}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\frac{{\mathrm{Kc}}_w{\mathrm{\ρ}}_w}{\mathrm{\λ}}\underset{s}{\∫\∫}{N}_i[\frac{{\∂H}_i}{\∂x}\·\frac{{\∂H}_j}{\∂x}\·\frac{{\∂H}_i}{\∂y}\·\frac{{\∂H}_j}{\∂y}]\left|J\right|\mathrm{dxdy};$

$G_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{g}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\frac{\mathrm{c\ρ}}{\mathrm{\λ}}\underset{s}{\∫\∫}{N}_i{N}_j\left|J\right|\mathrm{dxdy};$

$\text{|J|=}\left[\begin{array}{cc}\frac{\∂x}{\∂\mathrm{\ξ}}& \frac{\∂x}{\∂\mathrm{\η}}\\ \frac{\∂y}{\∂\mathrm{\ξ}}& \frac{\∂y}{\∂\mathrm{\η}}\end{array}\right]=\left[\begin{array}{c}\mathrm{\Σ}\frac{{\∂N}_i}{\∂\mathrm{\ξ}}{x}_i\mathrm{\Σ}\frac{{\∂N}_i}{\∂\mathrm{\ξ}}{y}_i\\ \mathrm{\Σ}\frac{\∂N_i}{\∂\mathrm{\η}}{x}_i\mathrm{\Σ}\frac{\∂N_i}{\∂\mathrm{\η}}{y}_i\end{array}\right],$ 是雅克比矩阵的行列式。

以下将结合具体实施例对本发明进一步说明。

假设无越流或地下水的垂向补给，即垂向由于渗流而引起的热量损失可以忽略，只考虑岩石介质的热传导作用。由于煤层厚度相对于地层尺寸较小，在模拟剖面温度分布时，燃空区高度的改变对燃空区上部地层的温度分布影响不大，可先假定燃空区高度为一定值。建立剖面二维热源模型，通过对比不同剖面形状，并将现场探测温度与模拟温度分布比较，可得到较合理的剖面形状以及温度分布如图2。

通过对工程实例的计算和分析得到长半轴a＝18m，即燃空区长度约36m，计算温度值与实测温度值较为接近，如图3。

假定煤层为均质、等厚、各向同性的无限体，根据燃空区的剖面形状建立煤层平面内等效热源模型。在探测区域除观测孔外，其余大部分钻孔孔底都位于煤层内，因此，在煤层平面内，多次试探，利用各孔煤层中的温度大小，可以插值得到煤层中的温度场分布。在仅考虑热传导的情况下，以探测到的孔底温度为热源，对探测区域内的温度场进行有限元数值模拟，并通过对比不同短轴大小的计算结果，得到的探测区域内煤层温度分布如图4、图5、图6。

当短半轴b＝3m时，赋予热源平均温度值95℃，利用comsol有限元数值软件模拟得到渗流影响下的煤层温度分布，不难发现图中的高温影响范围较实际情况偏小，且比较现场探测的温度值也偏低；当短半径b＝4m时，同样赋予热源温度值95℃，模拟得到的渗流影响下的温度分布，与现场探测温度值较为接近；当短半轴b＝5m时，模拟得到的温度分布，高温热源影响范围较大，各钻孔底部温度分布较实际情况偏高。因此可以断定，燃空区短半轴约为4m。

综上，燃空区长、短半轴分别为a＝18m，b＝4m，近似认为燃空区的宽度和高度c相等，因此燃空区体积应用联通试验计算得到的燃空区体积约为1500m³，两种方法的计算结果较为接近。

综上，本发明一种煤层地下燃空区的体积建模方法，依据国内外研究现状的同时，基于以往方法的可操作性以及现场地质条件的复杂多变性，从热传导基本理论出发，基于燃空区的高温与常规地下温度场的相对恒温，通过分析渗流对温度场的影响机理，推导出渗流作用下的温度场控制方程。引入贝塞尔函数，推导出燃烧通道周围温度场分布方程，并给出了有限元数值求解，借助现场温度探测资料，用有限元数值软件反演渗流作用下的温度场分布，并确定出燃烧通道的高温影响边界。实际运用过程中，通过调整燃空区的截面尺寸，获得与现场探测拟合性最佳的煤层温度分布，可以反演出煤炭地下燃烧燃空区形状及尺寸，进而计算出燃空区的近似体积。该方法能够较为准确地判定煤炭地下燃烧通道的形状及尺寸，比较其他计算手段，所提出的燃烧通道的燃空区热源模型具有一定的实用价值。

以上实施例仅为说明本发明的技术思想，不能以此限定本发明的保护范围，凡是按照本发明提出的技术思想，在技术方案基础上所做的任何改动，均落入本发明保护范围之内。

Claims (1)

1.一种煤层地下燃空区的体积建模方法，其特征在于包括如下步骤：

(1)根据岩体内部单元体热量平衡原理，在岩层平面内建立渗流作用下的一维非稳定温度场方程：

$\mathrm{c\ρ}\frac{\∂T}{\∂t}=-c_w{\mathrm{\ρ}}_w\frac{\∂\left(\mathrm{vT}\right)}{\∂x}+\frac{\∂}{\∂x}\left(\mathrm{\λ}\frac{\∂T}{\∂x}\right)$

假设岩体内部有一椭圆形热源，且该热源温度处处相等，渗流影响下的二维非稳定温度场控制方程如下式：

$\mathrm{c\ρ}\frac{\∂T}{\∂t}=q_0-c_w{\mathrm{\ρ}}_w[\frac{\∂\left[v_x\left(H\right)T\right]}{\∂x}+\frac{\∂\left[v_y\left(H\right)T\right]}{\∂x}]+\mathrm{\λ}(\frac{{\∂}^{2}T}{\∂x^2}+\frac{{\∂}^{2}T}{\∂y^2});$

其中：T—温度，℃；q₀—恒定热源单位面积发热功率，w/m²；c_w—水的比热，J/(kg·℃)；ρ_w—水的密度，kg/m³；c—岩石的比热，J/(kg·℃)；ρ—岩石密度，kg/m³；λ—岩石导热系数，J/(m·s·℃)；v—地下水在岩体中渗透流速；v_x、v_y—地下水在岩体中分别沿x、y方向的渗透流速，m/s；K—渗透系数，m/s；H—水头差；

对于稳定流来讲满足：

$\frac{{\∂}^{2}H}{\∂x^2}+\frac{{\∂}^{2}H}{\∂y^2}=0$

将 $v_x=K\frac{\∂H}{\∂x},v_y=K\frac{\∂H}{\∂y},$ 代入到上式得：

$\mathrm{\λ}(\frac{{\∂}^{2}T}{\∂x^2}+\frac{{\∂}^{2}T}{\∂y^2})-c_w{\mathrm{\ρ}}_w[v_x{T}_x+v_y{T}_y]+q_0=\mathrm{c\ρ}\frac{\∂T}{\∂t};$

(2)以形函数N_k为权函数，其中k＝1,2,…，将恒定热源量放于定解条件中，在二维平面内对渗流影响下的温度场方程用Galerkin加权余量法得到：

$\underset{s}{\∫\∫}{N}_k\{\mathrm{\λ}(\frac{{\∂}^{2}T}{\∂x^2}+\frac{{\∂}^{2}T}{\∂y^2})-c_w{\mathrm{\ρ}}_w(v_x{T}_x+v_y{T}_y)-\mathrm{c\ρ}\frac{\∂T}{\∂t}\}\mathrm{dxdy}-\underset{\mathrm{\Γ}}{\∫}{N}_k\frac{\∂T}{\∂n}\mathrm{d\Γ}=0$

上式中的第一项为在求解域内的面积分，第二项是满足边界条件的线积分，对于无需满足第三类边界条件的情况，此时认为材料表面放热系数为零；

(3)在时间域内，用向后差分法进行二维离散，得到渗流场影响下的温度场有限元求解方程：

$(L+P){\left\{T\right\}}_t+\frac{\left[G\right]}{\mathrm{\Δt}}({\left\{T\right\}}_t+{\left\{T\right\}}_{t+\mathrm{\Δt}})=0$

其中： $L_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{l}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\left\{\underset{s}{\∫\∫}\right[\frac{\∂N_i}{\∂x}\·\frac{\∂N_j}{\∂x}+\frac{\∂N_i}{\∂y}\·\frac{\∂N_j}{\∂y}\left]\right|J\left|\mathrm{d\ξd\η}\right\};$

$P_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{p}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\frac{K{c}_w{\mathrm{\ρ}}_w}{\mathrm{\λ}}\underset{s}{\∫\∫}{N}_i[\frac{\∂N_i}{\∂x}\·\frac{\∂N_j}{\∂x}+\frac{\∂N_i}{\∂y}\·\frac{\∂N_j}{\∂y}]\left|J\right|\mathrm{dxdy};$

$G_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}{g}_{\mathrm{ij}}=\underset{e}{\mathrm{\Σ}}\frac{\mathrm{c\ρ}}{\mathrm{\λ}}\underset{s}{\∫\∫}{N}_i{N}_j\left|J\right|\mathrm{dxdy};$

$\left|J\right|=\left|\begin{array}{cc}\frac{\∂x}{\∂\mathrm{\ξ}}& \frac{\∂x}{\∂\mathrm{\η}}\\ \frac{\∂y}{\∂\mathrm{\ξ}}& \frac{\∂y}{\∂\mathrm{\η}}\end{array}\right|=\left|\begin{array}{c}\mathrm{\Σ}\frac{\∂N_i}{\∂\mathrm{\ξ}}{x}_i\mathrm{\Σ}\frac{\∂N_i}{\∂\mathrm{\ξ}}{y}_i\\ \mathrm{\Σ}\frac{\∂N_i}{\∂\mathrm{\η}}{x}_i\mathrm{\Σ}\frac{\∂N_i}{\∂\mathrm{\η}}{y}_i\end{array}\right|,$ 是雅克比矩阵的行列式。