祁南矿东风井冻结温度场时空演化规律分析

黄诗清，荣传新，龙 伟，马昊辰，何骏珍

(安徽理工大学 土木建筑学院，安徽 淮南 232001)

人工冻结法[1-5]已广泛应用于地铁、隧道及煤矿等施工，该工法能够适应十分复杂的工程条件，能够大大提升煤矿井筒施工的安全性。随着煤矿冻结凿井深度的不断增加，地层条件也逐渐复杂，因此，对于不同地层不同土性冻结效果差异的研究必不可少。学者们从诸多方面对冻结法展开了研究分析，在数值计算方面，文献[6-11]以各个矿井为实际工程背景，基于相关地质参数，通过有限元计算软件建立冻结温度场数值计算模型，对多圈管冻结温度场的发展特性和冻结温度场的影响因素展开研究。在理论研究方面，刘波等[12]推导了竖向直排三管冻结壁温度分布计算公式，并推广到直排和多排冻结壁温度场分布计算，进而研究了斜井冻结壁温度场的发展规律；张涛等[13]通过分析冻结管内盐水流动的特点和状态，同时基于相似理论，得到了冻结管内盐水状态对温度场的影响规律；盛天宝等[14]通过对某矿多圈孔冻结壁温度场实测，分析了冻结前期存在浅部片帮、掘进速度慢以及冻结孔冷量浪费的原因。部分学者采用室内试验和数值计算相结合的形式，陈军浩等[15]利用多圈管冻结模型试验并结合数值模拟对冻结管有无偏斜2 种情况的冻结壁温度场的发展特性进行了对比分析；任建喜等[16]通过物理力学实验、现场实测以及数值模拟相结合的方法，分析了冻结法凿井期间冻结壁的受力机制和冻结壁温度场的分布规律；周盛全等[17]通过现场实测、室内试验以及数值模拟分析了冻结地层的热物理参数分布规律、冻结壁冻融规律及其力学特性；Wang Zhi 等[18]通过模型实验的方法研究了钙质黏土层冻结温度场的发展规律；Wang Bin 等[19]采用数值模拟对冻结温度场的发展情况进行预测，结果表明在采用局部差分冻结技术后，内排冻结孔形成的冻结壁向内扩展范围得到限制，井壁降温速率明显降低；Hu Xiangdong 等[20]研究了FSPR 的冻结特性，同时探索了冻土墙在冻结和开挖时期随时间发展的变化规律；Yao Zhishu 等[21]基于实际工程，采用人工冻结法对白垩纪地层凿井温度场进行数值模拟和试验分析。

前人对冻结温度场的研究已具有一定的广度和深度，但大部分仅从单一层位土体且不考虑冻结孔实际成孔位置偏斜来对冻结温度场的发展特性展开研究，而在实际工程中，冻结孔偏斜是十分复杂且不可忽略的问题。陈红蕾等[22]针对深冻结井筒温度场成孔弱界面，应用数值模拟软件分析了深冻结井筒冻结孔在3 个层位上实际成孔位置下的温度场弱界面参数；汪仁和等[23]在考虑了土层中水的相变潜热以及冻结温度随冻结时间变化的前提下，分析了冻结管偏斜和不偏斜2 种情况下的冻结壁发展特性；焦华喆等[24]通过分析600 m 深地层的钻孔偏斜、水的相变潜热以及不同温度下的导热系数、比热容等影响因素，对地层冻结温度场的发展情况进行了预测。笔者以安徽祁南矿东风井冻结法凿井为工程背景，利用大量实测数据从多方面对祁南矿表土层冻结温度场特性展开分析，同时考虑冻结孔实际成孔位置和不同埋深土体之间初始地温、比热容及导热系数的差异，研究细砂、钙质黏土及砂质黏土3 种不同埋深土体冻结温度场的分布与发展规律，基于建立的数值计算模型对冻结温度场进行预测，该研究旨在为相关冻结工程的设计施工提供参考依据。

1 工程概况

祁南煤矿位于安徽省宿州市埇桥区祁县镇境内，北距宿州市约23 km，南距蚌埠市约70 km。其东风井位于矿井南侧，距离工业广场约4.5 km，井筒地面高程+22.4 m，井口高程+24.0 m，井筒采用冻结法施工，井筒净直径6 m，冻结深度405 m，临时锁口7 m，冻结段井筒掘砌长度395 m，基岩段掘砌深度48.7 m(含井底水窝0.7 m)，马头门两侧各5.0 m。

冻结管相关设计参数见表1，冻结孔、测温孔以及水文孔的布置如图1 所示，本次计算主面路径(图1)2 种，一种计算路径同时通过主排孔和辅助孔，称主面1；另一种计算路径同时通过主排孔和防片帮孔，称主面2。冻结孔设计布置剖面如图2 所示，其中包括4 个测温孔，C1、C2 孔深405 m，C3 孔深390 m，C4 孔深280 m；3 个水文孔，S1 孔深36 m，S2 孔深130 m，S3 孔深230 m，单号主排孔共22 个，孔深385 m，双号主排孔共22 个，孔深405 m。

图1 冻结孔布置与主面路径Fig.1 Freezing holes layout and main surface paths

图2 冻结孔剖面图Fig.2 Section diagram of the frozen hole

表1 冻结管主要技术参数Table 1 Main technical parameters of frozen pipe

2 冻结壁温度场数值计算模型

2.1 基本理论

立井冻结温度场是一个具有相变、移动边界、内热源、边界条件复杂的不稳定三维导热问题。冻结壁横向尺寸远小于纵向尺寸，且冻结壁纵向的热传导较微弱。因此，在分析冻结壁温度场时，可以将三维冻结温度场简化为二维平面冻结温度场问题。由热物理学和冻土学理论，得出立井冻结温度场的控制微分方程[25]为：

式中：θn为冻结温度场中任意一点的温度，℃；t为冻结时间，d；r为冻结区域内任意一点到井筒圆心的距离，m；an为导温系数，m2/s，an=λn/(ρnCn)，λn、ρn和Cn分别为导热系数、密度和比热容，W/(m·K)、kg/m3、J/(kg·K)；n为 土体的状态，n=1为未冻土，n=2为冻土。

在冻结开始前，地层温度的初始条件为：

式中：θ0为土层的初始温度，℃。

在冻结过程中，冻结管壁与周围土层满足Dirichlet 边界条件，其表达式为：

式中：(xP,yP)为冻结管管壁的坐标；θc(t)为冻结管内的盐水温度，℃。

距离冻结区域无穷远处满足Dirichlet 边界条件，其表达式为：

2.2 土体热物理参数

依据室内土体热工试验，选择埋深218 m 钙质黏土层、埋深225 m 细砂层以及埋深259 m 砂质黏土层作为不同土性、不同埋深的研究层位，得到不同层位土体热物理参数(表2)。

表2 土体热物理参数Table 2 Soil thermal physical parameters

2.3 初值与边界条件

根据冻结前，现场实测的地层初始温度可知，钙质黏土层位、细砂层位以及砂质黏土层位的初始温度分别为21.50、21.63 和22.56℃。冻结管的边界条件取现场实测的盐水温度，如图3 所示。

图3 盐水温度趋势Fig.3 Brine temperature trend

2.4 数值计算模型的建立

祁南煤矿东风井的几何模型取半径为40 m，并假设为均质且各向均匀的土体，忽略地下水渗流对温度场的影响，忽略井筒地层以及冻结管的竖向传热，将冻结温度场简化为二维平面问题。

基于冻结孔的实际成孔位置来建立冻结温度场有限元计算模型，有限元模型采用的是三节点三角形的二维实体热单元来进行网格划分，对井筒周围土体的网格划分加密，使计算结果更准确，在远离冻结管的区域，由于温度梯度变化较小，单元梯度变化稀疏，剖面单元剖分不加密，数值模型一共划分为11 374 个单元，其中包括340 个顶点单元，518 个边界单元，最小单元质量0.482 6，温度场数值模型网格划分如图4 所示。

图4 温度场数值模型网格划分Fig.4 Meshing of temperature field numerical model

3 冻结壁温度场时空演化规律

3.1 温度场实测数据分析

根据190～300 m 埋深段内冻结145 d 中4 个测温孔的温度变化，绘制如图5 所示的实测温度三维图。C1 和C2 测温孔距离冻结管位置最远，降温速度最慢；细砂层位相比于钙质黏土层位以及砂质黏土层位具有更快的下降趋势。C3 测温孔位于主排孔和辅助孔之间，在2 排孔冷量叠加的影响下，降温速度较快，冻结40 d 左右，测温孔温度到达0℃左右，土体内的水开始结冰同时释放相变潜热，与冻结管释放的冷量相互抵消，因此，该时间段内出现短暂的温度基本不变的相变平台。C4 测温孔距离辅助孔最近，在多圈冻结孔冷量叠加的影响下，降温速度最快，冻结75 d 防片帮孔停冻，因此，该测温孔降温速度变缓，此后测温孔温度在一定范围内产生波动，其原因是在开挖到相应层位时，外井壁浇筑混凝土的同时产生混凝土水化热使测温孔温度短暂上升，随着混凝土内水化反应的减缓，其所释放的热量也在减少，冻结管内的温度又会重新下降。

图5 4 个测温孔实测数据变化三维图Fig.5 Three dimentional diagram of measured data of four temperature measuring holes

3.2 数值计算模型验证

图6 为砂质黏土层位的4 个测温孔实测结果与对应位置的数值计算结果对比图，4 个测温孔模拟与实测误差均在±2.0℃以内，模拟与实测结果具有较高的一致性，因此，通过数值计算对各个层位冻结温度场的发展情况进行预测是可行的。

图6 埋深259 m 砂质黏土层位实测结果与数值计算结果对比曲线Fig.6 Comparison curves between measured results and numerical calculation results of sandy clay with buried depth of 259 m

3.3 有效冻结壁平均温度和厚度

有效冻结壁指扣除入荒径冻土后所剩下的那部分冻结壁，且考虑冻结壁厚度的不均匀性，取不同方位的有效厚度，计算其平均值作为有效冻结壁厚度。冻结壁平均温度取冻结面上冻结壁温度在其面积上的加权平均值。

如图7 所示，冻结壁有效平均厚度与冻结时间呈正相关关系，即冻结时间越长，冻结壁厚度越大。在冻结期内，钙质黏土层位冻结壁厚度以0.010 5 m/d 的速度增长，细砂层位冻结壁厚度以0.012 8 m/d 的速度增长，砂质黏土层位冻结壁厚度以0.011 5 m/d 的速度增长，细砂层位的冻结壁发展速度要快于钙质黏土和砂质黏土，钙质黏土和砂质黏土的冻结壁发展速率大致相同。在冻结壁交圈的初期阶段其厚度的增加速度明显，原因是交圈初期冻结壁厚度由冻结管圈径以外以及冻结管圈径以内同时向两侧扩展，当冻结壁扩展至开挖荒径以后，外侧冻土与周围土层接触，使其冻结壁发展速度变缓。在相同冻结时间和相同冷量情况下，细砂层位所形成的冻结壁厚度最大，砂质黏土形成的冻结壁厚度低于细砂，钙质黏土所形成的冻结壁厚度最小。细砂、研质黏土及钙质黏土层位所形成的冻结壁最终厚度分别为7.15、6.88 和6.63 m，均达到冻结壁设计指标值6.2 m，满足设计需求。

图7 冻结壁有效厚度及平均温度与时间关系Fig.7 Relationship between effective frozen wall thickness and average temperature and time

冻结壁的平均温度随着冻结时间的变化划分为2个阶段，在冻结80～130 d 内，冻结壁平均温度随冻结时间快速下降，钙质黏土层位、细砂层位、砂质黏土层位的冻结壁平均温度下降速度分别为−0.032 50、−0.042 38、−0.039 25℃/d；从冻结速度来看，细砂层位的冻结速度要快于砂质黏土层位，钙质黏土层位冻结效果最差；从最终冻结效果来看，细砂层位冻结壁最终平均温度最低，砂质黏土层位次之，钙质黏土层位最高。细砂、砂质黏土及钙质黏土层位所形成的冻结壁在开挖时平均温度分别为−19.86、−18.74 和−18.22℃，均低于开挖时冻结壁平均温度设计指标值−15℃，满足设计需求。

产生以上2 种现象的主要原因是由于土体本身性质决定的，细砂的导热系数大于砂质黏土和钙质黏土，导热系数越大，冷量传递就越快，因此，冻结壁平均温度就下降得越快，同时冻结壁厚度增长得越快。

3.4 主面温度场时空变化

图8—图10 为不同层位冻结壁2 个主面温度分布情况。将冻结壁温度场沿径向划分为A、B、C 三个区域，在相同冻结时间前提下，B 区的冻结速度大于A 区，C 区速度最慢。A 区是最靠近井筒中心的位置，冻结管距离该区域较远，但A 区并未接触周围土体，且没有其他热源对其产生影响，因此，A 区冻结速度在3 个区域里面适中。B 区介于主排孔和辅助孔之间，在主排孔和辅助孔冷量叠加的影响下，周围土体温度下降十分迅速，同时位于排孔之间的土体未接触外部热源，故B 区内土体温度下降速度最快。C 区内土体直接与外部土体接触，外部土层源源不断向该区域内土体提供热源，因此，冻结速度最慢。对于不同层位的土体，在同一冻结时间下，主面都具有相同的规律，即细砂降温效果最优，砂质黏土低于细砂，钙质黏土最差。

图8 埋深 218 m 钙质黏土层位2 个主面的温度随时空变化关系曲线Fig.8 Temperature variation curves of two main surfaces of calcareous clay with buried depth of 218 m

图9 埋深 225 m 细砂层位2 个主面的温度随时空变化关系曲线Fig.9 Temperature variation curves of two main surfaces of the fine sand layer with buried depth of 225 m

图10 埋深259 m 砂质黏土层位2 个主面的温度随时空变化关系曲线Fig.10 Temperature variation curves of two main surfaces of sandy clay with buried depth of 259 m

3.5 井帮温度

巷道开挖到埋深218、225、259 m 3 个层位的时间分别是冻结116、120、130 d，因此，本文分别提取模型中对应冻结天数的井帮温度与现场实测数据进行对比，井帮模拟温度和实测数据见表3，利用COMSOL Multiphysics 软件分别提取埋深218、225、259 m 3 个层位的井帮平均温度，为−6.51、−8.81、−8.10℃，井帮模拟平均温度与实测平均温度相差均在1℃以内，根据《煤矿冻结法开凿立井工程技术规范》，井筒垂深在150～250 m 的掘进段高内的黏土层井帮温度应在−4～−8℃，因此，模拟和实测结果均符合规范要求。

表3 3个层位的井帮温度计算结果与实测结果对比Table 3 Temperature comparison between measured results and numerical calulation results of hole side wall in three layers

4 结 论

a.数值模拟结果与实测结果的变化规律基本一致，表明了采用数值模拟来演化煤矿立井冻结壁温度场的可靠性和可行性。

b.测温孔实测数据指出，距离冻结孔越近的测温孔降温速度越快。在开挖至相应层位时，受混凝土水化热的影响，距离井壁越近的测温孔产生温度变化幅度越大。

c.在相同冻结时间条件下，冻结壁有效平均温度和平均厚度数值模拟结果均表明：细砂层位冻结效果最优，砂质黏土次之，钙质黏土最差。数值计算与现场实测均表明，冻结壁平均温度均低于−18℃，冻结壁有效厚度均达到6.6 m 以上，深部表土层(200 m 以下)开挖时井帮温度均在−4℃以下，满足设计要求的相应指标值。

d.冻结孔沿径向由内到外将冻结温度场划分为3个区域。同一冻结时间时，辅助孔与主排孔之间的B 区降温速度最快，井筒内A 区降温速度次之，主排孔外的C 区降温速度最慢。