湛江湾海域温排水温升分布及热环境容量模拟

孙 琰，许源兴，韩泽文，李 卓，何露雪，徐广珺,3，刘大召,3

（1.广东海洋大学电子与信息工程学院，广东 湛江 524088；
2.中海油信息科技有限公司天津分公司，天津 300450；
3.广东省海洋遥感与信息技术工程技术中心，广东 湛江 524088）

随着我国区域经济的快速发展，各行各业的用电需求与日俱增，近年来沿海地区大力兴建电厂，使得河口海岸地区的电厂布局紧张[1]。多数电厂采用直流冷却系统，大量冷却水直接流入海洋，造成受纳水体温度上升，改变水体环境质量，从而产生不可忽视的潜在生态环境问题。

在对电厂排热影响范围及程度的研究中，数值模拟作为常用方法之一得到了广泛发展与应用，通过计算得到海水流场与高温水的分布特征，从而对温排水进入水体后的漂移扩散规律进行定量描述。McGuirk 等[2]采用二维深度平均模型计算得到温排水在排放口附近的温升分布；
Hamrick 等[3]利用ECDC 建立三维数值模型模拟Peach Bottom 核电站温排水对附近水域水温的影响；
汪求顺[4]等通过建立温排水数值模型，分析杭州湾温排水的热扩散过程，总结该处海域潮差与温升变化的关系；
白玉川等[5]结合渤海湾围填海模拟结果，分析不同围填海条件下渤海湾的温排水造成的温升扩散分布规律；
李绍武等[6]采用MIKE 数值模型与杜比公式结合的计算方法，分析防波堤对电厂附近海域温升分布的影响。在温排水漂移扩散规律的研究基础方面，有学者通过数值模型结合周边环境敏感区的特点，预测分析电厂排热口不同布置方式对温升范围的影响，对不同排水方式与排放位置进行优化筛选，给出最优排放方案与选址建议[7-8]。此外，数值模拟方法也应用于环境容量的研究中，通过综合考虑受纳海域水质标准、近岸生态环境等因素计算出受污染水域环境容量总量与剩余环境容量，为海滨电厂排热量确定、水质目标管理与水质功能区划提供技术支持[9-12]。

本研究基于Delft3D 数值模型，综合模拟湛江湾内现存的多处电厂温排水造成的温升与热量积累效应，利用试算法计算各排热口附近海域的热环境容量总量与剩余热环境容量，对湛江湾海域电厂热排放进行全面统筹规划，为温排水指标的确定提供科学依据。

1.1 研究区自然条件概况

研究区为湛江湾（图1），湛江湾西北临湛江市，南北长约15 km，东西宽约24 km，属于半封闭型海湾，其入口小、内腹大，南西北三面由东海岛、湛江市、南三岛环绕组成，是一狭长天然良好水域[13]。湛江湾纳潮面积约264.1 km2，港湾纳潮量约5.4×108m3，湾内径流弱、潮流强。湛江湾海域的潮汐属于不正规半日潮，两次高低潮的潮差较大，根据湛江港多年实测资料统计，年平均高潮位3.04 m，平均低潮位0.87 m，多年平均潮差为2.16 m，平均涨潮与落潮历时分别为7 h 与5.5 h，涨潮历时大于落潮历时，落潮时的流速大于涨潮时流速。

图1 电厂、验证点位置示意Fig.1 Geographical location of power plants and verification points

1.2 湛江湾海域环境功能区及温排水排放分布

湛江湾海域近岸建有多处港口、码头、养殖区且存有多个海洋环境保护区，湾内的水质条件对生态环境产生重要的影响[14]。对污染物排放的限制主要依据，是海洋功能区划与所对应的海水水质标准。依照《广东省海洋功能区划(2011‒2020)》对湛江湾的功能区划与所执行的水质标准，湛江港保留区按规划要求为“保留现状”，根据广东海洋大学2020 年9 月所开展湛江湾水质调查航次中16 个站点的水质实测数据，分析化验表明该区为四类水体，后续将按照四类水质标准计算）。根据《海水水质标准》(GB3097-1997)中规定，对于第一、二类海域海水水质，人为造成的海水温升夏季不超过当时当地1 ℃，其他季节不超过2 ℃；
对于第三、四类海域海水水质，人为造成的海水温升不超过当时当地4 ℃的标准确定各湛江湾各功能区所对应的水体温升限制。根据收集的湛江市规划与湛江湾海洋工程环评报告等[15-17]，湛江湾内电厂温排水排放位置、排放方式、排放流量与温升数据见表1，其中京信电厂采用深水排放方式，其余两电厂排热口分别设在浅层和表层，处于湛江湾内湾的湛江电厂，排放流量最大，达到142 m3/s，三个电厂温升在8～9 ℃之间。湛江湾作为吸收三处电厂余热的受纳水体，水体的温升对湛江湾生物种群的生长繁殖与群落结构产生一定的影响。针对湛江湾的电厂排热布局现状，计算研究区域海域的温升分布以及热环境容量，有利于生态环境的保护与资源的合理利用。

表1 各电厂的温排水排放数据Table 1 Data of thermal discharge of the three power plants

1.3 Delft3D数学模型

本研究采用Delft3D 模型的FLOW 模块，计算模拟湛江湾海域水动力环境。Delft3D 是一个可用于河口、港湾及海洋的数学模型，对于潮流场、污染物输移扩散的模拟已被应用于许多研究中[18-21]。FLOW 模块在Navier-Stokes 方程的基础上，基于静水压强假定和Boussinesq 近似，采用循环隐式进程ADI（Alternating Direction Implicit）方法进行离散求解非线性浅水方程。在正交曲线坐标系(ξ,η)下，沿水深积分的连续性方程为：

其中，ζ为潮位，d为水深，t为时间为曲线坐标系与直角坐标系之间转换系数；
u、v分别为ξ和η方向上的速度分量；
Q为单位面积的源或汇流量。

模型在ξ和η方向的动量方程为：

其中，f为科氏力；
ρ0为水体密度；
Pξ和Pη为ξ和η方向水压力梯度；
Fξ和Fη为ξ和η方向紊动量通量；
Mξ和Mη为源（汇）项的动量分量。

热量输运方程：

其中，c为温度，ω为垂向速度，σ为垂向坐标，DH与DV为水平和垂向扩散系数，λd为一阶衰减过程。模型中通过对k-ε紊流模式，基于紊动扩散原理与实测值，对扩散系数进行校正[22]，最终得出DH取1～10 m2/s，DV取0.001 m2/s。S为单位面积的源和汇项，其包含水的排放或抽出以及自由表面的热量交换：

其中，qin与qout分别为单位排放源与单位取水源，单位为1/s，cin为排放温度，Qtot为通过水体自由表面与大气的热量交换，有些研究中采用综合散热系数来表征水气之间的热量交换过程，但使用经验公式计算后精度有限[23]，本研究选用更为精确的海洋热通量模型，通过引入太阳的垂直辐射强度、相对湿度与云覆盖率，计算背景辐射、蒸发以及对流引起的热量损失，此外，当空气与水体温度发生自由对流时，该模型也计算了潜热与感热的自由对流[24-25]，计算公式为：

其中，Qsn为短波辐射，Qan为长波辐射，Qbr为反向散射，Qev为蒸发热通量，Qco为对流热通量。研究所用气象数据均来源于NCEP的逐月再分析数据。

1.4 模型网格与参数

基于Delft3D-FLOW 模块建立的模型的计算域，包括湛江湾及其近海水域约5.65 × 103km2，其东北-西南方向约为116.2 km，西北-东南方向约为68 km。为了模型边界处的计算稳定，外海开边界设置为弧形。模型采用正交曲线网格，根据研究区海域的地形与水动力环境特征划分网格，使得网格方向与岸线方向总体垂直，并在湛江湾内重点关注水域进行网格加密处理，整个计算域包含10 246 个网格与6 406 个网格节点，网格分辨率为30～1 800 m。模型垂向分为10 层，从表层至底层每层占10%比例。外海开边界采用俄勒冈州立大学中国海潮汐模型所提供的8 个主要分潮（M2，S2，N2，K2，K1，O1，P1，Q1）调和常数驱动潮位。模型采用冷启动，初始温度场设置为湛江湾夏季海水平均温度28.9 ℃，京信电厂、湛江钢铁电厂与湛江电厂排热温度分别为36.9、37.9、36.9 ℃，按照电厂的实际排热深度，在模型中相对应层深设置排放源。水陆边界采用绝热条件。经过率定计算域底床糙率取0.015～0.022，水平与垂向紊动黏滞系数分别取10 m2/s 与1×10-6m2/s。计算时间步长为60 s，计算时段为7～9 月，模型输出的时间分辨率为1 h。模型的计算域网格见图2。

图2 模型计算域网格Fig.2 Model grids of the calculation area

2.1 模型验证

为验证Delft3D模型的模拟效果，将模型计算的潮位、潮流与实测数据进行对比验证，验证模拟时段分为大潮、小潮两个阶段，大潮验证时间为2016年3月8 日12 时至3 月9 日12 时，小潮验证时间为2016年3月16日16时至3月17日16时，研究区内L1、L2两个潮位站点与S1、S2 两个潮流站点分布见图1。图3和图4分别表示大、小潮期间的L1、L2点的潮位验证结果和大潮期间S1、S2潮流验证结果。

图3 潮位验证Fig.3 Verification of tidal level

图4 大潮流速流向验证（S1、S2）Fig.4 Verification of tidal current speed and direction during spring tide(S1、S2)

为了进一步对模型计算结果评估，采用平均绝对误差(MAE)、均方根误差（RMSE）、相关系数（r）以及Willmott[26]提出的统计学方法模型技术评分(Skill)进行定量评价。

模型技术评分：

其中，Xi为模拟值，Yi为观测值，N为时间序列为模拟值的平均值为观测值的平均值。MAE、RMSE 值越接近于0，代表模型计算的效果越好；
r越接近于1，代表模拟计算的值与观测值越接近；
模型技术评分Skill>0.65，表示模型比较符合实测值，Skill的值越低，表明计算与实测越不相符。

在验证时间段内，L1、L2 潮位验证的MAE 为0.023、0.05，RMSE 为0.26、0.19，r为0.98、0.99，模型技术评分为0.98、0.99；
S1、S2的流速验证的MAE为0.003、0.03，RMSE 为0.11、0.12，r为0.80、0.82，模型技术评分为0.89、0.90；
S1、S2 流向验证的MAE 为11.37、23.52，RMSE 为59.82、35.50，r为0.79、0.94，模型技术评分为0.89、0.95。

综上所述，本研究建立的潮流数学模型，在潮位、潮流方面实测与模拟计算值拟合度较好，模拟计算值与实测值数据的变化趋于一致，可以认为模型能较准确地反映湛江湾海域潮流环境。

2.2 模拟结果分析

2.2.1 潮流场模拟分析 图5 为模型计算的湛江湾及其近海海域，涨急流场与落急流场的分布图，通过计算分析发现该海域的涨、落潮方向大致与东海岛、南三岛东侧岸线垂直，涨潮时来自外海的潮波通过模型开边界进入计算域，海水大致呈东南向西北方向传播，落潮时潮流沿原路返回，流向外海。

图5 计算域涨急与落急流场分布Fig.5 The current of the calculation area at the moment of fastest flood and fastest ebb

湛江湾内海域海流受潮汐作用较为明显，在湛江湾独特狭长地形影响下呈现明显往复流特征，涨、落潮期间海水从湛江湾口航道汇入或流出，涨潮时海流从湛江湾口以及南三岛北侧水道汇入湾内，在湾口到湾顶大致呈东南-西北流向；
落潮时湛江湾内北部水道与西侧海流在汇聚后一并向东南方向流出湾外。流速方面，湾口区域为流速高值区，其流速最大可达1.4 m/s；
海流进入湾内后，过水断面增大使得流速下降，基本表现为航道深水区平均流速大于其他海域，浅滩、岛屿以及近岸地区流速最小，总体上涨潮平均流速略小于落潮平均流速。

2.2.2 温升平面分布 目前我国对温排水的监管与控制，多采用“一事一议”的方法来确定每一个滨海电厂温排水对环境的影响，这就可能忽略了多个排热口之间的互相叠加影响，本研究综合三个电厂的排热数据，计算得出湛江湾内温排水造成的温升分布，主要集中在东海岛北侧与湛江湾内湾河道，图6给出湛江湾大小潮期间涨憩与落憩四个特征时刻，湛江湾海域表层的温升分布图。

在东海岛北侧海域往复流的带动下，在此的京信电厂温升区沿岸线呈带状分布，京信电厂的排热对其西侧的东头山岛产生了不同程度的影响，并产生一定的热量积累，在大潮涨憩时对东头山岛影响的范围最广（图6(a)）。

位于湛江湾口的湛江钢铁电厂排热口附近海域，大潮时潮动力强劲，温排水得到充分掺混并且能随着涨落潮更远，这使得钢铁电厂与京信电厂之间的热量积累更加明显，两电厂之间在大潮时产生的温升大于0.5 ℃。在小潮潮型时，水体掺混能力弱，高温水扩散距离短，0.5 ℃温升范围主要集中在排热口附近（图6(c-d)）。

图6 特征时刻温升分布（温升：0.5～4 ℃）Fig.6 Temperature rise distribution at characteristic time（Temperature rise:0.5～4 ℃）

湛江湾上游的湛江电厂，排热量较大且为浅层排放，一个潮周期内温升带分布变化规律比较明显。涨潮时温排水沿河道向北转移至浅水区，0.5 ℃温升线最远到达湾顶，在落潮时温排水随潮流向南转移，在大潮落憩时0.5 ℃温升线向南最远扩散约14.6 km（图6(b)）。值得注意的是湾顶区域属五里山港海洋保护区，其执行温升不超过1 ℃的水质标准，但在大潮涨憩时刻2 ℃温升线深入保护区约1.4 km（图6(a)），已远超出限定功能区划水质限制。

在温升面积方面，温排水在各排热口附近产生了不同程度的影响，不同潮型0.5 ℃的温升包络面积大潮略大于小潮，大潮时涨、落憩的0.5 ℃温升包络面积分别为73.3、73.21 km2，小潮涨、落憩0.5 ℃温升包络面积分别为72.6、67.7 km2；
3.5 ℃的温升包络面积小潮大于大潮：大潮时涨、落憩的3.5 ℃温升包络面积分别为0.96、1.4 km2，小潮涨、落憩3.5 ℃温升包络面积分别为4.1、3.6 km2。总体来看，大潮高温升（3.5 ℃）面积小，低温升（0.5 ℃）面积大，而小潮时与之相反，这是大潮的水动力强于小潮，使得大潮时水体掺混程度较高造成的。

2.3 热环境容量计算

2.3.1 热环境容量计算方法 环境容量是环境的一种属性，指环境所能容纳某种特定活动或活动频率的能力。环境容量的大小即为特定海域自净能力强弱的指标，是海水自净能力综合表现的定量描述[27]。

温排水进入受纳水体后，与环境水体发生掺混稀释，并随潮流输移扩散将热量带向远区。热污染相较于其他污染物质没有浓度的概念，所以在热环境容量的计算中，对热污染更重视的是温升范围的控制。试算法是目前计算热环境容量的较为常用方法。假设某污染源的排热量W（单位：m3·℃/s），计算其全潮最大温升分布，然后根据温升影响范围是否满足环境评价指标来调节W，如温升影响符合环境标准，则可以继续增大W，直至满足评价指标的最大热排放量W0，即为热环境容量。根据热环境容量的计算结果，可以对附近海域的水热环境进行评估，从而达到热排放总量控制的目的。本研究根据《海水水质标准》和《广东省海洋功能区划(2011‒2020)》的相关规定，确定湛江湾外的湛江—珠海近海农渔业区执行温升不超过1 ℃的一类水质标准；
东海岛北侧工业城镇用海区、湛江港保留区与内湾水道执行温升不超过4 ℃的四类水质标准；
特呈岛海洋保护区与内湾水道北侧五里山港海洋保护区执行温升不超过1 ℃的二类水质标准。

此外，一般污染物在刚排出到受纳水体当中时，会使得一定范围内水体温度显著增加，在此范围外污染物受到混合、稀释，浓度迅速降低，以此范围内的水体温度来研究热环境容量是不准确的，这就需要引入“混合区”这一概念，来代表排热口附近的允许水质超标的区域。混合区设定需要全面综合的考虑，若混合区设定范围过大，会忽略一定的水体处理能力；
若范围设定过小，则将会过度限制污染物的排放[28]。针对湛江湾水域情况，根据《污水海洋处理工程污染控制标准（GB18486—2001）》：污水排放到面积小于600 km2海湾时，混合区范围面积A由此来判定：

其中，A0为计算至港口位置的海湾面积，单位为m2。

湛江湾模型计算域面积约为269.13 km2，小于600 km2，故适用式（10）计算混合区面积。通过计算得出混合区面积约为1.35 km2，混合区半径约为1 162 m。但是由于湛江湾内湾水道较窄，1 162 m的混合区半径可能占据整个水道，又根据国家环境保护总局《近岸海域环境功能区管理办法》，确定混合区范围时要综合考虑当地的水动力条件等，湛江湾内湾河道潮流、水交换情况与湛江湾其他海域相比较为独特，故本研究在内湾河道不设置稀释混合区。

2.3.2 热环境容量总量计算 本研究以三个电厂现有热排放为起算点，采用模型试算法逐步改变W，直至温升包络线边界处的温升达到水质限制标准，即为所求热环境容量。通过36 次试算得到湛江电厂排热口附近海域热环境容量为2 952(m3·℃)/s，湛江钢铁电厂排热口附近海域热环境容量为1 870(m3·℃)/s，东海岛北侧京信电厂排热口附近海域热环境容量为4 745 (m3·℃)/s。东海岛北侧京信电厂附近海域潮动力较强，且四类水质标准的功能要求使得该处海域的热环境容量最大；
而位于京信电厂东侧湛江湾口的湛江钢铁电厂排热口，虽然附近海域水交换率大，自净能力高，但同时要满足湾外湛江-珠海近海农渔业区所执行的一类水质标准，使得该处海域的热环境容量受到了极大限制；
位于内湾上游的湛江电厂附近海域水交换能力弱且在北侧五里山港海洋保护区二类水质标准的制约下总热环境容量最小。

2.3.3 剩余热环境容量 环境容量总量与现状污染物排海通量的差，即为剩余环境容量。当剩余环境容量大于0 时，代表区域海域理论上还能容纳更多该种污染物；
当剩余环境容量小于等于0时，表示污染物当下排放通量已经超过允许排放量，应当削减排放。结合表1 温排水排海流量，总热环境容量减去各电厂已使用热环境容量，得出三处电厂排热口附近海域的剩余热环境容量（图7）。其中东海岛北侧京信电厂排热量占用其附近海域热环境容量的33.6%，剩余容量为66.4%；
湛江钢铁电厂现有排放占用其附近海域62.8%，还有37.2%的剩余环境容量；
湛江电厂附近海域已使用的热环境容量已超出其最大容量，剩余热环境容量为-77.6%。综上所述，京信电厂附近海域环境容量还有较大剩余；
湛江钢铁电厂近海热环境容量总量较小，在目前排温情况下环境余量较少，后续开发利用需注意；
湛江湾内湾上游的湛江电厂的温排水排热量已经远超过该海域的热环境容量，在其附近海域容易产生热量积累，需要改进排污工艺流程或减排才能减少热污染对生态环境造成的损害。

图7 各电厂温排水热环境容量对比Fig.7 Comparison of heat environmental capacity of each power plant

本研究针对湛江湾潮汐水域多个电厂高温水排海现状，应用Delft3D 水动力与热通量模块，结合电厂排水工程的实际排放工况，建立能反映电厂温排水热力特性的三维数学模型，绘制了特征时刻的温升分布图，并结合现有规范区划，使用试算法计算了热环境容量，综合分析得到以下结论：

1）三处电厂排热口附近均有一定的热量积累，3.5 ℃的温升包络面积在大潮涨、落憩时分别为0.96、1.4 km2；
在小潮涨、落憩时分别为4.1、3.6 km2。

2）潮动力对半封闭海湾湛江湾中温排水的扩散与分布起主导作用。在大潮时高温升（3.5 ℃）面积小，低温升（0.5 ℃）面积大，小潮时反之。

3）在现有的排放布局下，湛江湾内湾水道北端的五里山港海洋保护区，在大潮涨憩时刻2 ℃温升线深入保护区约1.4 km，超出了水质限定标准；
东海岛北侧的东头山岛周边海域受温排水作用明显，有一定的热量积累。

4）计算得出湛江湾各电厂附近海域的热环境容量：京信电厂>湛江电厂>湛江钢铁电厂，分别为4 745、2 952和1 870(m3·℃)/s。进一步计算各电厂附近海域的剩余热环境容量，京信电厂附近海域剩余热环境容量最大（66.4%），而湛江电厂的热排放已经超标77.6%，对该海域生态环境产生了热污染影响。

因此，湛江电厂可通过削减排热量，或通过增加通风冷却塔、循环水池等优改进排水工艺，降低电厂排热对海域生态环境的影响；
东海岛沿岸两电厂热环境容量仍有剩余，可适当加以开发利用，但要注意排热量过大会使两电厂热污染产生叠加，避免温升带在潮流作用下到达特呈岛海洋二类水质保护区和湛江湾外的湛江-珠海近海农渔业一类水质区。