深井回灌式水源热泵井群运行的地下含水层传、蓄热性能模拟研究

１　问题的提出
　　
　　深井回灌式水源热泵技术作为一种有益于环境保护和可持续的冷热源形式，在国内外空调工程界已经得到了越来越多的应用^{[1] [2]}，^[3]给出了其基本原理与相关技术分析。这一系统方式利用温度全年相对恒定的地下水作为水源热泵的水源，通过建造抽水及回灌井群，实现夏季抽冷水、灌热水，冬季抽热水、灌冷水的这一全年角色轮换的运行过程，地下含水层内部的热量或冷量被提取、蓄存和转移。井群是深井回灌式水源热泵系统的一个关键组成部分，其正常运行与否决定了应用水源热泵系统工程的成败，井群的设计布局应当是慎之又慎的关键环节。目前国内进行此类工程的井群设计和施工过程中，系统方案的可行性判据基本取决于单井出水量是否满足要求，以及能否实现良好的人工回灌。然而在进行该类工程井群部分的可行性分析和设计中，还需要考虑以下几方面的问题：
　　（1）当地含水层中的能量蓄存、转移过程。
　　应用深井回灌方式，需要在设计阶段知道当地含水层的能量提供能力有多大，也就是系统可负担的建筑容积率极限是多
少。当建筑物全年冷热负荷不均匀时，系统对于含水层温度的常年影响效果如何，是否会造成含水层"背景温度"逐年降低或升高，从而导致系统运行失败，相关的应对策略如何制定？
　　（2）"热贯通"影响的避免
　　由于回灌水与原始含水层温度存在的差异，在导热和对流等作用下，回灌井水"温度锋面"会导致近抽水井出水温度有不同程度的升高或降低，通常称为"热贯通"现象。如何确定适宜的井间距，如何确定井群的布局，避免"热贯通"的影响，是设计人员关心的主要问题。对于高密度住宅小区或城区商用建筑应用深井回灌式水源热泵系统来说，由于可利用建筑用地的面积限
制，如何优化井群布局及其各自对应的抽水或回灌角色，最大限度地避免"热贯通"的不利影响是尤为关键的。
　　（3）水文地质条件的影响作用
　　该问题的核心是如何考虑速度相对较大的当地地下水流动的存在对于地下含水层温度场和井群布局的影响，如何在建筑物冷全年热负荷不均匀的情况下，利用自然地下水流场的存在，合理地优化各井的抽水、回灌角色和轮换方式，从而实现对于能源的最优利用。
　　（4）地面机组和管道系统形式和运行模式的影响作用
　　利用"小流量、大温差"的系统运行方式，能够实现对于含水层蓄能的最大利用，同时减少对于地下水资源最小程度的开采利用^[3]。在"小流量、大温差"和传统的"大流量、小温差"两种运行工况下，地下含水量水层温度场全年变化过程，以及所导致对应的井群布局差异如何，也是研究设计人员所关心的问题。
　　以上这些问题都需要寻求对于井群部分的含水层水热运动过程适用的分析工具，研究井群当地地下含水层的水热运动与水文地质条件、环境气象因素和工程措施之间的关系，为该类系统的相关设计与分析提供的依据，以推动集中式水源热泵机组应用的进一步推广。
　　在此工程背景下，笔者通过比较目前工程学术界流行的含水层流动传热模拟程序，选择利用了美国地质调查局开发的地下水流动、传热、传质三维有限差分模拟程序HST3D，对一典型双井承压含水层进行了全年温度场和流场模拟，对该程序应用本问题的功能性和适用性作出评价，指出其需要完善之处。
　　
2 含水层中的水热运动及相关数值模拟
　　
　　2．1 含水层中的水热运动
　　
　　自然界含水层中的地下水流动一般满足达西定律，而含水层内部的传热过程包括：
　　（1） 地下水的对流换热过程；
　　（2） 地下水的导热过程；
　　（3） 固体骨架的导热过程；
　　（4） 由于通过多孔介质孔隙的不同流动通道液体的机械混合造成的局部热弥散，以及由于不同地质成分构造混合所造成的宏观热弥散过程；
　　（5） 地下水与固体骨架之间的传热。
　　在含水层传热中，当含水层骨架颗粒较小和流体流动雷诺数较低时，可以假定流体温度场与固体骨架的温度场时一致的，因此可以不考虑上面第5项的热量传递。
　　第4项由含水层多孔介质热弥散机理所造成的热量传递，对于含水层内部的整体传热过程存在着不可忽略的影响，特别是单井周边含水层中由于相对较高的地下水孔隙流速，热弥散的影响作用尤为突出。
　　
　　2．2 含水层水热运动数值模拟
　　含水层水热运动的建模和模拟计算工作，在含水层季节性热蓄能等相关领域已经进行了较长时间的深入研究。含水层季节性热蓄能的应用和研究，早于上世纪七十年代中期已经相继在我国、北欧、北美等地区和国家开展。美国加州大学的Lawrence Berkeley Laboratory (LBL)建立了相应的单井含水层蓄能有限差分数值模拟程序CCC，并通过该程序对现场实验进行了模拟和分析^[4]。国内陈兆祥^[5]和薛禹群^[6]等亦完成了相关模拟计算和现场实验工作。国际能源组织于1993年至2000年期间所完成的IEA ANNEX 8： Implementing Underground Thermal Energy Storage，其中一个子项目就是关于地热蓄能的设计分析工具应用与评价，G?ran Helstr?m ^[7]在其为该子项目所作的报告中列出了适用于含水层水热运动分析现行软件，包括AST、TWOW、SUTRA、Tradikon、HST3D等程序。Chiasson^[8]的中亦列出了适用于含水层水热运动分析的现行软件，除以上所列以外，还包括SWIFT、AQUA3D、FEFLOW等。下面将对HST3D程序作原理和功能介绍。
　　
3 HST3D简介
　　
　　HST3D^[9][10]是Three-Dimensional flow，Heat, and Solute Transport model 的简称，它是美国地质调查局（USGS）于80年代末开发的一套开放型研究用程序。HST3D采用控制容积的能量平衡法对三维流动、传热和传质微分方程进行离散求解，能够实现饱和含水层中流动、传热和传质过程的非稳态模拟，可用于饱和地下含水层相关流动、传热和传质问题的模拟，包括热田和土壤热、海水入侵、放射性核废料填埋等问题。HST3D具有很强的实用性，其性能包括程序设计结构的模块化、离散方法的简单化和求解方法的多样化、允许采用多种坐标系及不等距等，其开放性的模块化结构信纸科研人员可以根据需要添加、修改或删除相应的模块。
　　HST3D所求解的流动、传热以及物性方程^[11]分别如下：
　　饱和含水层的流动微分方程：
　　　　　 　　　　　　　　　　　　　　　　　　 （1）
　　饱和含水层的传热微分方程：
　　　　　 　　　　　　　　　　　　　　　　　　 （2）
　　假定密度ρ为压力和温度的函数，其液体物性方程：
　　　　　 　　　　　　　　　　　　　　（3）
　　HST3D对于离散方程的系数矩阵的求解方法包括：（1）三对角直接求解法（2）逐次超松弛迭代法（3）基于红黑排序的通用共轭梯度法（4）基于D4Z排序的通用共轭梯度法。
　　HST3D能够处理第一、第二和第三类边界条件，能够处理点源和面源问题，在三维网格坐标方向允许设置不同的土壤传导参数及容积参数。能够处理承压含水层问题，以及存在自由水面的潜水含水层问题。
　　我们选择HST3D作为问题计算分析工具，为进一步评价该程序应用于本问题的功能性和适用性，利用 HST3D对一典型近似工况下的双井承压含水层的全年温度场和流场进行了模拟。
　　
4 双井承压含水层模拟
　　
　　选取区域为长300m×宽200m×厚30m的具有上下不透水层的双井承压含水量水层为计算模型（见图1），中部承压含水层以及上下不透水层的厚度均为10m，相关水力热力参数见表1。模型区域中央为一抽水井和一回灌井，两井相距100m，两井均为完整井（透水井壁空越整个承压含水层）。模型的初始温度为15℃，ABCD、EFGH边界面为15℃恒温边界条件，ABFE、CDHG边界面为定水头边界条件，ADHE、BCGF边界为不透水边界。模型设定为无自然水头条件。
　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　图1 承压含不层计算模型示意图
　　
　　　　　　模拟模型的含水层水力热力参数 　　　　　　　　　　　　　　　　　　　　　　　表1

　　为了尽可能模拟水源热泵系统的全年"大温差，小流量"的运行工况，首先进行持续50天的夏季工况运行，抽水及回灌流量均为1200t/d，回灌水温恒为25℃；然后为持续50天的过渡季工况，两井停止运行；最后为持续50天的冬季运行，抽水及回灌井轮换角色，流量均为1200t/d，回灌水温恒为6℃。需要强调的是，由于HST3D程序的输入功能限制，本算例不得不采用固定流量和固定回灌温度。
　　该算例的模拟结果如下：
　　由图2可以看出，夏季工况期间，回灌热水锋面已经到达抽水井，出现"热贯通"现象，同时部分回灌热量以导热为主的方式向上下不透水层传递。由图3可以看出，过渡季工况期间，热量传递过程以导热为主，大部分夏季回灌热量在含水层内部实现"跨季节"储存。由图4可以看出，冬季工况期间，回灌冷水锋面同样已经到达抽水井，但是由于夏季回灌热水的存在，在一定程度上缓解了冷水锋面对于抽水井出水温度的影响。
　　　　　　
　　　　　　　　　　图2 25℃连续回灌50天（夏季工况）含水量水层中心剖面温度分布，A井回灌，B井抽水
　　　　　　
　　　　　　　　　　　图3 停止运行连续50天（过渡季工况）含水层中心剖面温度分布，A和B井停止运行
　　　　　　
　　　　　　　　　　图4 6℃冷水回灌连续50天（冬季工况）含水层中心剖面温度分布， A井抽水，B井回灌
　　
　　由图6可以看出，夏季工况期间，抽水井出水温度逐步上升，由原始含水层温度15℃升至约19℃；冬季工况运行开始阶段，由于夏季回灌热量在含水层中的蓄存，抽水井出水温度远高于含水层原始水温，达24℃，可以看出在冬季运行期间实现了部分夏季回灌热量的"热回收利用"。
　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　图5 双井运行期间含水层中心平面流速分布示意图
　　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　图6 冬夏季抽水井出水水温变化曲线
　　
5　HST3D适用性评价
　　
　　笔者认为，作为适用于深井回灌式水源热泵系统井群部分含水层水热运动的模拟分析的工具，需要能够实现以下几方面的功能：（1）能够处理非稳态问题，能够读入动态的边界条件参数（2）能够实现对于井群的参数设置（3）能够处理热边界条件（4）能够实现压力场（水头）的计算（包括单井水头和远端边界的水头影响）（5）能够反映符合工程实际的真实的物理过程，能实现地面系统与井群两部分的联合运行工况分析。
　　通过利用HST3D进行双井承压含水层的模拟，可以看出HST3D可以满足以上所列前4项功能，能够实现一定水文地质条件下含水层水热运动的非稳态模拟，能够给出在一定的井群运行工况下能量在含水层中转移、蓄存的过程分析，能够给出抽水井的动态温度变化。
　　但是HST3D还不能够实现地面以上部分（热泵机组）和地面以下部分（井群及共周边土壤）的全年联合运行分析，也就是说不能够根据建筑物全年动态负荷变化判断回灌温度和水量，从而对于符合工程实际上的过程进行模拟。由于HST3D具有一定的开放性和可拓展性，进一步的工作将为添加相关的建筑物负荷模块，完善其作为深井回灌式水源热泵系统井群运行的地下含水层传、蓄热性能模拟研究的计算分析工具。