淡澳河口——大亚湾三维潮流及海水入侵模拟

谢培，张亚辉，乔飞*

1.环境基准与风险评估国家重点实验室，中国环境科学研究院

2.中国环境科学研究院环境分析技术测试中心

河流是陆源物质和污染物质进入海洋的重要通道，陆地进入海洋的物质85%由河流搬运[1]。河口区是水动力复杂的活跃水文生态系统和陆海相互作用的集中地带，由于大多沿岸河口区人类活动频繁，河口水体储存了大量污染物质，在涨落潮不稳定的环境特性下，水体污染物质不断振荡难以扩散。

大亚湾是位于珠江口东侧的亚热带半封闭性海湾，近几十年随着沿岸人口增长、经济及相关产业的快速发展，大亚湾水体环境富营养化问题突出[2]。在大亚湾12条入海河流中，最大的河流输入是西北部的淡澳河，已有研究表明，大亚湾入海河流中淡澳河入海污染物通量占比最高[3]。淡澳河为感潮河段，水流除受上游径流影响外，还受海洋潮流的强烈影响，不同潮流特性对河口区影响差异大，河口内段受边界约束潮波变形，因此明确该区域的潮流场及水动力特性，对于控制水环境污染具有极其重要的意义。

目前数值模拟是较为普遍且有效的环境建模方法，由于天然水体是三维的，三维数值模拟能全面反映水环境系统流场的基本特征。在三维数值模拟模型通用性上，国内鲜有模型能够实现复杂边界准确拟合以及水流分层可视化精细模拟，当前应用最为普遍的三维数值模拟模型包括美国的环境流体动力学（EFDC）模型和FVCOM模型、丹麦的Mike系列模型以及荷兰的Delft3D模型等[4]。陈青毅等[5]利用三维数值模拟河口水库的水流流动和水体交换，制定了合理的水库引排水方案；
姜恒志[6]应用EFDC模型建立了太湖三维水动力数学模型，分析了太湖流场特性；
齐珺等[7]建立三维水动力模型考察了长江水系武汉段河床的冲淤变化。尽管国内外学者已在典型河流湖库就水动力特征开展了一些研究[5-8]，但鲜有对淡澳河下游感潮河段、河口及其毗邻海区水体交换的整体探究，有关淡澳河流域水环境问题探究多集中于海水入侵带来的灾害，如土壤盐渍化[9]，以及陆地-海洋相互作用对地下水排泄的影响[10-11]等方面。为系统了解该水域复杂的水动力和咸淡水交互作用，本研究构建了淡澳河下游感潮河段、河口及其毗邻海区为一体的三维水动力水质模型，旨在揭示径流潮流作用下的水动力过程和内在影响规律。此外，考虑到充分认识不同径流量与潮流交互作用引发河口水环境的变化，以系统全面地了解淡澳河口特性，笔者定量估算了不同水期盐水入侵特征、咸淡水交互与水质的响应关系，以期为淡澳河流域水环境治理提供依据。

1.1 研究区域

淡澳河为分洪渠与响水河汇合后至白寿湾出海口段，于 114°33"57"E、22°43"29"N 河口处汇入大亚湾，是大亚湾入海河流中流量最大且水质最差的河流。淡澳河河长约14 km，下游河口处宽度为500 m，上游河宽约为200 m，为较典型的上窄下宽喇叭状河口。淡澳河河段由于海水的顶托作用，污染物不易向外扩散，河段水质较差。淡澳河上游的淡澳分洪渠为人工分洪渠，早期主要为缓解淡水河及其下游的洪水压力，近年来随着惠阳区的发展，大量工业废水和生活污水通过淡澳分洪渠和淡澳河进入大亚湾。

1.2 研究方法

1.2.1 模型选择

模拟区域包含淡澳河感潮河段及整个河口区，河道与河口在横向尺度上存在很大的差异，并且区域水动力过程复杂，对模型适用性的要求较高。EFDC模型垂向上采用σ坐标变换，能较好拟合近岸复杂的岸线和地形，在河流、水库、湖泊、河口等多类水体都得到广泛应用，并体现了很好的适用性。由于淡澳河口水动力条件复杂，采用长期应用和改进的河口海岸三维数值模型，该模型已长期应用于河口水动力过程和盐水入侵等方向的研究，并取得了众多成果[6-8,12-17]。

1.2.2 模型构建与验证

1.2.2.1 网格设置

模型采用水平曲线非正交网格，范围包括整个淡澳河河口、大亚湾和临近海区，上游流量边界设定为响水河与淡澳河交汇处上游约600 m处，下游为外海开边界。由于模拟范围较广，河道、河口和海区各区域水体空间尺度变化较大，为尽可能满足不同区域计算精度同时保证较高的计算效率，模型采用变尺度网格系统，河口内网格分辨率为15 m×15 m～40 m×80 m不等，河口外网格较稀疏，分辨率为80 m×80 m～125 m×200 m不等。垂向采用σ坐标，均匀分为5层。考虑到淡澳河河口区域浅滩较多，模型选用动边界处理，运用干湿判别法处理潮滩移动边界，临界水深取0.2 m。

1.2.2.2 条件设置

模型上游设置为流量型边界，下游设置为水位边界，模型上下游边界、姚田桥和虎爪断面（国控点位）位置分布如图1所示。本研究模拟2019年1月1日——12月31日淡澳河口——大亚湾区水动力水质过程，上边界采用日流量过程，下边界由潮位驱动，其潮流特性主要受来自太平洋的潮波制约，潮位数据资料（由各分潮调和常数合成得到）取自全球潮汐数值模式NAOTIDE计算结果（http://www.miz.nao.ac.jp/staffs/nao99/index_En.html）。模式岸线采用2019年海图数字化资料及2019年水深实测资料。水质指标主要选择COD、氨氮、TP等，水质浓度采用实测值。

图1 网格边界条件分布Fig.1 Grid division and boundary conditions of research region

此外，为提高模型计算效率，本次模拟采用动步长，基准步长设置为0.1 s，实际步长由模型根据实时克朗数自动调整，经统计平均时间步长约为1 s，平水期和枯水期稍长，丰水期稍短。

1.3 参数率定验证

1.3.1 主要参数设置

模型设置水平涡黏系数为Ax=Ay=100 m2/s，扩散系数为1×10-5m2/s，粗糙高度是表征河流底板对水体阻力的变量，相当于曼宁公式中的糙率，采用单元划分法[18]求解糙率，最终取值为0.001～0.01 m。

1.3.2 模型验证

为保证模型的可靠性和适用性，水动力水质模型验证采用3组不同水期资料：1月（枯水期）、6月（丰水期）和10月（平水期）。主要水动力水质参数验证站位为虎爪断面，计算结果表明模拟值与实测值均呈现较好的相关性（表1）。水位和盐度模拟与实测值表现为强相关，相对误差小于5%，主要水质指标中氨氮相对误差小于10%，COD和TP相对误差略高（在可接受范围之内），主要原因可能在于COD和TP输入量统计不够，导致局部断面模拟浓度出现一定的误差。

表1 水动力水质参数模拟误差统计Table 1 Error analysis of simulated hydrodynamic water quality parameters

淡澳河水体模拟纵向分为5层，分别提取底层（k=1）、中部（k=3）和表层（k=5）盐度数据，通过不同水层含盐量计算淡水与盐水比例。淡澳河上游含盐量约为0‰，中游含盐量约为15‰，下游含盐量约为35‰，河流入海区盐度高于河流上中游盐度，符合入海河流盐度分布规律。以淡澳河虎爪断面月内不同水层盐度分布为例（图2），计算得到该断面表层、中部和底层的淡盐水比例。水体底层含盐量约为10.8‰，淡盐水比例约为2∶1；
水体中部含盐量约为6‰，淡盐水比例约为7∶1；
水体表层含盐量约为1.7‰，淡盐水比例约为 16∶1。

图2 水体垂向盐度分配Fig.2 Distribution of vertical salinity of water body

2.1 潮流场

大亚湾潮流性质以不规则半日混合潮型为主，由于受地形影响，外海潮波传至大亚湾内变形较大，以致潮汐日现象显著，潮波主导着水流运动，潮流动力作用较弱。湾内潮流基本是顺水道（主航道）的往复流[19]，涨潮平均流速为 0.065 m/s，落潮为0.073 m/s，落潮流速比涨潮流速大，与吴仁豪等[20]的研究结果一致。水动力模拟结果显示，淡澳河虎爪断面的表层平均流速为0.153 m/s，底层平均流速为0.094 m/s，表层流速约为底层流速的1.6倍。涨潮流向以NW（320°）为主，落潮流向以 SE（100°）为主，海浪以涌浪为主，风浪为辅，平均波高为0.2 m，常浪向为SE。港湾内缺少长期波浪观测资料，短期观测资料显示，6级风时，港内水域波高仅0.6 m，12级台风时，波高仅1.3 m。

大亚湾月内潮汐变化如图3所示。从图3可以看出，每月8～10 d为日潮，20～22 d为半日潮。不同潮型日潮（17日）和半日潮（11日）的潮流场分布如图4所示。日潮涨急流速为0～0.963 m/s，其中表层涨急最大流速为1.019 m/s，底层涨急最大流速为0.878 m/s；
半日潮涨急流速约为0～0.918 m/s，其中表层涨急最大流速为1.003 m/s，底层涨急最大流速为0.767 m/s，不同水层均表现为日潮涨急流速大于半日潮涨急流速。

图3 大亚湾潮汐变化趋势Fig.3 Trend of tide change in Daya Bay

图4 日潮和半日潮涨急流域潮流场分布Fig.4 Distribution of tidal current at diurnal and semidiurnal tides during typical moment

考虑到不同潮汐类型对河流水质的影响，本研究分析了不同潮汐类型下的潮流场分布，由图5可知，日潮高潮时，表层海水入侵到虎爪断面下游约400 m处，底层海水入侵到虎爪断面上游约1 100 m处；
半日潮高潮时，表层海水基本停留湾内，底层海水入侵到虎爪断面上游约700 m处。总体上表现为日潮和半日潮的高潮时刻底层海水可入侵到虎爪断面上游700～1 100 m处，但表层海水并不能到达虎爪断面上游。

图5 日潮和半日潮流域表、底层潮流场分布Fig.5 Distribution of tidal current at diurnal and semidiurnal tides on the surface and bottom

2.2 盐水入侵

通常河流淡水的盐度在0.5‰以下，河水与陆架海水的混合水体盐度为0.5‰～30‰，陆架海水的盐度为30‰以上，因为底层盐度梯度更大，界面混合影响小，可采用底层盐度等值线作为河海水团分界线。模拟得到淡澳河上游姚田桥至虎爪断面水体的底层盐度分布如图6所示。由图6可知，0.5‰盐度线分布在虎爪断面上游，距离姚田桥约1 490～1 520 m处，在河口海岸水环境整治时可依据分界线对淡水和海水进行分区域治理。此外，由于盐度分布受径潮流等影响，季节差异大，分界线较为复杂，但水文要素往往有周期性规律，研究基于2019年水文资料初步划分了河海水团分界线，考虑到分界线的稳定性和合理性可考虑取多年平均值来划分。

图6 淡澳河上游姚田桥至虎爪断面底层水体盐度分布Fig.6 Salinity distribution of bottom water from Yaotianqiao to Huzhao section in the upstream of Dan"ao River

淡澳河上、下游边界处不同水层的含盐量基本一致，其中上游边界处含盐量约为0‰，下游边界处含盐量约为35‰。中游地区由于海水入侵呈现明显的盐度分层流，研究以虎爪断面为代表断面，通过分析不同水期下水体表层和底层盐度分布（图7），确定水体淡盐水分配比例。

图7 中游河段不同水期水体盐度分布Fig.7 Distribution of water salinity in the midstream in different water periods

枯水期时，表层和底层含盐量约为10.4‰和20.7‰，经计算得到表层和底层咸淡水占比分别约为1∶2.5和10∶7；
平水期时，表层和底层含盐量约为1.7‰和10.8‰，计算得到表层和底层咸淡水占比分别约为1∶16和1∶2；
丰水期时，表层和底层含盐量约为0.3‰和5.0‰，计算得到表层和底层咸淡水占比分别约为1∶100和1∶6。由咸淡水占比可见，丰水期虎爪断面水环境质量主要由淡澳河上游来水主导，枯水期时，虎爪断面水质情况由淡澳河上游来水和海水共同作用，其中表层受上游来水影响大，底层受海水入侵影响大。

2.3 咸淡水交互对水质的影响

最新研究表明[21]，大亚湾海域营养盐绝大部分通过河流输入，其中氮通量输入9.4×105kg/a，磷通量输入1.7×105kg/a，淡澳河是大亚湾入海污染物输入的主要途径，淡澳河输入氮通量占76.7%，磷占82.4%。为进一步证明上游输入对淡澳河口产生了明显影响[22]，摸清上游来水流量与河流水质的响应关系，以2019年8月实测水文水质数据为基准，即上游来水流量为6.0 m3/s，同时依据实际水文情况设定枯水期小流量为1.0 m3/s和丰水期大流量为10.0 m3/s。

由模拟结果可知（图8），枯水期潮期间污染物浓度变化不明显，虎爪断面COD降幅为39%～42%，氨氮浓度降幅为27%～30%，TP浓度降幅为13%～15%；
丰水期小潮期间污染物浓度增加幅度最为明显，虎爪断面COD增幅为1%～12%，氨氮浓度增幅为11%～14%，TP浓度增幅为6%～7%。由以上数据可知，上游来水流量变化直接影响虎爪断面水质，流量越大对虎爪断面的水质影响越大，以小潮期间最为显著，建议进行相关工程优化调度，以减小小潮期间的下泄流量。

图8 上游来水变化对虎爪断面水质影响分析Fig.8 Analysis of the influence of upstream water on the water quality of Huzhao section

(1)大亚湾潮流性质以不规则半日混合潮型为主，每月约8～10 d为日潮，20～22 d为半日潮，落潮流速高于涨潮流速，日潮涨急流速略高于半日潮涨急流速。由不同潮汐类型对淡澳河潮流场的影响可知，日潮、半日潮的高潮时表层海水并不能到达虎爪断面，但底层海水可入侵到虎爪断面上游700～1 100 m处。

(2)基于2019年水文资料初步划分河海水团交界线，位于虎爪断面上游距姚田桥约1 490～1 520 m处，可为河口海岸带淡水和海水分类整治提供依据。盐水入侵表明，淡澳河水体含盐量枯水期＞平水期＞丰水期；
以虎爪断面为例，枯水期含盐量高达16‰，底层海水含量高于淡水，枯水期表层海水含盐量较平水期增加明显，相当于平水期的6倍。

(3)上游来水流量与水质响应模拟研究表明，来水流量直接影响虎爪断面水质，其中，枯水期污染物浓度降低15%～40%，而丰水期污染物浓度增加10%左右，以小潮期间增加幅度最为显著。基于此，可建议进行相关工程优化调度，以减小小潮期间的下泄流量。