基于MIKE3的丰枯水期防城湾水交换能力研究*

2020-08-02周争桥邱文博李冠宇

广西科学 2020年3期

周争桥，夏维,邱文博,李冠宇

(1.交通运输部南海航海保障中心广州海事测绘中心，广东广州 510320; 2.国家海洋技术中心，天津 300111)

0 引言

防城湾地处北回归线以南的低纬度地区，属于亚热带海洋性季风气候，冬季温和，夏季多雨，季风明显，受灾害性天气影响较明显。防城湾分为东湾和西湾，两湾被陆地隔开，通过狭长水道连通，湾区被主城环抱，人口密集，工业发达，每年有大量的污水排入海湾。因防城河的存在，目前污水口多放置在西湾。水交换能力是评价海湾环境容量的重要指标，表征了海湾自身的物理自净能力，海湾的水质情况与水交换能力的强弱直接相关[1]。在径流和潮动力等条件下，开展防城港东、西湾海域水交换能力的研究，形成海湾交换能力的科学认识，以期制订合理的对策，降低对海湾生态环境的影响，对防城港东、西湾未来的合理开发利用有着重要的科学意义。

箱式模型[2-3]是最早的水交换能力的研究方法，随着数值计算模型的发展，对流-扩散模型[4]和Lagrange质点追踪[5]逐渐成为研究海湾水交换的主要方法。高抒等[2]基于狭长海湾多箱物理模型研究了象山港的水交换机制，该模型按垂向充分混合时间小于1 d和超过1 d两种情况分别运行，结果均表明象山港是一个水体交换缓慢的海湾。赵亮等[6]基于ECOM水动力模型，采用Lagrange质点追踪法研究了胶州湾的整体和分区域水交换能力，指出深水区水交换时间为7 d，部分海区水交换时间最长可达2个月。李希彬等[7]应用三维水动力和水交换数值模型，得出渤海湾水体半交换周期为323 d，海湾西北部海域和南部海域的水交换率较低。陈振华等[1]基于 POM 水动力模式，考虑漫滩过程，得出钦州湾的整体平均水体半交换时间约为18 d。朱金龙等[8]利用二维对流-扩散模型，研究了岸线变迁影响下的芝罘湾水交换能力变化过程。曹雪峰等[9]基于无结构三角网格的FVCOM海洋模式，在考虑径流影响环境下，建立了钦州湾海域水动力数值模型，指出内湾茅尾海的水体交换主要受制于潮汐和径流，而外湾钦州湾的水体交换受制于潮汐和波浪。陶磊等[10]建立了基于欧拉观点的渤海区域水龄模型，通过对比潮汐单独作用下与潮汐和季风共同作用下渤海湾内水龄的时空分布特征，研究了这两个因素在渤海湾水交换中的作用，指出季风是促进渤海湾水交换的重要因素。目前，对流-扩散模型和Lagrange质点追踪是研究水交换时被广泛接受和采用的方法，王聪等[11]利用Lagrange质点追踪和ECOM模式研究了大亚湾的水交换能力，分析了两种方法结果差异产生的原因。

针对防城湾水交换的研究并不多，陈宪云等[12]采用传统方法，利用形态因子面积及潮差因子计算防城港东湾2008年和2012年的纳潮量，对比2008年，2012年防城港东湾纳潮减少量占总纳潮量的5.5%，导致该海域海水交换能力减弱。蒋磊明等[13]从纳潮量、潮流场、余流和波浪4个因素分析防城港钢铁项目的围填海工程对附近海域水动力环境的影响，认为围填海工程所引起的海流、波浪、纳潮量等因素的变化均较小，纳潮量最大减少8%。上述研究为本研究提供了良好的研究基础，为了对防城湾水交换能力有全面的认识，本文基于无结构的三角网格建立了高分辨率的 MIKE3三维水动力数值模型，采用Lagrange质点追踪方法，选取水体半交换周期和水体交换律作为评价指标，开展防城湾水交换能力研究，以期为防城湾的科学管理决策提供支持。

1 材料与方法

1.1 模式介绍

采用丹麦水力学研究所的三维水动力模型MIKE3开展本研究。该模型相比MIKE21更适合深水水体的模拟研究，被广泛用于河口、海岸、海洋的研究中[14-17]。MIKE3 Flow Model采用标准Galerkin 有限元法进行水平空间离散，在时间上，采用显式迎风差分格式离散动量方程与输运方程。MIKE3模型是基于雷诺平均化的Navier-Stokes方程，满足Boussinesq假设和静水压力分布。

1.2 数值模式配置

(1)计算区域设置

根据北部湾海域的潮波运动规律，在确保计算区域边界效应对防城湾无影响的前提下，进行计算范围的确定。开边界为越南云屯县-中国北海市大风江口，网格在防城湾、防城河口和东兴湾进行了加密，计算区域及网格见图1。计算区东西宽约150 km，南北长约75 km，总面积约为50 km2，外海的计算网格步长为1—2 km，在防城河流域周围海域网格最小分辨率约为10 m，整个计算区域内水平方向上共有节点数31 879个，网格数59 364个。垂向上采用sigma坐标，垂向共分为4个σ层。

图1 计算区域网格及水深示意图Fig.1 Diagram of grid and water depth in calculation area

(2)水深和岸界

外海水深数据采用美国国家地学测量中心提供的Digital Bathymetric Database Version 5.2(DBDB5)原始数据集，近岸水深采用电子海图以及测量资料，在海图水深利用过程中，通过沿海海洋站潮位与最低理论基准面关系对海图水深进行订正，确保其准确性。

采用海图中岸界，并利用海岸线勘测资料和Google Earth最新岸线信息进行订正。

(3)计算时间步长和底床糙率

模型计算时间步长根据CFL条件进行动态调整，确保模型计算稳定进行，保证计算的CFL数小于0.8，最小时间步长0.1 s，最大时间步长120 s。

(4)模型初始条件

网格节点的水位与网格单元的流速都取0。潮汐在开边界的振幅、迟角和潮流流速等要素使用TPXO9数据。TPXO9是2018年最新发布的版本，在中国近海有较高的准确度，共包括8个基本分潮(M2、S2、N2、K2、K1、O1、P1、Q1)，2个长周期分潮(Mf和Mm)，以及M4、Ms4、Mn4、2N2、S1等分潮。本文在开边界通过潮位输入了M2、S2、N2、K2、K1、O1、P1、Q1共8个分潮。

(5)外部强迫条件设置

风速和防城河流量设置根据多年统计资料进行确定。根据防城河中游的长歧水文站 1956—2004 年49 a资料分析，全流域多年平均径流量为2.105×1010m3，多年平均径流深2 353.0 mm。丰水期(4—10月)流量可达99.85 m3/s，枯水期(11月—次年3月)流量约为丰水期的33%，流量可达33.28 m3/s。根据防城港市气象局提供的1994—2014年风速资料，丰水期多年平均风速为2.7 m/s，枯水期多年平均风速为3.7 m/s。

(6)工况设置

主要研究了两种条件下防城港东湾和西湾水体交换情况：1)丰水期+ SSW向多年平均风(2.7 m/s)；2)枯水期+NNE向多年平均风(3.7 m/s)。

1.3 评价指标

水交换是指水体通过对流和扩散等物理过程与周围水体相互混合，在海湾、近岸海域的水交换研究中，一般是通过定义各种时间尺度来描述水体交换能力的强弱或交换速度快慢，并使用半交换时间、交换率和存留时间等描述其时间尺度。水体交换率是指湾外水体通过潮汐作用进出河口或海湾，每次进入或流出湾内的外海水占整个进潮量或退潮量的比率。而水交换半周期是指在潮汐和径流的共同作用下，湾内水体交换出50% 到湾外时的周期数，在一定程度上表征了海湾的水交换能力，在钦州湾[1]、胶州湾[6]、丁字湾[18]等海湾均有所应用。

水体交换半周期在不同的潮周期下出现不同的变化与当时的潮汐现象、径流作用、风的作用、温度盐度变化等因素有关，因此要确定防城港海域水体交换的周期不能仅根据某一两次潮周期下的计算来确定，而应该是建立在多个潮周期的计算上，统计出各海域的平均水体交换率和平均水体交换半周期。

本文以海湾口门处作为分界线，选取大潮落潮中间时刻释放保守物质，保守物质浓度海湾内设置为1，湾外和径流设置为0。由于潮流的动力作用，防城湾内保守物质的平均浓度必然会随涨、落潮过程出现潮周期性的变化，也会出现大小潮的半月周期性变化，本文利用PL64TAP低通滤波器[19]滤掉潮流部分，分别计算丰枯水期防城湾东、西湾保守物质浓度变化。本文对多个潮周期的水体交换率和水体交换半周期计算取平均值。

1.4 模型验证

为了验证模型的准确性，利用水动力模型模拟了该海域的潮位和潮流状况，潮位和海流的实测站位布设见图2，经纬度坐标和观测时间如表1。

图2 潮位和潮流验证点位置图Fig.2 Location map of verification points of tide level and tidal current

表1 潮位和潮流验证点坐标Table 1 Location of verification points of tide level and tidal current

2 结果与分析

2.1 模型验证结果

(1)潮位验证

模型输出了仙人洞临时潮位站位置2019年9月29日—2019年9月30日潮位，将潮位数据和模型输出的潮位数据订正至同一基面，绘制潮位过程曲线(图3)。结果显示，模型可以较好地模拟潮位过程。

图3 仙人洞验潮站潮位实测模拟过程曲线Fig.3 Simulation process curve of tidal level measurement at Xianrendong tidal station

(2)潮流验证

采用4个潮流实测站位表层流速资料对模拟结果进行验证，观测时间为2019年9月29日—2019年9月30日。潮流验证曲线见图4。验证结果表明，对应观测点潮流模拟结果与实测资料基本吻合，能够较好地反映项目周边海域潮流状况。

图4 各站位表层潮流验证曲线Fig.4 Surface current verification curves in each stations

潮位和潮流模拟结果显示，本文建立的水动力模型可以真实反映防城湾实际的水体流动特点，可作为该海域水交换研究的模型。

2.2 流场特征分析

模拟结果(图5)表明海域潮流性质为不正规日潮，涨潮时，海水从外海进入近海各个海湾，落潮时，海水向南退入外海，落潮流流速大于涨潮流流速，涨潮流流向为偏北向，落潮流流向则为偏南向，在岸边的速度方向与岸线较为一致。在涨潮过程中，外海潮波传入防城港湾，流向都指向口门，在潮波传播过程中，受地形约束强度较大，到达渔万岛南端分为两支，一支沿牛头岭向西北方向流去，另一支沿暗埠口江向东北方向流去。流速在牛头岭深槽附近流到最大，过了牛头岭和企沙半岛后，潮流进入东、西湾内湾，由于水深变浅，水流阻力增大，流速逐渐减小。无论在外湾还是内湾，潮流都受岸线和深槽的影响，流向有顺着深槽的趋势。落潮是涨潮的逆过程，流场情况类似。在低潮时，渔万岛东北和西北面海域流速较小，较大面积的滩涂露出。西湾狭窄水道处涨落潮流速可达1 m/s以上，东湾张罗潮最大流速为0.8—0.9 m/s。

图5 大潮涨急(左)和落急(右)时刻流场Fig.5 Flow field of high (left) and low (right) tide during spring period

2.3 水体交换率

丰枯水期防城湾东湾、西湾保守物质浓度变化见图6。对保守物质随时间变化分析可知，在潮汐径流和季风的共同作用下，两个海湾区域的保守物质浓度均以较快的速度减小，西湾河口区域的保守物质浓度迅速减小，保守物质浓度减小速度明显大于防城湾东湾。由于天文潮的动力作用，两个海湾湾内平均浓度必然会随涨、落潮过程出现潮周期性的变化，也会出现大小潮的半月周期性变化。

图6 东、西湾丰水期和枯水期保守物质浓度变化曲线Fig.6 The variation curve of the concentration of the conservative substance in the East Bay and the West Bay during the wet and dry period

本文对丰枯水期不同时间段的保守粒子浓度进行了分析(图7-8)，可知：东湾的保守物质浓度的空间分布呈现由湾口向湾内增加，主水道保守物质浓度明显大于两侧浅水，说明主水道的水体交换能力较强，随着时间的增加，保守物质浓度减小速率降低。西湾保守物质浓度的空间分布呈现由河口向湾外增加，冲淡水主水道保守物质浓度明显大于两侧浅水，说明主水道的水体交换能力较强，随着时间的增加，保守物质浓度减小速率降低。

图7 防城湾丰水期保守物质浓度分布Fig.7 Distribution of conservative substance concentration in Fangcheng Bay during the wet period

2.4 水体交换时间

经计算，防城湾西湾、东湾水体半交换时间如表2所示，防城湾西湾丰水期水体半交换周期约为1.7 d，枯水期水体半交换周期约为3.6 d，东湾丰水期水体半交换周期约为15.2 d，枯水期水体半交换周期约为27.8 d。从半交换周期来看，丰水期半交换周期明显短于枯水期，径流对污染物扩散起到加速作用，作用较季风明显。