基于数学模型的乐清湾环境容量算法应用

谢中宇，刘涵心，陈界仁

（1.国家海洋局第二海洋研究所工程海洋学重点实验室，杭州国海海洋工程勘测设计研究院，浙江杭州 310012；2.河海大学水利水电学院，江苏南京 210098）

基于数学模型的乐清湾环境容量算法应用

谢中宇1，刘涵心2，陈界仁2

（1.国家海洋局第二海洋研究所工程海洋学重点实验室，杭州国海海洋工程勘测设计研究院，浙江杭州 310012；2.河海大学水利水电学院，江苏南京 210098）

在实测资料分析的基础上建立了乐清湾潮流水环境数学模型,对乐清湾区域水动力环境及污染物扩散输移进行数值模拟，以COD为代表，计算分析了乐清湾在大麦屿港区规划中满足海洋功能区划水质标准下的COD环境容量。结果表明：乐清湾完成一次水交换时间为90 d，排污附近海区完成一次水交换为30 d。大潮期72 h排污口排放污染物COD影响范围大致为北至乐清湾底，西南至瓯江北口，南至洞头峡，东至鹿西岛以东15 km的海域，湾内COD浓度为1.231～1.232 mg/ L。允许混合区的COD剩余环境容量大潮期为8 148.99～8 164.98 t/a；小潮期为8 144.64～8 162.75 t/a。大麦屿港区规划实施后满足乐清湾COD环境容量要求。

环境容量；污染物扩散；水体交换；数学模型；乐清湾

乐清湾位于浙江省东南部沿海，瓯江口北侧，三面环陆属典型半封闭性海湾，岸线总长约184.7 km，海湾形似葫芦，口门宽约21 km，纵深达40 km，按地形地貌形态分为内湾、中湾和外湾三大部分[1]，岸线以下海湾总面积463.6 km2。乐清湾海区属半封闭性强潮海湾，控制本海区潮波运动的是以M2分潮为主的东海前进潮波系统，潮波由洞头洋经湾口诸海峡传入湾内。根据海域实测资料统计，湾口海域最高、最低潮位分别为4.30 m和-3.81 m；平均高、低潮位分别为2.77 m和-2.25 m，量值与坎门站比较接近。

由于土地资源开发、围塘养殖、航运等需要,台州港大麦屿港区规划于玉环半岛西侧岸段建设包括普竹、连屿、大岩头在内的三个围填海作业区及附属码头，填海面积375.1 hm2，新建码头泊位28个。为研究该工程后乐清湾区域特定污染水的环境容量，利用计算分析结果科学分配和控制排放总量，本文基于MIKE21数学模型，建立了包括乐清湾、隘顽湾、温州湾在内的大范围水动力模型，利用实测水文观测资料对模型进行验证，在此基础上利用资料分析和数值模拟，了解乐清湾的水动力环境。其次，在水动力模型基础上建立区域水交换模式，并计算码头工程区污染物扩散过程及分析影响，从而估算区域内环境容量。

1 潮流数学模型的建立及验证

采用丹麦水力研究所MIKE 21（2007）模型建立目标海域的潮流数学模型，运用实测潮位、流速资料进行模型验证。

1.1 潮流数学模型

数学模型控制方程为：

式中：ζ为潮位，h为水深，p、q分别为x、y方向上的垂线平均单宽流量，C为谢才系数为糙率系数，H=h+ζ，n=0.022，g为重力加速度g=9.81 m/s2，f为风应力系数，V、Vx、Vy分别为风速及其在x、y方向的分量，柯氏力Ω=2 wsinφ，w为地转角速度，φ为纬度取28.1°(计算区域平均维度)，pa为大气压力，ρw为海水密度，x、y为直角坐标，t为时间，τxx、τxy、τyy分别为剪切应力分量，u、v为水平流速分量，紊动粘性系数采用Smagorinsky公式计算

为求解上述控制方程，必须要有合适的定解条件，所需的定解条件包括初始条件和边界条件。

在海域中的开边界上，海面的水位根据沿岸潮位站来给定。MIKE 21 HD采用交替方向隐格式(ADI)求解二维浅水潮波方程，方程矩阵采用双消除法（Double Sweep）算法求解，该格式具有二阶精度[2]。

模型对码头群桩的处理采用MIKE21 Pier结构，其原理为：

式中：τp为等价剪切应力，F为单个桥桩的拖曳力，n为桥桩的密度，Δx、Δy为网格间距。

总体拖曳力为群桩阻力和底摩擦力之和：

式中：CD为拖曳系数，ρ为海水密度，Be为桥桩的等效宽度，He为桥桩在水中的高度，V为流速。

1.2 模型配置

模拟排污口位于乐清湾口门，潮动力环境较强，为了能更好地模拟乐清湾大范围及局部潮流场情况，采用网格实时耦合的方法：耦合计算区分为大、中、小三层同步计算的区域，具体情况如图1，网格的Y向均为正北向向东偏转18°，网格精度分别为270 m、90 m、30 m，时间步长10 s，共计矩形网格数745 800个。

计算区域的大区北起椒江南岸-大陈岛一线、南至鳌江南岸-南麂列岛，计算大区开边界上的水位，由浪玑山、下大陈、南麂和石坪4个潮位站的水位得到，潮位站具体位置如图1。开边界处的水位同时包含11个主要分潮（K1、O1、P1、Q1、M2、S2、N2、K2、M4、MS4和M6）。温州三大河口的内边界由石仙妇、瑞安、鳌江三个岸站的水位控制（含11个分潮）；另外考虑到湾内径流对模型精度的影响，在模型中清江、楚门河等乐清湾主要河流入湾口处加强径流源，流量0.032 6 m3/s[7]，模拟期内不间断释放。乐海码头、坎门和洞头为3个潮位验证站，S1-8为8个潮流验证站（国家海洋局第二海洋研究所为甬台温复线乐清湾大桥工程实测水文站位，观测内容包括流速、流向、水深和悬沙浓度，观测时间为26～28 h[3]）。模型涉及地形主要采用本海区历史海图（13640、13710、13715、13741、13770、13771、13781）及业主单位提供的工程前沿测深图。

图1 模型配置图Fig.1 The model configuration

1.3 潮位、流速验证

图2给出了潮位验证曲线，由图可见，模拟所得的三个潮位与观测结果吻合良好，高、低潮位误差基本小于10 cm，相位基本一致，其中的差异一方面可能是由网格离散点与潮位站之间的位置差异所造成的，另一方面，可能是由于实际水深与离散网格水深取值的差异所致。

图2 潮位验证图Fig.2 Tidal level validation

分析计算与实测潮流流速可知，模拟所得的流场基本与实际情况相符，计算值和实测流速相对误差小于15%；总体看来，模拟所得的8个连续站的潮流与实测流吻合得较好，区域内流场模拟计算结果基本反映了该海域潮流和潮波的实际变化，流场模拟验证基本达到了后续预测计算的要求[4]。

2 水环境数值模拟及环境容量估算

为分析排污对海域水环境的影响，在潮流数学模型基础上建立污染物扩散数学模型，计算乐清湾水体交换变化、污染物扩散特性，以COD污染物为水质预测因子，计算分析排放典型污染物对海域水质、水环境容量的影响。

图3 大、小潮潮流流速、流向验证图Fig.3 Tidal current velocity and direction validation of spring-neap tide

2.1 海湾水交换及水体更新周期

水交换状况对环境容量有着极重要的影响，既代表了外来物质的输入能力，也反映了湾内物质的稀释能力。为模拟乐清湾的水交换状况，先假设湾内均匀分布着浓度为1的溶解态保守性物质，外海水的保守物质浓度为0，先运行潮波运动模式，待潮流稳定后开始运行溶解态保守物质的对流-扩散模式，因而在某一时刻湾内保守物质的浓度c可代表湾内水经过一段时间对外扩散后剩余的浓度，即1-c可表示湾内水体被外海置换比率。

以溶解态的保守性物质作为湾内水的示踪剂，建立乐清湾海域对流-扩散型的水交换数值模式[5]，湾内水示踪剂的控制方程为：

其中，(u，v)深度平均流速在笛卡尔坐标(x，y)方向的分量；c保守性溶解态湾内水的示踪剂浓度：k=kh+ks+kt，k扩散系数，kh垂向结构的余环流引起的水平输运的扩散系数，ks垂向剪切引起的水平输运的扩散系数，kt紊流扩散系数，这里均取模型默认稳态系数值。需要指出的是径流对于封闭海湾水交换具有一定的促进作用，乐清湾流域有山溪性河流30余条，主要有大荆溪、清江、楚门河等，多年平均径流量约103 000 m3。因此与水动力模型相同的，在模型中清江、楚门河等主要河流入湾口处加径流源强，流量0.032 6 m3/s，模拟期内不间断释放，径流边界保守物质浓度为0。

运用建立的乐清湾海域污染物对流-扩散数值模型模拟水体交换过程，根据运行10～90 d后乐清湾表层保守物质浓度分布情况可知，受半封闭狭长地形影响，从湾口向湾底，乐清湾水域中的水交换由快变慢，且差别很大，具体模拟湾内各区域的海水置换率变化情况见表1。一般认为海域水体交换置换比率达到70%即基本完成一次水交换，根据水交换模拟结果，排污口所在的玉环半岛西侧海域，完成一次水交换的时间不到30 d。模拟结果如图4所示。

图4 水交换模拟30 d保守物质分布Fig.4 Distribution of conservative material modelled in 30 days

表1 模拟湾内各区域海水置换率变化（%）Tab.1 The replacement rate of the sea area(%)

2.2 海域污染物扩散数值模拟

污水扩散模拟采用Mike21AD（2007）模型进行计算，模拟预设排污口位于乐清湾口门，距岸边250 m，不间断连续排放，根据港区控制性详规[8]拟定排污口日排放量为8 000 m3/d。污水扩散计算以COD作为水质预测评价因子研究所排废水对周边环境的影响，排放源强为30 mg/L，在污染物总量不变的前提下，排污口源强统一按流量0.1 m3/s加入。根据实测资料得到海域COD本底值分布状况，计算可得COD平均本底值为1.23 mg/L，在计算时将海域COD本底值分布状况加入模型作为污水扩散计算初始浓度。同时将上述排污口源强加入模型连续运行365 d，得到将一年的污染物全部排放以后的海域状态，再分别进行大潮和小潮期污染物扩散过程的模拟，并分析区域内COD污水最大增量分布情况，可知：

（1）大潮期：72 h内COD污水影响区域大致为北至乐清湾底，西南至瓯江北口，南至洞头峡，东至鹿西岛以东15km海域，整个乐清湾内的污水浓度在1.231～1.232 mg/L之间，其中江岩山以北的区域污水浓度均大于1.231 5 mg/L。

（2）小潮期：72 h内COD污水影响区域大致为北至乐清湾底，西南至瓯江北口，南至洞头岛西岸，东至鹿西岛以东7 km海域，整个乐清湾内的污水浓度在1.231～1.232 5 mg/L之间，其中玉环电厂厂址以北的区域污水浓度均大于1.232 mg/L。

关于排污口初始稀释度的预测，根据潮汐潮流数值模拟的结果，排污口附近流速取大潮平均流速为0.81 m/ s，取排污口流量为8 000 m3/d、单管排放、21.1 m的平均水深及海水密度1.025 t/m3，计算得出排污口初始稀释度约为1 241.5，可见由于排污口处湾口主槽内，水深流急，排污口附近海域污水扩散条件优良。

图5 大潮期COD污水最大增量分布Fig.5 Maximum increment of COD in Spring tide period

表2 污水最大浓度增量等值线包络面积（km2）Tab.2 The envelope area of the maximum concentration increment of sewage

2.3 工程海域环境容量估算

根据《中国环境保护标准汇编-环境质量与污染物排放2003》中有关污水海洋处置工程污染控制标准（GB18486-2001）的章节[6]，污水海洋处置工程污染物的混合区规定：若污水排往开敞海域或面积≥600 km2(以理论深度基准面为准)的海湾及广阔河口，允许混合区范围：Aa≤3.0 km2。若污水排往面积＜600 km2的海湾，混合区面积必须取按以下两种方法计算所得值（Aa）的小者：

式中L为扩展器长度,m。

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

对于重点海域和敏感海域，划分污水海洋处置工程污染物的混合区时还需要考虑排放点所在海域的水流交换条件、海洋水生生态等。

根据《中国海湾志》记载，乐清湾岸线以下海湾总面积为463.6 km2，其中滩涂面积（岸线至理论深度基准面)约为220.8 km2，约占海湾总面积48%；水域面积（理论深度基准面以下）约为242.8 km2，占海域总面积52%左右，平均水深10 m左右[7]。因此按照规定根据式（7）计算，得到允许混合区面积Aa为0.6 km2。

由温州市海洋环境功能区划，玉环半岛连屿至横址山前沿为四类环境功能区，乐清湾内的其他水域为二类环境功能区，考虑到污染物扩散模拟的结果，污水在牛头颈排海后将影响包括湾底养殖区域在内的乐清湾内大部分海域，因此本规划按照二类水质标准进行环境容量计算。二类水质标准为COD≤3.0 mg/L，根据评价海域的水环境质量现状监测结果，COD的海域平均值为1.23 mg/L，将目标浓度减去水质现状值，在二类海域内排污区的COD污染物允许浓度增量为1.77 mg/L。

由上述允许混合区面积，将其概化成以排污口为圆心，面积0.6 km2，半径437 m的圆型区域，讨论排污口COD排放浓度30 mg/L连续排放的情况下该区域内的环境容量。根据水体交换过程模拟的结果，允许混合区所在海域完成一次水交换的时间小于30 d；污染物扩散模式运算的起算点为水动力模块基础上加入海域平均值和排污口源强连续运行1 a后的海域状态，因此可以认为，起算时海域内污染物浓度达到动态平衡状态。统计污染物扩散模式的计算结果，排污管道建成运营后，允许混合区内的COD大潮期平均浓度增量为0.003 mg/L，最大浓度增量为0.006 5 mg/L；小潮期平均浓度增量为0.003 5 mg/L，最大浓度增量为0.007 4 mg/L，将污染物允许浓度增量乘以允许混合区的水体体积，可以得到区域内的环境容量见表3。这里由排污口流量和污水浓度得排污管排海COD污水87.6 t/a，通过计算可以得到本规划排污对允许混合区环境容量的贡献值。

表3 剩余环境容量估算及污水浓度增量对其贡献率Tab.3 Estimation of residual environmental capacity and the contribution rate of sewage concentration increment

3 结语

海洋环境容量作为海洋环境管理的重要手段，对其计算能够对海洋水质环境安全保护提供参考依据。本文以乐清湾大麦屿港区规划污染物排放项目为例，通过二维水质数学模型计算该项目对区域COD环境容量的贡献值，为后期水质环境管理工作及同类型海洋环境容量计算提供技术借鉴。

排污项目所在乐清湾为浙东南典型半封闭强潮流海湾，通过水体交换模拟可知湾口处的水体在潮汐作用下的交换周期不到30 d。湾口处排放污染物扩散条件优良，影响范围基本包括整个乐清湾海域，浓度增量除项目前沿均较小。根据水质模型的预测结果，结合项目所处海洋功能区划的水质标准，计算得到本项目允许混合区的COD污水剩余环境容量大潮期为8 148.99～8 164.98 t/a；小潮期为8 144.64～8 162.75 t/a。

[1]黄秀清.乐清湾海洋环境容量及污染物总量控制研究[M].北京:海洋出版社,2011:16-17.

[2]李 峰,谢中宇.台州港大麦屿港区规划数值模拟计算及潮间带生物现状调查报告[R].2012.

[3]国家海洋局第二海洋研究所.浙江省甬台温高速公路复线乐清湾大桥及接线工程数学模型试验报告[R]//海洋水文及工程可行性专题研究.2011.

[4]中华人民共和国交通部.JTJ/T 233-98海岸与河口潮流泥沙模拟技术规程[S].1998.

[5]董礼先,苏纪兰.象山港水交换数值研究Ⅱ.模型应用和水交换研究[J].海洋与湖沼,1999,30(5):465-470.

[6]国家环境保护总局,国家质量监督检验检疫总局.GB 18486-2001污水海洋处置工程污染控制标准[S].2001.

[7]谢钦春,冯应俊,等.中国海湾志第六分册浙江省南部海湾[M].北京:海洋出版社,1993.

[8]交通部规划研究院.台州港大麦屿港区控制性详细规划报告[R]//工程可行性专题研究.2006.

Application of Mathematical Model in the Yueqing Bay Environmental Capacity Algorithm

SHAO Zhu-feng，LIU Han-xin2，CHEN Jie-ren2
(1.The Second Institute of Oceanography,Soa,Hangzhou 310012;2.College of Water Conservancy and Hydropower Engineering,Hehai University,Nanjing 210098,China)

Based on the analysis of the measured data,the hydrodynamic and environment model is established,according to the numerical model of Yueqing bay area hydrodynamic environment and pollutant transport,calculation and analysis of the environmental capacity of typical bays have been carried out in the marine functional zoning.The calculated results is shown that:the water exchange time is 90 days in Yueqing bay, the water exchange time is 30 days near the sewage sea.The COD concentration is 1.231-1.232 mg/L at the inner bay in the tide period.The residual environmental capacity of COD is 8 144.64-8 162.75 t/ain spring tide period.The residual environmental capacity of COD is 8 148.99-8 164.98 t/ain neap tide period.The environmental capacity of Yueqing bays meet the requirements after the introduction of Damaiyu port planning.

environmental capacity;pollutant dispersion;water exchange;mathematical model;Yueqing bay

X26

A

1008－830X(2016)05－0430－06

2016-07-19

谢中宇（1983-），男，浙江苍南人，工程师，研究方向：海洋工程数学建模与环境风险.E-mail:xiezhongyu23@163.com