Page 123 - 2022年第53卷第5期
P. 123
３ 数值模型
３．１ 概念模型 以工程场地为中心向四周各扩展 １０００ｍ作为计算边界，模型底部以第Ⅲ承压含水层
底部作为底部边界。根据土体的性质，模型垂向上共剖分 ５层，分别为潜水含水层、第Ⅰ弱透水层、
第Ⅰ承压含水层、第Ⅲ弱透水层和第Ⅲ承压含水层。由于计算边界与抽、灌井距离远，结合南通地区
水文地质条件、渗透系数参考值和影响半径经验公式，地下水源热泵系统的运行过程中抽、灌地下水
影响范围小于 １０００ｍ，地 下水源 热泵 系统 运 行对数值模 型边 界影 响 很 小。因 此，四周 边界 处 理 为
Ｄｉｒｉｃｈｌｅｔ边界，即水力学上处理为定水头边界、热力学上处理为恒温边界、水化学上处理为恒定浓度
边界；根据计算域水力学特征，模型顶面处理为自由面边界，底部处理为 Ｎｅｕｍａｎｎ边界，即水力学上
处理为隔水边界、热力学上处理为热通量边界、水化学上处理为零通量边界。
本项目利用第Ⅰ承压含水层的浅层地热能，需要抽取、回灌第Ⅰ承压含水层的地下水，由于各层
存在水位差，各含水层之间存在水力联系，地下水流特征概化为三维非稳定流。地下水源热泵系统运
行过程中，假定每个节点的温度同时代表地下水和土颗粒温度，忽略地下水密度随温度发生的变化。
３．２ 数学模型 各含水层之间的黏性土层作为弱透水层和含水层一起直接参与计算，地下水流为三维
非稳定流。取坐标轴方向与各向异性介质主方向一致，地下水系统渗流控制方程如下 ［１５］ ：
ｈ
ｈ
ｈ
( ) ( ) ( ) ＋ Ｗ＝ Ｓ ｈ （１）
＋
＋
Ｋ
Ｋ
Ｋ
ｘ ｘｘ ｘ ｙ ｙｙ ｙ ｚ ｚｚ ｚ ｓ ｔ
式中：Ｋ 、Ｋ 、Ｋ 为各方向水的渗透系数，ｍ?ｄ；ｈ为地下水位，ｍ；Ｗ为源汇项，１?ｄ；Ｓ 为含水层
ｘｘ ｙｙ ｚｚ ｓ
贮水率，１?ｍ。
假设地下水和含水介质骨架的热动平衡是瞬时完成的，即含水介质骨架与周围地下水具有相同的
温度，并忽略由于温度差引起水的密度不一样而引起的上下自然对流的影响。地下水系统热量运移控
制方程如下 ［１６ － １８］ ：
（ｎｕＴ） （ｎｕＴ） （ｎｕＴ）
Ｔ
Ｔ
Ｔ
( ) ( ) ( ) [ ｘ ｙ ｚ ] ＋ Ｑ ＝ ｃ Ｔ
＋
＋
－ ｃ
ｘ λ ｘ ｘ ｙ λ ｙ ｙ ｚ λ ｚ ｚ ｗ ｘ ＋ ｙ ＋ ｚ ｃ ｓ ｔ （２）
ｕｕ
＝ ＝ ｕ ＋ （ － ） ｉ ｊ （３）
ｉｊ
λ ｉｊ λ ′ ＋ λ ｐ β Ｔ δ ｉｊ β Ｌ β Ｔ ＋ λ ｐ
ｕ
为各方向水的热动力弥散系数，Ｗ?（ｍ·℃）；Ｔ为岩土体温度，℃；ｃ、ｃ分别为水
式中：λ ｘ 、λ ｙ 、λ ｚ ｗ ｓ
３
３
和土体的体积热容量，Ｊ?（ｍ·℃）；ｕ为多孔介质中地下水平均流速，ｍ?ｄ；Ｑ为热的源汇项，Ｊ?（ｍ·ｄ）；
ｃ
为 Ｋｒｏｎｅｃｋｅｒ记号；ｕ、ｕ为地下水平均流速分量，ｍ?ｄ；
β Ｌ 、β Ｔ 分别为纵向和横向热弥散度，ｍ；δ ｉｊ ｉ ｊ
分别为热机械弥散系数和热传导系数，Ｗ?（ｍ·℃）。
ｉｊ
λ ′、λ ｐ
根据多孔介质流体动力学，溶质迁移的过程中满足质量守恒定律，本次研究针对稀物质在水中迁
移过程，可假定流体 的密 度为 定 值。取 坐 标 轴 方 向 与 主 渗 透 方 向 一 致，则 地 下 水 溶 质 运 移 控 制 方
程为 ［１９ － ２１］ ：
（ｕｃ） （ｕｃ） （ｕｃ）
ｃ
ｃ
ｃ
( ) ( ) ( ) [ ｘ ｙ ｚ ] ＋ Ｉ ＝ ｃ （４）
Ｄ
Ｄ
＋
Ｄ
－
＋
ｘ ｘｘ ｘ ｙ ｙｙ ｙ ｚ ｚｚ ｚ ｘ ＋ ｙ ＋ ｚ ｔ
ｕｕ
ｉ ｊ
ｕ ＋ （ － ） ＋ Ｄ （５）
ｄ
ｉｊ
Ｄ ＝ Ｄ′ ＋ Ｄ ＝ α Ｔ δ ｉｊ α Ｌ α Ｔ ｄ
ｉｊ
ｕ
２
式中：Ｄ 、Ｄ 、Ｄ 为各方向的水动力弥散系数分量，ｍ ?ｄ；ｃ为地下水中溶质的质量浓度，ｍｇ?Ｌ；
ｚｚ
ｘｘ
ｙｙ
分别为纵向和横向溶质弥散度，ｍ；Ｄ′为机械弥散
ｉｊ
ｎ为孔隙度；Ｉ为溶质源汇项，ｍｇ?（Ｌ·ｄ）；α Ｌ 、α Ｔ
２
２
系数，ｍ ?ｄ；Ｄ 为分子有效扩散系数，ｍ ?ｄ。
ｄ
研究表明：当地下水流速较快时，机械弥散系数远大于分子有效扩散系数，因此忽略分子扩散系
数的影响 ［２２］ 。
地下水平均流速分量分别表示为：
— ６ ２ ３ —
