天津地区新近系馆陶组热储层，除在潘庄凸起、双窑凸起核部缺失外普遍存在，具有埋深浅、温度适中、开发经济成本低的特点[1]。伴随地热开采量的不断增加，地下水位持续下降，孔隙性热储开发对深部热储特征有无影响，是否会引起水化学性质改变、地下温度场和地下水压力的下降，进而引发资源枯竭、地面沉降等环境问题，一直是社会关注的焦点。为了地热资源的可持续开发利用和保护，天津市开展了多年地热回灌工作[2-3]，随着回灌工作的不断持续，大部分热储水位下降幅度减缓或出现回升。同时，也发现不同区域地热井回灌能力差别很大，制约了一些地区对馆陶组热储资源的储量补偿。数据显示，大城凸起、武清凹陷、板桥凹陷、北塘凹陷地热井单井回灌量分别为34～85 m3/h、45～120 m3/h、60～109 m3/h、62～100 m3/h[4]。如何有效地提高馆陶组热储回灌能力，业已成为地热资源开发保护的焦点。赵苏民[5]、鲍卫和[6]、陈莹[7]等，从地热井成井工艺角度研究分析，认为过滤器可以减少水阻、增大过水性能有利增加回灌量；林建旺[8]、宗振海[9]、张平平[10]等认为在回灌初期，物理堵塞、化学堵塞和微生物堵塞不是影响衰减的主导因素，主要原因是冷锋面前移过程中，与热储层作用产生的一系列阻力影响所致；林黎[11]、王光辉[12]等认为储层水文地质特征、回灌流体物性特征以及灌井与采井间距、回灌方式等影响回灌能力；杨吉龙[13]认为加强信息化管理水平，引入第三方监管机制，实现地热资源实时化监管和精细化管理，可加强地热资源可持续开发利用。本文通过大量回灌试验数据及水岩试验数据，通过对地热流体特定离子浓度、下伏地层岩性、砂岩微观孔隙度分析，对天津市新近系馆陶组热储层回灌的影响因素进行定性分析，加深了对新近系馆陶组回灌影响因素的认识。
1 区域地质背景
天津地区在构造凸起区双窑凸起和潘庄凸起的核部，馆陶组缺失[14]，在缺失带两侧因沉积环境不同，导致馆陶组岩性、地下水径流、水化学类型等具有明显的差异性，形成完全不同的两套馆陶组砂岩热储。冀中坳陷区馆陶组热储，以辫状河沉积环境为主，武清凹陷为冀中北部地区的沉积中心，热储层沉积厚度最大，但埋深较浅，沉积旋回明显，整个剖面呈粗-细-粗三分性。受沉积环境影响，该区馆陶组地热流体处于深循环的地质环境，在平面上地下径流较慢，垂向上受岩性的影响，与浅层流体交换较小，该区域整理热储渗透系数较小。黄骅坳陷馆陶组热储为河流冲、洪积相沉积，该构造区域馆陶组砂岩热储渗透系数明显大于冀中坳陷，地下径流强烈。
虽然不同构造区域呈现出不同的热储特征，但沧县隆起、冀中坳陷、黄骅坳陷的馆陶组沉积物源相似，根据岩性分析，在馆陶组沉积时期，中北部发育冲积扇辨状河沉积，构成燕山山前冲积平原部分，而南部远离燕山物源影响范围，主要发育低弯度河流沉积[15]。
2 回灌试验
2.1 样本收集
区域内共收集整理2010 年至今54 处回灌井回灌试验数据[16]，其中包括射孔成井47 处，滤水管、筛管成井7处，主要分布在武清凹陷及北塘凹陷的西南侧，大城凸起与板桥凹陷数据相对较少。为了保证回灌数据有效性和综合研究性，结合目前馆陶组成井技术水平，选择射孔成井工艺的回灌井作为研究对象，共47处。回灌试验井位置见图1。
图1 馆陶组热储层分布图
Fig.1 Distribution map of thermal reservoir of Guantao Formation
2.2 参数计算
为了更加准确进行回灌能力研究，将不同构造单元，回灌试验获得的、不同温度的观测水位统一换算到25℃下的校准水位埋深，消除温度效应，获取更加科学合理的试验参数，可以更好地开展热储回灌能力评价工作。
利用回灌试验数据资料，采用裘布衣公式，求取回灌井水文地质参数[17]。假设地热井所在区域水文地质条件符合地下水向承压完整井运动的规律。同时，在数据资料整理过程中，考虑到温度对热储压力的影响，需对温度进行统一校正，消除地热井开发利用过程中井筒效应的影响。温度的传导扩散特性，认为地热井流体内温度与深度呈一维线性关系。校正水位由公式1、2计算：
回灌井的单位回灌量q回的计算公式为：
回灌井的注水渗透系数K注、导水系数T注，采用式4、5来计算：
式中：h校—校正后动水位埋深（m）；
H—热储层取水段中点的埋深（m）；
h静、h动—试验静水位、动水位埋深（m）；
ρ测、ρ热—回灌流体、热储温度密度（kg/m3）；
ρ校—校正后温度对应的密度（kg/m3）。
q回—单位回灌量（m3/h·m）；
Q注—试验测得稳定回灌量（m3/d）；
s注—回灌试验时稳定动水位与试验前静水位之差（m）；
K注—校正后（25℃）注水渗透系数（m/d）；
T注—校正后（25℃）注水导水系数（m2/d）；
M—热储有效厚度垂直厚度（m）；
R注—回灌影响半径（m）；
rw—回灌井热储段井半径（m）。
2.3 回灌效果分析
武清凹陷：试验数据4组。试验数据表明（表1）：WQ-23B井回灌能力最大，稳定最大回灌量达到120m3/h，由于该井位于构造单元的分界处，可能受构造的影响不具有代表性；其余3眼回灌井在武清凹陷上沿纵向分布，试验数据表明3眼回灌井回灌量自北向南逐渐衰减，3组回灌试验回灌量60～70m3/h[18]。
表1 回灌井试验数据计算结果表
Table 1 Calculation results of reinjection well test data
序号构造位置1 2 3 4 5 6 7 8 9 10板桥凹陷11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28北塘凹陷大城凸起武清凹陷井号TG-43 TG-29B TG-55B TGR-42D DG-28B TG-37B TG-46B TG-31B TGR-41D TG-61B TG-60B TGR-40D TG-62 TG-63B TG-21B TG-19B JH-26B WQ-06B WQ-36B WQ-38B XQ-25B BC-11B BC-13B BC-14B WQ-37B WQ-25B WQ-28B WQ-44校正后水位/m 106.81 111.11 101.27 135.32 130.2 87.67 82.06 128.62 121.03 131.33 129.18 124.00 108.26 108.70 127.70 125.54 54.40 111.97 109.76 121.83 88.14 115.27 118.46 122.74 139.08 54.87 131.69 113.63单位回灌量/（m³/h·m）0.71 0.70 1.09 1.21 2.00 2.82 7.13 4.85 2.92 0.94 0.98 2.76 1.05 4.28 2.02 2.07 1.30 1.31 1.39 0.78 0.77 1.50 2.02 1.58 1.91 1.04 1.44 1.30注水渗透系数/（m/d）0.28 0.18 0.44 0.36 0.51 0.38 3.15 1.48 0.79 0.21 0.3 0.73 0.57 0.91 0.56 0.47 0.35 0.33 0.44 0.2 0.59 0.49 0.56 0.53 0.36 0.28 0.39 0.11
大城凸起：试验数据10组。大城凸起现有回灌井主要集中在大成凸起与武清凹陷交界处，除XQ-25B 和JH-26B 井位于构造的南部区域外，其余8 眼回灌井与WQ-23B 井类似，分布在大城凸起和武清凹陷的分界线附近。试验结果表明：南部区域地热井最大回灌量为45-50m3/h，其余8组构造边界的回灌试验数据差异性较小，平均灌量为65～75 m3/h。
基于以上试验数据，馆陶缺失带以西的两个构造单元，馆陶组热储具有如下特征：在构造交界处，地热回灌能力受深部断裂及地层缺失影响，回灌能力较好，其他区域自南向北呈衰减之势。
北塘凹陷：试验数据13组。与冀中凹陷馆陶组地热井相比，黄骅坳陷回灌能力较好，北塘凹陷上的13组回灌试验平均回灌量80～85 m3/h。
板桥凹陷：试验数据8组。其中DG-07和DG-17井为石油改造地热井，受成井工艺影响，回灌能力不具备代表性，回灌效果60～80 m3/h。其余6眼地热井均在板桥凹陷北部区域，靠近海河断裂，回灌试验平均灌量75 m3/h，在本次试验中，热储回灌能力最强。
通过地热井回灌试验，北塘凹陷的馆陶组地热井回灌量80～85 m3/h，热储条件较好，回灌能力最强；大城凸起和板桥凹陷馆陶组热储回灌能力差异较大，但整体回灌效果较好，平均回灌量65～75 m3/h，武清凹陷馆陶组热储回灌量普遍低于70 m3/h，在本次样本中回灌能力一般。
为进一步分析不同构造单元的回灌能力，选取单位回灌量及注水渗透系数作为评价核心因素，对本次试验数据进行校正计算，剔除异常点，数据分析对比见表1，4个构造单元回灌能力排序为北塘凹陷>板桥凹陷>大成凸起>武清凹陷[16]。
3 回灌影响因素分析
3.1 砂岩胶结程度
砂砾岩的胶结成分主要包括泥质物质、粘土矿物、铁氧化物等，其中泥质物质是影响其胶结程度的主要因素，不同来源、不同沉积环境的砂砾岩胶结成分和胶结程度有着显著的差异，同时受砂岩孔隙度、饱和度、压实度等的影响。所以为了降低影响因素，选择同一构造单元、埋深相近的样本进行分析。
北塘凹陷、板桥凹陷、大城凸起、武清凹陷馆陶组底板埋深分别为1 800～2 200 m、1 800～2 100 m、1 300～1 800 m、1 700～2 700 m，埋深不同导致馆陶组砂岩胶结能力不同，进而影响不同构造单元上回灌能力的差异[19]。本次样本选择武清凹陷与板桥凹陷进行砂岩岩性测试，并进行对比分析（表2）。
表2 样品干密度测定统计表
Table 2 Statistical table for determination of the dry samples density
样品编号干密度（mg/L）1#1.73 2#1.77 3#1.73 4#1.59 5#1.77 6#1.65
试验过程：选取武清凹陷1#、2#、3#砂岩样本、板桥凹陷4#、5#、6#砂岩样本，依据试验技术流程，进行矿化度浸出试验，每隔12 h或24 h测量浸出液体的矿化度，试验持续时间180 h，测试精度±1mg/100 ml[16]。
结果分析：如图2 在武清凹陷样品试验结果，TDS排序为2#>1#>3#，与样本泥质含量检测结果正向耦合。在板桥凹陷样品试验结果中，TDS 排序为4#>6#>5#，与样本泥质含量检测结果正向耦合，但与武清凹陷样本测试结果不同的是，6#溶出速度前期较慢，但随着时间的延续，最终TDS与4#样本趋于一致，说明该样本虽其胶结程度较好，但其样本中泥质含量占比较大，溶解比例较大；5#样本TDS 最低，对比泥质含量测试结果，该样本值介于4#和6#之间，TDS 溶出较低主要由于样本干密度测定值最大，说明其致密程度高于4#和6#，溶出率较低。基于本次砂岩浸出试验，可以得出初步结论，在本次试验周期内，砂岩的干密度和泥质含量是决定砂岩样品溶出物质的快慢程度和TDS 浓度的两个核心因素，样品干密度与砂岩溶出物质速度呈反比，干密度越大TDS浓度越小；泥质含量与溶出物质速度呈正比，泥质含量越高TDS浓度越高。
图2 砂岩浸出试验对比图
Fig.2 Contrast map of Sandstone leaching test
综上，武清凹陷砂岩溶出物质速度相较板桥凹陷砂岩较慢，说明其胶结程度较高，直接影响孔隙型热储的裂隙发育度。通过回灌试验与砂岩测试结果表明，武清凹陷胶结能力大于板桥凹陷，从区域上表现为武清凹陷回灌能力弱于板桥凹陷。
3.2 泥质含量
依据公式4、公式5，影响回灌量大小的主要因素是渗透率，渗透率表征了岩石空隙中流体的流动能力。馆陶组砂岩热储渗透率主要决定于空隙微观结构，如喉道大小、孔隙大小及颗粒大小等，这些因素均受地层有效应力和地热流体运移的影响，而砂岩的胶结程度和含泥量正是影响储层有效应力和流体流动的重要影响因素，即对渗透系数影响最大的就是砂层中泥质含量，宏观上表现为砂岩储层中泥质含量与回灌量呈反比关系。
本次对馆陶组砂岩热储岩样进行取样测试，样本数量8组，根据测试结果绘制砂岩样本孔隙率与泥质含量拟合曲线（图3），结合前人研究成果分析二者直接的拟合关系，其中包括的天津大港和河南周口[20]3个样本数据，以及渤海湾盆地东濮凹陷沙河街组[21]4个样本数据。本次试验数据拟合曲线可以看出，随着砂岩中泥质含量的增加，溶出孔隙率呈自然对数递减趋势，通过数据拟合，建立Φ=49.6-9.6ln（C）回归方程。
图3 砂岩孔隙度与泥质含量关系图
Fig.3 Relationship between porosity and shale content of sandstone
相比之下，文献[20]回归方程为：
文献[21]回归方程为：
通过比对3 组测定数据，回归方程函数形式类似，本次曲线正好拟合在前2组模型之间，所以认为本次样本分析结果合理可靠，砂岩热储孔隙度随着泥质含量的增大而减少，但泥质含量增加导致的孔隙率减少的趋势呈衰减之势，不再是核心影响因素。
为了进一步研究砂岩热储中泥质含量对孔隙度的影响，对本次砂岩样品进行扫描孔隙率测定，与地球物理测井过程中获得的表观孔隙率相比，并计算出相对偏差（DE）[16]。
由图4 可以看出，孔隙度的相对偏差值（DE）与泥质含量（C）之间具有典型的一次线性关系DE=4.6267C-20.114（R2=0.9905）。随着泥质含量的增加，孔隙度偏差由2.7%逐渐增大至197%。从该拟合曲线可以推断如下结论：1）砂岩热储中泥质成分主要由非常微小的矿物颗粒组成，受限于CT扫描的分辨率，对样本内部极其微细的孔隙无法测得；2）并不是所有的孔隙都具有导水性，如微小闭孔和黏土矿物层间孔隙，这些孔隙一方面不参与地热流体的运移，另一方面受流体的表面张力与泥土矿物的吸水性影响，阻碍地热流体的运移，起到负面影响。基于以上结论，也解释了为什么砂岩热储成井时，地球物理测井显示孔隙度较好，但地热井涌水或回灌效果不佳。所以当要对热储地质条件进行特定研究时，不能仅凭测井试验得到的表观孔隙率来判断砂岩的渗透性，要综合考虑影响因素，否则可能出现测试结果与实际情况相反的情况。例如本次研究，通过表观孔隙度选取样品武清凹陷孔隙度大于板桥凹陷，但根据岩心测试3D板桥凹陷孔隙度大于武清凹陷，与二者回灌能力强弱相吻合，呈正相关关系。所以想要定性定量地分析泥质含量对回灌的影响，必须进行必要的岩样测试，不能仅仅根据测井解释结果判断馆陶组回灌能力。
图4 砂岩孔隙率计算偏差随泥质含量变化图
Fig.4 Variation of porosity calculation deviation of sandstone with shale content
3.3 水动力对回灌的影响
一般情况下水动力条件的强弱可以用地热流体中离子浓度表征，即离子浓度大的地方，表示水流缓滞，水力联系较弱；反之离子浓度小的地方，水动力条件强烈[22-24]。
馆陶组热储地热流体从北到南水化学类型由简单的HCO3-Na型逐渐变为复杂的Cl·SO4-Na型。缺失带西侧，馆陶组热储地热流体从北到南总矿化度由686.7 mg/L逐渐增大到2 200 mg/L，与地热流体径流的方向基本一致；缺失带东侧也呈同样的变化趋势，总矿化度由北部的483.4 mg/L 逐渐增大到南部的1 925 mg/L。
取样过程：为更好地分析离子浓度对回灌的影响，进行馆陶组地热流体水质取样。武清凹陷、大城凸起、北塘凹陷、板桥凹陷分别取样5组、8组、7组、9组，送天津市地质矿产测试中心进行地热水质全分析，获得结果绘制离子浓度与矿化度关系曲线（图5、图6）。
图5 氯离子与矿化度关系曲线图
Fig.5 curve of relationship between chloride ion and mineralization degree
图6 钙离子与矿化度关系曲线图
Fig.6 curve of relationship between calcium ion and mineralization degree
测试结果分析：氯离子一般来自沉积岩中岩盐、氯化物的溶解，岩浆岩含氯磷灰石、方钠石的风化溶解、火山喷发物的溶滤以及海水氯离子的溶解度很大，随着矿化度的增大，钠离子浓度增大，而且钠离子在水中最为稳定，不易被吸附，在水中迁移能力强，一般浓度变化趋势线与地下水径流方向一致。
由图5我们可以看出矿化度与氯离子线性关系拟合度很高。同时，我们得到各构造单元地热流体中氯离子浓度排序:板桥凹陷>北塘凹陷>大城凸起>武清凹陷，这个顺序基本与回灌能力一致。由图6，我们可以看出钙离子浓度在冀中凹陷（<10 mg/L）与黄骅凹陷（>20 mg/L）具有明显的分区性，在冀中凹陷中大城凸起钙离子浓度略低于武清凹陷，黄骅凹陷中北塘凹陷钙离子浓度明显低于板桥凹陷。
综上，钙离子越低表明地下径流溶滤作用越强，氯离子浓度越低表明地下水运移越快。二者与馆陶组热储层单位回灌量大小呈负相关关系，即钙离子、氯离子低的区域，更有利于储层回灌，表征在区域上不同构造单元回灌量大小关系为板桥凹陷>北塘凹陷>大城凸起>武清凹陷。
3.4 下伏地层对回灌能力的影响
武清凹陷下伏地层为古近系，但同处冀中凹陷的大城凸起北部下伏地层为中生界，南部地区为石炭-二叠系，本次选取这三个区域地热井，分析下伏地层对回灌能力是否存在明显影响[25-26]。
数据统计：下伏地层为石炭二叠系馆陶组JN-09、JN-21 地热井单位回灌量为0.87～3.73 m³/h·m，平均为1.77 m³/h·m；下伏地层为中生界单位回灌量分别为0.77～1.78 m³/h·m，平均为1.41 m³/h·m；下伏地层为古近系时单位回灌量分别为0.36～1.91 m³/h·m，平均为1.2 m³/h·m。
原因分析：大城凸起位于沧县隆起西侧，为一个走向NNE、倾向NWW的单斜断块构造，由于古近系的缺失，新近系馆陶组下伏地层为中生界或石炭系基岩地层。研究表明[27]，下伏地层岩层倾角陡缓程度不同，其受力机制也不同，会对后期新沉积的新近系馆陶组成岩裂隙产生影响。基岩裂隙的发育，导致竖向越流补给增强，使得新近系馆陶组与中生界或石炭系地层之间有一定水力联系，进而影响上覆馆陶组的回灌效果。同时，下伏地层为石炭系的馆陶组储层的区域靠近沧县隆起核部，地层倾角更大，基岩更破碎，回灌试验结果表明回灌效果更好。而武清凹陷下伏地层为古近系，岩性以厚层泥岩为主，在某种意义上阻断了馆陶组与下伏地层间的垂向越流联系，导致武清凹陷馆陶组回灌效果较弱。
4 结论
影响馆陶组热储回灌能力的因素是复杂的，本文通过对水岩试验、回灌试验等数据的分析，初步认为影响馆陶组砂岩热储回灌的主导因素是孔隙度和泥质含量。微观上，泥质含量的变化，引起砂岩的表观孔隙度和扫描孔隙率的差异性，所以不能单一通过测井孔隙度来定性回灌能力的大小，热储砂岩的胶结性、矿物质溶出快慢程度也决定了回灌效果；宏观上，下伏地层的岩性构造，影响热储孔隙度发育和各水层之间的水力联系，进而影响地热流体的运移能力，可以通过离子浓度变化表征出来。本文通过本次样本试验，对馆陶组回灌能力进行分析研究，得出如下结论：
（1）多年回灌数据表明，天津地区馆陶组回灌能力从大到小为：塘沽凹陷>板桥凹陷>大城凸起>武清凹陷。
（2）砂岩热储泥质含量是影响回灌能力的主要因素，其直接影响砂岩3D孔隙度，通过样本测试，武清与滨海新区馆陶组砂岩热储孔隙度与泥质含量满足Φ=49.6-9.6ln(C)的对数关系，当泥质含量超过30%，孔隙度陡降，且孔隙度偏差增大。
（3）馆陶组砂岩热储孔隙度发育与水岩反应有关，砂岩胶结度越大，越不利于砂岩热储孔隙发育，钙、氯离子浓度低的区域回灌效果更好，表征的地热流体溶滤作用和迁移能力越强。
（4）馆陶组砂岩储层下伏地层为基岩时，相较古近系回灌效果更好，这主要与层间越流与垂向构造裂隙发育有关。
[1]地质矿产部地质专报.天津市区域地质志[M].北京:地质出版社，1982：205-230.
[2]殷肖肖，闫佳贤，李虎，等.2021年度天津市地热资源开发利用动态监测报告[R].天津：天津地热勘查开发设计院，2021.
[3]殷肖肖，赵苏民，蔡芸，等.近三十年天津市地热大规模开发热储动态特征研究[J/OL].地质学报:1-17.
[4]高新智，赵艳婷，杨宝美，等.天津市孔隙型砂岩热储垂直向对比回灌参数研究[R].天津：天津地热勘查开发设计院，
[5]赵苏民，孙宝成，马忠平，等.天津地区孔隙型回灌井井身结构与完井技术[A].地温资源与地源热泵技术应用论文集（第三集）[C].2009:6.
[6]鲍卫和.天津地区新近系回灌井成井工艺研究[R].天津市帅源科技有限公司，2008.
[7]陈莹，王攀科，吴烨等.河南兰考地区地热回灌影响因素分析及对策[J].钻探工程，2022，49（06）:146-152.
[8]林建旺，赵苏民.天津地区馆陶组热储回灌量衰减原因探讨[J].水文地质工程地质，2010，37（05）:133-136.
[9]宗振海，王平，王冰，等.稳态测温方法在天津地热动态监测中的应用[J].地质调查与研究，2018，41（03）:228-231.
[10]张平平.山东省孔隙热储地热回灌堵塞机理研究[R].山东省地质矿产勘查开发局第二水文地质工程地质大队（山东省鲁北地质工程勘察院），2019.
[11]林黎，王连成，赵苏民，等.天津地区孔隙型热储层地热流体回灌影响因素探讨[J].水文地质工程地质，2008，（06）:125-128.
[12]王光辉，赵娜，唐永香，等.孔隙型地热资源回灌模式研究—以天津市滨海新区为例[J].地质调查与研究，2014，37（02）:155-160.
[13]杨吉龙，汪大明，牛文超，等.天津地热资源开发利用前景及存在问题[J].华北地质，2022，45（03）:1-6.
[14]田树信，翟子梅.天津市岩石地层[M].武汉:中国地质大学出版社，1996.
[15]赵苏民，孙宝成，林黎.沉积盆地型地热田勘查开发与利用[M].北京：地质出版社，2013.
[16]高新智.天津市孔隙型砂岩热储回灌能力评价模型研究[D].中国地质大学（北京），2018
[17]薛禹群.地下水动力学[M].北京：地质出版社，1997.
[18]沈健.天津市东丽湖地表水热储回灌技术研究[D].中国地质大学（北京），2015.
[19]宋美钰，刘杰.天津地热流体水化学特征及同位素分析[J].地质调查与研究，2018，41（02）:138-144.
[20]王华军，于萌，赵苏民，等.新近系地热钻孔砂岩孔隙微观分析及渗透特性对比[J].太阳能学报，2019，40（06）:1790-1796.
[21]王瑞飞，何润华，苏道敏，等.深层砂岩油藏储层孔喉特征参数及预测模型[J].地球物理学进展，2016，31（05）:2160-2165.
[22]孙宝成，曾梅香，林黎，等.天津地热对井回灌系统中的同位素示踪技术[J].地质调查与研究，2005，（03）:187-192.
[23]张芬娜，贾志，李虎，等.天津塘沽地区馆陶组地热流体演变特征研究[J].地质调查与研究，2016，39（02）:149-152.
[24]靳宝珍，杨永江，李俊峰，等.天津静海县史家庄一带馆陶组地热流体质量及实用性[J].地质调查与研究，2007，（03）:219-223.
[25]阮传侠.天津地区雾迷山组热储地热回灌研究[D].中国地质大学（北京），2018.
[26]江国胜，赵苏民，李虎，等.天津市奥陶系热储赋水规律研究[R].天津：天津地热勘查开发设计院，2015：9-12.
[27]庞菊梅，庞忠和，孔彦龙，等.岩溶热储井间连通性的示踪研究[A].中国科学院地质与地球物理研究所2014 年度（第14届）学术年会论文汇编[C].2015.
[28]王连成.天津市新近系馆陶组地热流体回灌研究[D].中国地质大学（北京），2014.
[29]周鑫.沉积盆地孔隙型地下热水回灌堵塞机理研究[D].长安大学，2013.
[30]申小龙，李浩，王伟，等.关中西部眉县城区地热井回灌试验Q-p 曲线方程探究[J].煤田地质与勘探，2022，50（06）:165-174.