我国地热资源丰富，约占全球资源总量的8%^[1]。地热资源最常用的开采方法是直接采水取热，但这种方法的适用范围有限^[2-3]。1980年，R. N. Horne^[4]首次提出井下同轴闭式换热器的概念，建立了二维流动传热数值模型；K. Morita等人^[5-6]提出了绝热中心管的概念，并试验验证了其可行性；在此基础上，研究人员分析了工程参数、地质参数^{[2, 7-8]}及保温材料^[9-10]对产能的影响，现场应用中取得了较好的效果^[11-12]。

石油工程及地热领域广泛应用模糊数学方法。赵振峰等人^[13]采用模糊集合理论，进行了体积压裂水平井的产量预测；朱晓东等人^[14]建立了分层模糊系统，用来预测石油钻井参数，并对预测模型进行了训练和优化；局培等人^[15-16]采用模糊综合评判法，进行了PDC钻头优选及欠平衡钻井方式优选；赵春森等人^[17]采用模糊综合评判法，评价了注水井分层注水效果；朱兆群等人^[18]利用改进的模糊–灰色综合评判方法，对储层进行了定量评价；方舟等人^[19]利用模糊数学方法，对海上稠油热采效果进行了评价；张伟等人^[20]建立了适用于高温地热EGS开发方案的综合评价方法，用来优选目标地层的注采方案。

闭式系统主要利用循环工质从地热储层提取热量，循环工质对换热性能影响较大；但是目前评价循环工质换热性能的方法仅通过单一因素评价，无法全面评价循环工质的换热性能。为此，笔者首先建立了三维流动传热瞬态模型，利用夏威夷HGP-A井的现场数据进行验证；然后采用层次分析法，确定出口温度、换热功率、循环压耗及性能系数等4个指标的权重；最后基于模糊综合评判法，建立了综合评价模型，对9种循环工质的换热性能进行综合评价优选。

1.   井下同轴闭式地热系统数值模型建立

1.1   模型建立及假设

井下同轴闭式地热系统主要包括地热储层中的一口直井和井筒中呈同轴位置关系的套管与中心管（见图1），其中中心管一般使用保温材料^[9]。循环方式主要有正循环和反循环2种方式。前人已经证明反循环的方式具有更高的换热效率^[6]，即将循环工质从地面通过高压泵注入环空，然后通过中心管返回地面。由于储层和井筒内流体存在温度差异，环空中的循环工质会从井壁获取热量，加热后的循环工质流经热交换器被利用。井筒和储层之间由套管和水泥分隔，套管和水泥的导热系数影响换热效率^[2]。并在其内部发生热传导。循环工质与周围储层之间存在较大温差，井筒与储层之间会产生热对流和热传导，这可以弥补井筒周围的热损耗，并且地下水的流动可以加快这一过程。

图  1  井下同轴闭式地热系统的组成及换热过程

Figure  1.  Composition and heat exchange process of the CBHE closed-loop geothermal system

下载: 全尺寸图片 幻灯片

基于华北地区地质条件，建立了耦合井筒和储层三维非稳态数值模型，模型进行以下假设：1）套管、水泥和岩石等固体部分是均质的、具有各向同性；2）固体部分的热物理性质恒定且与温度无关；3）岩石是水饱和介质；4）地下水沿水平方向流动且满足达西定律，储层中传热满足局部热平衡理论；5）利用非等温管流模型对环空和内管中流体的流动传热进行描述，忽略循环工质的径向流动；6）考虑华北地区的高静水压力和相对较小的地热梯度，不考虑循环工质的相变。

1.2   控制方程

地热储层中，考虑局部热平衡，储层岩石与地热流体之间的热交换可表示为：

式中：${\left( {\rho {C_{\rm{p}}}} \right)_{{\rm{eff}}}}$为有效体积热容，J/（K·m³）；T_r为储层流体温度，K；ρ_g为储层流体的密度，kg/m³；C_p,g为储层流体热容，J/（kg·K）；${v_{\rm{g}}}$为储层流体流速，m/s；q为热源项，W；${\lambda _{{\rm{eff}}}}$为有效导热系数，W/（m·K）；ϕ为储层的孔隙度；λ_s为储层岩石导热系数，W/（m·K）；λ_g为储层流体导热系数，W/（m·K）；ρ_s为储层岩石的密度，kg/m³；C_p,s为岩石热容，J/（kg·K）。

一般用达西定律来描述地热储层中的流体流动，控制方程如下：

式中：$K$为储层渗透率，m²；$\mu$为储层流体黏度，Pa·s；z为竖直方向距离，m；g为重力加速度，m/s²。

采用非等温管流来描述井筒内环空和中心管中循环工质的流动和传热过程，相关的控制方程如下：

式中：A_p为管横截面积，m²；f_D为达西摩擦因子，由Churchill模型^[21]确定；ρ_f为工质密度，kg/m³；C_p,f为工质等压热容，J/(kg·K)；v_f为工质流速，m/s；p为中心管中的压力，Pa；d_p为管内径，m；λ_f为工质导热系数，W/（m·K）；Q_wall为地热流体通过换热器壁传递给循环工质的热量，W/m；${T_{\rm{ext}}}$和${T_{\rm{f}}}$分别为中心管外部和内部的温度，K；$Z$为管壁的周长，m；${\left( {hZ} \right)_{\rm{eff}}}$为等效传热系数，W/（m·K）；$h$为对流换热系数，W/（m²·K）；${T_{\rm{ext}}}$为管道外的温度，K，可以通过求解方程（1）—（5）得到。

h和${\left( {hZ} \right)_{\rm{eff}}}$的计算公式为：

式中：r_i和r_o分别为管内径和管外径，m；h_int为内部对流换热系数，W/（m²·K）；h_ext为外部对流换热系数，W/（m²·K）；λ_p为管导热系数，W/（m·K）；$Nu$为努塞尔数，可通过Gnielinski方程^[22]得到：

式中：Re和Pr分别为雷诺数和普朗特数。

根据上述方程，可以将井筒和地热储层中的流动传热通过${Q_{\rm{wall}}}$和${T_{\rm{ext}}}$进行耦合。

1.3   模型验证

为了验证模型的准确性，将夏威夷HGP-A井的生产数据^[22]与模型的模拟结果进行对比。1991年，K. Morita等人在美国夏威夷Kapoho地区的HGP-A井进行了现场试验^[23]。该井位于基拉韦厄东部裂谷带，距基拉韦厄火山口42 km，完钻井深1 962.00 m，后来一直关井；改造时下入封隔器，试验井段为0～876.50 m，井身结构如图2所示。

图  2  HGP-A井的井身结构

Figure  2.  Casing program of Well HGP-A

下载: 全尺寸图片 幻灯片

该井井底初始温度为110 ℃，地层平均孔隙度为30%，地层导热系数为1.60 W/（m·℃），水泥导热系数为0.99 W/（m·℃），套管导热系数为46.10 W/（m·℃），保温内管导热系数为0.06 W/（m·℃）。试验时，排量保持4.8 m³/h，入口温度保持30 ℃。首先，将模型的储层参数、运行参数及井身结构调整与HGP-A 井相同，进行数值模拟；然后，将前7 d的生产数据与模拟结果进行对比，结果如图3所示。从图3可以看出，数据吻合较好，说明建立的模型准确可靠。

图  3  模拟结果与夏威夷HGP-A井生产数据的对比^[23]

Figure  3.  Comparison between simulation results and production data of Well HGP-A^[23]

下载: 全尺寸图片 幻灯片

2.   综合评价模型

2.1   换热性能评价指标选取

结合现场实际和文献调研，选取以下4个评价指标，评价井下同轴闭式换热系统的换热性能。

1）出口温度T。评价循环工质换热性能最直接的因素是出口温度，表征循环工质通过在环空和中心管循环换热之后达到的温度。

2）换热功率P_out。换热功率是评价循环工质换热性能最重要的参数，表示循环工质从地热储层中提取热量的能力：

式中：P_out为系统换热功率，kW；T_in和T_out分别为循环工质入口处和出口处的温度，K；C_p,in和C_p,out分别为循环工质入口和出口处的恒压热容，J/（kg·K）；Q_in和Q_out分别为循环工质入口和出口处的质量流率，kg/s。

3）循环压耗p_c。循环压耗表示工质在循环换热过程中的压力损失，影响注入泵的选择，从而影响生产成本。计算式为：

式中：p_out和p_in为分别为出口压力和注入压力，MPa。

4）性能系数（coefficient of performance，COP）。性能系数表示单位能耗的换热功率，代表循环工质在系统损失单位能耗时所能得到的热功率。计算公式为：

式中：C_p为性能系数；P_con为系统能耗，kW。

2.2   隶属函数确定

不同的性能评价指标单位不同，因此无法直接进行比较，利用比例转换法，对属性指标参数归一化。

对于越大越优型的收益类属性指标，其隶属函数为：

对于越小越优型的支出类属性指标，其隶属函数为：

式中：${a_{ij}}$表示第j个特征样本中第i个元素，${a_{j,\min }}$表示第j个特征样本中的最小值；${a_{j,\max }}$表示第j个特征样本中的最大值；i=1, 2,…, n；j=1, 2,…, m。

由此可得隶属函数归一化矩阵B：

2.3   权重确定

对出口温度、取热功率、循环压耗及性能系数等4个指标进行层次分析，按照“1—9标度法”构造判断矩阵，结果如表1所示。

表  1  指标判断矩阵

Table  1.  Index judgment matrix

指标	出口温度（B1）	换热功率（B2）	循环压耗（B3）	性能系数（B4）
出口温度（B1）	1	1/8	1/2	1/4
取热功率（B2）	8	1	4	2
循环压耗（B3）	2	1/4	1	1/2
性能系数（B4）	4	1/2	2	1

下载: 导出CSV 

| 显示表格

计算得到综合矩阵的最大特征值为4，对应的特征向量为[–0.1085　 –0.8677　 –0.2169　 –0.4339]，归一化得[0.0667　0.5333　0.1333　0.2667]。

对矩阵进行一致性检验。一致性指标$I_{\rm{C}}$为：

式中：${\lambda _{\max }}$和$n$别为判断矩阵的最大特征根和阶数。

查表得，平均随机一致性指标I_R=0.9。随机一致性比率$R_{\rm{C}}$为：

采用MATLAB，计算$I_{\rm{C}}$=0，$R_{\rm{C}}$=0<0.1，说明满足一致性要求，权重分配是合理的。即权重向量W=[0.0667　0.5333　0.1333　0.2667]。

2.4   综合评价模型

根据上面利用层次分析法求出的各个指标的权重矩阵 W和各个指标的隶属函数矩阵 B，计算各个样本的综合评价值E。E值越大，表明循环工质换热性能越好；反之，则表明换热性能越差。

3.   数值模拟及综合评价结果

3.1   数值模拟结果

选取水、CO₂及7种有机工质等9种工质进行对比分析，9种工质的基本性能参数如表2所示。通过美国国家标准与技术研究院（NIST）的工质物性查询平台REFPROP，得到不同工质的热物性参数，包括密度、恒压热容、黏度、导热系数和比热率。考虑循环工质的压力变化（1～20 MPa）和温度变化（273.15～373.15 K），利用MATLAB软件建立工质真实物性数据集，涵盖了井内温度压力变化范围，将数据集输入根据COMSOL建立的模型进行运算。

表  2  循环工质基本性能参数

Table  2.  Basic parameters of the circulating working fluids

工质	化学名称	临界压力/MPa	临界温度/℃
水		22.12	374.15
CO2	二氧化碳	7.38	31.20
戊烷	正戊烷	3.37	196.40
R134a	1,1,1,2-四氟乙烷	4.06	101.06
R152a	1,1-二氟乙烷	4.52	113.26
R227ea	1,1,1,2,3,3,3-七氟丙烷	2.93	101.65
R245fa	1,1,1,3,3-五氟丙烷	3.65	154.01
1234ze	1,3,3,3-四氟丙烯	3.64	109.37
R600a	异丁烷	3.66	134.98

下载: 导出CSV 

| 显示表格

模型采用河北雄安新区西柳-1井基础数据作为参考。该井是一口开采地热流体的直井，后来由于产能过低被废弃。该地热井完钻井深2 530.00 m，在目的层段井深1 800.00和1 810.00 m下入裸眼封隔器进行封隔。由现场测试数据得知，地温梯度为2.7 ℃/100 m，井底初始温度为57 ℃，地层平均孔隙度为20.0%，储层基质平均渗透率为200 mD。地层导热系数为3.0 W/（m·℃），水泥环导热系数为0.7 W/（m·℃），油/套管导热系数为43.75 W/（m·℃）。地下水流速1.27×10^–7 m/s。设定循环工质在井内的初始温度与地层温度相同，储层边界处无流动，模型的上下边界为绝热表面。从现场条件来看，入口温度和流量分别保持在10 ℃和20 m³/h；中国北方采暖期一般为4个月左右，模型的运行时间设定为120 d。分别计算运行120 d后的4个换热性能评价指标，结果如图4所示。

图  4  不同循环工质的换热性能指标计算结果

Figure  4.  Calculation results of the indexes for heat exchange performance of different circulating working fluids

下载: 全尺寸图片 幻灯片

从图4可以看出，各种循环工质在不同的换热性能评价指标上表现各异。R600a具有最高的出口温度及最小的循环压耗，但是它的换热功率却最低，其出口温度比最低的水高出56.92%，循环压耗比最高的R227ea低59.35%，而换热功率比最高的CO₂低32.61%。CO₂除了换热功率最大外，性能系数也最大，其性能系数比最低的R227ea高出125.56%。但根据以上结果，无法对不同循环工质的换热性能进行统一且全面的评价，因此，需要建立综合评价方法对其进行评价，来优选最优循环工质。

3.2   综合评价结果

下面利用建立的综合评价方法评价9种循环工质，4种换热性能评价指标初始值构成原始决策特征矩阵：

利用指标参数归一化求隶属函数的方法，对决策矩阵进行处理。其中，出口温度、换热功率和性能系数属于收益类指标，可用式（16）计算；循环压耗属于支出类指标，可采用式（17）计算。按照以上方法，利用MATLAB编程对矩阵进行归一化处理，得到归一化矩阵：

利用层次分析法求得的权重，对模糊综合评判矩阵进行加权，得到综合评价结果。由式（21）计算可得：

将综合评价结果转换为百分制，得到不同循环工质换热性能的综合评价值。根据综合评价值，可将循环工质划分为3类：优（综合评价值≥90）、较优（70≤综合评价值<90）和一般（综合评价值<70），不同工质的综合评价值及等级如表3所示。CO₂综合评价值高达95.12，换热性能最好，认为是最优循环工质，在后续现场应用中可以进行验证并推广使用；R600a、戊烷、水和R152a的综合评价值均大于70，也是不错的选择；虽然水的换热性能表现并不特别突出，但其成本最低，所以现场应用最为普遍。

表  3  不同循环工质的综合评价值及等级

Table  3.  Comprehensive evaluation scores and grades of different working fluids

循环工质	综合评价值	等级
CO2	95.12	优
R600a	81.04	较优
戊烷	77.54
水	75.95
R152a	73.47
R134a	68.83	一般
R1234ze	67.87
R227ea	64.90
R245fa	64.75

下载: 导出CSV 

| 显示表格

4.   结论与建议

1）基于现场基础数据，利用有限元软件COMSOL，建立了井下同轴闭式地热系统三维非稳态流动传热模型，利用夏威夷HGP-A井生产数据验证了模型准确性。选取了4个换热性能评价指标，采用层次分析法确定评价指标的权重，建立了循环工质换热性能综合评价方法。

2）首先确定了9种循环工质，通过NIST工质物性查询平台REFPROP得到工质真实物性；然后输入数值模型进行计算，得到不同循环工质的4个评价指标初始值；最后利用综合评价方法评价了9种循环工质的换热性能，评价结果表明CO₂的综合换热性能最好，认为是最优循环工质。

3）研究结果为优选井下同轴闭式地热系统循环工质提供了新的思路，可以采用CO₂作为循环工质来提高系统的换热性能；但是目前尚未分析各种循环工质在实际应用中的换热效果，还需进一步通过现场试验验证不同循环工质的换热性能。