铁山港作为广西北部湾经济区中转运输枢纽港区，是我国西南地区便捷的出海通道，也是我国与东盟国家海上贸易的重要口岸，具有优越的区位优势。为满足铁山港湾内船舶通行和临港工业快速发展的需求，北海市在铁山港区3.5万吨级进港航道基础上拓宽浚深，总长约23.6 km，分10万吨级、5万吨级、1万吨级和5 000吨级4段进行设计。由于施工情况复杂等因素，目前仅完成Ⅰ、Ⅱ标段。Ⅲ标段东侧多生长有红树林，红树林生态系统具有丰富的生物多样性和极高的生产力，同时发挥着诸多其他生态系统所不能替代的重要功能，研究Ⅲ标段疏浚工程引起的水动力改变对周围红树林的影响具有重要意义。水动力环境变化通常通过模型来模拟。辛文杰^[1]、陈波等^[2]通过数值模拟建立铁山港的潮流模型；施华斌等^[3]利用二维浅水模型建立广西沿海的水动力模型，分析河流对近岸水动力的影响；王丽娜等^[4]利用三维浅海水动力模型建立潮流场模型；李小维等^[5]、Yang等^[6]通过历史数据和数值模型研究铁山港的水动力特征；谢洁^[7]利用数值模型研究铁山港的潮流和泥沙输运，计算得到规划港区和航道的泥沙回淤强度及回淤总量。随着铁山港航道工程的建设发展，其Ⅲ标段疏浚工程会对周边海域的水动力产生一定影响，而水动力环境是泥沙输运、海洋生态环境变化的动力学基础^[8]，确定水动力环境变化对铁山港海域环境保护、制定发展规划具有一定参考意义。因此，本研究基于MIKE数值模型，同时为使模型计算更加稳定、输出结果更加精确，在光滩和红树林区域采用不同的曼宁系数^{[9, 10]}，详细分析铁山港航道Ⅲ标段疏浚工程引起的航道及红树林区域水动力改变情况，为相关建设和规划制定提供参考。

1 材料与方法 1.1 铁山港航道Ⅲ标段疏浚工程

北海铁山港航道三期工程Ⅲ标段位于北海铁山港区石头埠作业区至雷田作业区以东海域，起点为三期航道工程设计G点，终点为雷田作业区口门处K点，各标段详细信息见图 1和表 1。北海港铁山港区航道Ⅲ标段疏浚工程施工总工期约12个月，拟疏浚航道总长约9.78 km。全航段采用8 m³(北段)和13 m³(南段)抓斗挖泥船配泥驳开挖，由抓斗船和泥驳清礁，总疏浚量约为3.869 3×10⁶ m³，炸礁量约为3.657×10⁵ m³。

工程区附近植物主要为红树植物。在低潮潮位附近分布的红树植物有白骨壤Avicennia marina，高潮潮位附近分布有桐花树Aegiceras corniculatum、卤蕨Acrostichum aureum、老鼠簕Acanthus ilicifolius和半红树植物如黄槿Hibiscus tiliaceus、苦郎树Clerodendrum inerme，岸边一带生长有厚藤Ipomoea pes-caprae、海马齿Sesuvium portulacastrum、南方碱蓬Suaeda australis等。根据周边其他近岸工程调查可知，铁山港附近海域外源性高岭土在红树林区域的堆积和黏附会妨害白骨壤的呼吸及光合作用，导致白骨壤死亡^{[11, 12]}。高岭土的悬沙浓度除了与泥沙源强相关外，海洋水动力环境也对泥沙浓度有着重要影响，因此研究疏浚工程造成的水动力环境改变，对红树林保护具有间接参考价值^[13]。

1.2 模型的建立 1.2.1 模型介绍

采用平面二维数值模型MIKE21研究目标海域的潮流场。MIKE21采用标准Galerkin有限元法进行水平空间离散，采用显式迎风差分格式离散动量方程与输运方程^[14]研究工程造成的泥沙冲淤状况。

连续方程：

$ \frac{\partial \zeta}{\partial t}+\frac{\partial h u}{\partial x}+\frac{\partial h v}{\partial y}=0, $

(1)

向动量方程：

$ \begin{array}{r} \frac{\partial u}{\partial t}+u \frac{\partial u}{\partial x}+v \frac{\partial u}{\partial y}-f v=-g \frac{\partial \zeta}{\partial x}- \\ \frac{g u \sqrt{u^{2}+v^{2}}}{c^{2} h}+\frac{\partial}{\partial x}\left(N_{x} \frac{\partial u}{\partial x}\right)+\frac{\partial}{\partial y}\left(N_{y} \frac{\partial u}{\partial y}\right), \end{array} $

(2)

向动量方程：

$ \begin{gathered} \frac{\partial v}{\partial t}+u \frac{\partial v}{\partial x}+v \frac{\partial v}{\partial y}+f u=-g \frac{\partial \zeta}{\partial x}- \\ \frac{g v \sqrt{u^{2}+v^{2}}}{c^{2} h}+\frac{\partial}{\partial x}\left(N_{x} \frac{\partial v}{\partial x}\right)+\frac{\partial}{\partial y}\left(N_{y} \frac{\partial v}{\partial y}\right), \end{gathered} $

(3)

式中：t为时间(s)；x, y为原点0置于某一水平基面的直角坐标系坐标；u、v分别为流速矢量V沿x、y方向的分量(m/s)；ζ为相对于xy坐标平面的水位(m)；h=d+ζ，为总水深(m)；d为相对于xy坐标平面的水深(m)；N_x、N_y分别为x、y向水流紊动黏性系数(m²/s)；f为科氏参量(deg)；g为重力加速度(m/s²)；c=Mh^1/6，为谢才系数；M为曼宁糙率系数。

1.2.2 初始条件

$ \begin{aligned} &\left.\zeta(x, y, t)\right|_{t=0}=\zeta_{0}(x, y), \\ &\left.u(x, y, t)\right|_{t=0}=u_{0}(x, y), \\ &\left.v(x, y, t)\right|_{t=0}=v_{0}(x, y), \end{aligned} $

(4)

式中：ζ₀、u₀、v₀分别为ζ、u、v的初始值。

1.2.3 边界条件

流速法向量为0：

$ \mathit{\boldsymbol{V}} \cdot \mathit{\boldsymbol{n}} = 0, $

(5)

式中: n为固边界法向单位矢量。

开边界水位ζ和流速V可采用已知水位ζ^*(x, y, t)或流速V^* (x, y, t)控制：

$ \left.\zeta(x, y, t)\right|_{\Gamma}=\zeta^{*}(x, y, t) \text { (水位) }, $

(6)

$ \left.\boldsymbol{V}(x, y, t)\right|_{\Gamma}=\boldsymbol{V}^{*}(x, y, t)(\text { 流速 })。$

(7)

1.2.4 模型设置

本研究计算区域如图 2所示。铁山港附近海域岛屿众多，岸线曲折，两岸有较多红树林分布。因此采用非结构三角形网格可以较好地贴合自然岸线，提高计算精度和计算效率，并便于各种工程情况的准确布置。用动边界的方法对干、湿网格进行处理，在疏浚工程区附近进行网格加密。计算空间步长在工程海域精确到10 m，网格单元74 004个，网格节点39 307个。

工程附近的海域水深数据来自于中国人民解放军海军司令部航海保证部制作的海图数据(C1416710铁山港及附近)以及航道附近测量的水深数据。对于外海水深，本项目采用美国国家地学测量中心提供的DBDB5 (Digital Bathymetric Database Version 5.2)(https://shoreline.noaa.gov/data/datasheets/wvs.html)原始数据集，通过双线性插值方法插值到网格点上，工程前后水深见图 3。

潮汐在开边界的振幅、迟角和潮流流速等要素利用MATLAB中的T_tide工具包^[15]获得，采用M₂、S₂、N₂、K₂、K₁、O₁、P₁、Q₁ 8个主要分潮。模型计算时间为2021年4月1-23日，合计23 d，计算步长根据Courant-Friedrich Levy (CFL) number笛卡尔坐标下的浅水方程式进行动态调整，确保模型计算稳定。最大时间步长为120 s，最小时间步长0.1 s。底床糙率通过曼宁系数进行控制。因为红树林的高度和分布密度等因素显著改变底床糙率，为准确刻画红树林区域的流场时空分布，底床糙率通过空间变化的曼宁系数所形成的曼宁系数场进行控制。根据倪海祥等^[16]的研究及MIKE21中曼宁数定义，光滩曼宁系数设置为林地的10倍，因此本次计算曼宁系数n取4-40。

1.3 模型验证

为验证模型的准确性，利用水动力模型模拟工程海域的潮位和潮流状况。将模拟结果与1个临时潮位站(WL)、6个潮流观测站(1#-6#)测量结果进行对比。观测站位置见图 4，水位模拟对比曲线见图 5，潮流模拟对比曲线见图 6。

对比验证结果表明：潮位、潮流模拟结果与实测潮位、潮流资料基本吻合，能够较好地反映项目周边海域水位和潮流状况。

1.4 研究内容

前述验证结果表明文章叙述的模型可以用于工程前后水动力改变情况的研究，为其提供数值模拟技术。为更直观地比较疏浚工程前后的变化情况，选取航道9个断面27个特征点的大潮期流速(图 7)进行对比分析。同时，为研究红树林区域水动力改变情况，选取铁山港东岸和西岸红树林内共16个特征点的流速进行对比分析(图 7)。

2 结果与分析 2.1 疏浚工程前后潮流场分析

根据大潮期涨落急时刻流场图(图 8, 9)分析可知，由北部湾而来的潮流沿着西南向上溯。涨潮流以偏北向流为主，落潮流以偏南向流为主，涨潮流在英罗港西侧分成两股，一股向东进入英罗港，另一股折向西北上溯进入铁山港内湾；落潮时潮流则相反。铁山港深槽区域流速较大，岸边流速较小，岸边有浅滩出露。

工程前，工程附近海域涨急时刻流速为9.2-50.4 cm/s，大致呈南高北低分布，最大流速出现在Ⅲ段南端附近；落急时刻，工程海域整体流速相比涨急时刻增大，同涨急类似，流速呈现南高北低分布，最大流速出现在工程南部海域，为68.2 cm/s(图 10)。

工程后，因疏浚造成水深增加，因此航道内流速呈现降低趋势。涨急时刻，航道内流速降低幅度为0-12 cm/s，航道南部周边区域降低值较大，超过10 cm/s，最大为11.7 cm/s。落急时刻，航道内流速降低幅度为0-5 cm/s，周边流速降低相对涨急时刻幅度减少，不足5 cm/s (图 10)。

航道大潮期流速变化如表 2所示。涨急时刻，各特征点工程后流速较工程前变化无规律，其中以4号特征点流速增加最多，最多增加6.3 cm/s；23号特征点流速降低最多，最多降低11.7 cm/s。落急时刻，4号特征点流速增加最多，最多增加3.7 cm/s；24号特征点流速降低最多，最多降低9.8 cm/s。

2.2 疏浚工程对红树林区域的影响

红树林区域大潮期流速变化如表 3所示。涨急时刻，各特征点工程后流速较工程前大部分呈现降低现象，其中以H号特征点流速降低最多，最多降低4.3 cm/s；N号特征点流速增加最多，最多增加值为0.7 cm/s。落急时刻，L号特征点流速增加最多，最多增加0.1 cm/s；M号特征点流速降低最多，最多降低0.7 cm/s。工程完成后，因疏浚造成海域水深增加，进而影响红树林区域的流速情况，致使区域流速较工程前降低，对红树林区域起到间接保护作用。

3 结论

因铁山港区Ⅲ航道疏浚工程的开展，该区域水动力环境发生改变。数值模拟结果显示, 涨急和落急时刻，各特征点工程后流速较工程前变化无规律；北部浅水区域因水深增加，造成流速增大；而南部深水区域，出现流速降低的现象，这主要和截面过流通量有关。

对红树林区域而言，涨急和落急时刻，各特征点工程后流速较工程前大部分呈现降低现象，其中以H号特征点流速降低最多，最多不超过5.0 cm/s；个别点出现流速增大的现象，但是增加量不足1 cm/s。由此可知，水动力的改变对红树林产生的影响较小。

跟前人研究结果^{[17, 18]}相比，本研究在数值模型中，考虑光滩和红树林曼宁系数的差别，使模拟结果和实测结果更加接近，但是因缺乏红树林区域实测流速资料，因此模型的效果无法定量评价。同时，本研究在航道附近设置多条断面，分别在航道左侧、中间和右侧设置特征点，研究航道疏浚对航道内部及附近流场的改变，对比结果更加全面，可以明显看出因航道疏浚造成的水深增加，进而影响流速分布的情况。