“达维”台风影响下连云港航道浮泥模拟研究

李可萌，张庆河*，高山红，张娜，，杨华

（1.天津大学水利工程仿真与国家重点实验室，天津 300072；2.中国海洋大学物理海洋教育部重点实验室，山东 青岛 266100；3.交通运输部天津水运工程科学研究所工程泥沙交通行业重点实验室，天津 300456）

0 引言

连云港海岸为黏性细颗粒泥沙组成的淤泥质海岸，在极端天气如台风过程作用下，可能形成较强烈的泥沙运动，并在航道或港池中形成浮泥，甚至影响船舶的正常航行[1]。因此，研究台风影响下连云港海域的泥沙运动和航道中浮泥变化过程对于连云港深水航道建设具有重要意义。本文利用大气模型WRF模拟风场，利用SWAN风浪模型和FVCOM水动力泥沙模型建立波流泥沙耦合模型模拟海域水动力和泥沙运动，并将二维浮泥运动和固结模型嵌入FVCOM中，形成三维黏性泥沙和浮泥运动模型，对2012年“达维”台风影响下连云港海域的泥沙运动和航道中浮泥的形成与变化过程进行模拟，并和现场实测数据进行比较，以讨论强风浪过程中浮泥形成对连云港航道的影响。

1 浮泥运动模型介绍

连云港台风期间浮泥运动模型由以下几部分构成：包含数据同化模块WRFDA的大气模型WRF，波浪模型SWAN，水动力模型FVCOM，以及悬沙运动和浮泥运动及固结模型。

WRF模型是得到广泛应用的中尺度大气模型[2]，可以模拟台风过程风场。WRF模型中WRFDA模块是可以单独安装编译单独使用的模块，目前支持三维变分数据同化、四维变分数据同化、集合卡尔曼变换数据同化以及混合资料同化，通过卫星遥感或实测数据同化可以提高风场模拟精度。

波浪模型SWAN是由荷兰Delft理工大学开发的第三代海浪模式。SWAN模型可以模拟风浪生成和近岸区传播的多种物理过程，可描述台风浪在近岸区的变形[3]。

FVCOM是由美国马萨诸塞大学海洋科技研究院和伍兹霍尔海洋研究所联合开发，主要用于三维水动力泥沙数值模拟。模型采用非结构三角形网格，可适用于复杂岸线情况[4]。浮泥模型选用杨小宸所建浮泥运动数学模型[5-6]，并嵌于FVCOM模型中。

上述模型实现了近岸波流耦合，其数据传递利用MCT实现[7]，模型框图如图1所示。

图1 浮泥模型数据交换示意图Fig.1 Data exchangeof fluid mud model

2 2012年“达维”台风过程风浪流与泥沙运动模拟

2012年第10号超级台风“达维”于北京时间2012年7月28日20时，在日本东京东南方约1 330 km的西北太平洋洋面上生成。北京时间7月31日早晨，在西北太平洋洋面上加强为强热带风暴。北京时间8月1日08时，在日本九州岛东南部海面加强为台风。下面采用本文前述模型对该台风过程风场、波浪场、水动力、悬沙和浮泥过程进行模拟。

2.1 模型设置

2.1.1 WRF模型设置

WRF模型采用大小模型嵌套计算，其中大模型网格空间精度为9 km，网格数为300×200，大模型的纬度范围为27.04°N—42.25°N，经度范围为116.03°E—143.97°E。小模型网格空间精度为3 km，网格数为390×390，小模型的纬度范围为29.02°N—40.81°N，经度范围为116.00°E—129.42°E。小模型的计算地形范围见图2。WRF模型垂向采用σ坐标，分为34层。WRF模型初始场采用美国国家环境预报中心NCEP网站提供的6 h一次、精度1°×1°的FNL再分析数据（https://rda.ucar.edu/datasets/ds083.2/），初始温度采用全球海洋观测系统（NEAR_GOOS）日平均海表温度数据（MGDSST）（http：//near-goos1.jodc.go.jp/）。

图2 波浪小模型示意图Fig.2 Subdomain of wave model

为了提高风场模拟精度，本文采用NOAA-18单颗卫星的AMSUA亮温资料（https://rda.ucar.edu/datasets/ds735.0）进行数据同化，利用CRTM快速辐射传输模式处理AMSU-A亮温资料。AMSU-A全功率微波扫描辐射计一共包含15个探测通道，本文选择5～9通道，同化资料时间窗口设定为同化时刻的±1 h。

2.1.2 SWAN模型设置

SWAN模型采用水平向大小模型嵌套计算，小模型纬度范围为34.4°N—35.4°N，经度范围为119.3°E—120.1°E，网格空间步长最大为1 500 m，最小为20 m，地形范围见图2。风场模型嵌套区域的计算结果为波浪模拟提供风场，波浪大模型计算结果为小模型提供波浪谱边界条件，大、小模型计算时间步长均为5 min。

2.1.3 FVCOM模型和泥沙模型设置

水动力模型采用与波浪模型相同的大小区域网格进行嵌套计算。由风场模型嵌套区域的计算结果为潮流模拟提供风场，大模型边界水位由Mike软件系统[8]提供，小模型边界水位由大模型计算结果提供。大模型外模计算时间步长为1.0 s，内模计算时间步长为10 s。小模型外模计算时间步长为0.2 s，内模时间步长为2 s，垂向平均分为10层，底部粗糙高度取为0.000 1 m。

连云港泥沙中值粒径d50取为0.004 mm，密度取为2 650 kg/m3。冲刷模型[9]见式（1）。

式中：E0为冲刷系数，取为0.001 5 kg/（m2·s）-1；τb为床面剪切应力；τc为临界冲刷应力，取为[10]1.20 N/m2；a1、a2、a3、a4分别取为-0.144、0.904、-0.823、0.204。淤积模型[11]见式（2）。

式中：τd为临界淤积应力，取为0.07 N/m2；C为质量含沙量。黏性泥沙沉速公式[12]见式（3）。

式中：a、b、m、n分别取为1.05、0.58、2.37和4[5]。浮泥的流变参数根据实验来确定[1]，流变方程见式（4）。

式中：τB为宾汉体临界切应力；η为刚体系数。τB和η的经验关系见式（5）。

式中：Sv为体积含沙浓度；α1取为0.188 4；α2取为0.000 54；β1取为2.186 5；β2取为2.069 6。

SWAN模型和FVCOM模型通过MCT进行数据交换，数据交换时间间隔为30 min。

2.2 模型模拟结果与实测结果的比较

为了验证模型的合理性，收集研究区域测站气象、水文资料与模拟结果进行比较，“达维”台风过程气象和水文测站位置标注在图2中。

2.2.1 风场模拟结果

图3显示了同化模拟得到的“达维”台风路径与中央气象台确定的台风路径比较情况（以下均指北京时间）。

图3 台风路径结果对比图Fig.3 Comparison of tracks of Typhoon Damrey

图4 为气象测站A处模拟风速、风向与实测结果的对比结果。由图可知，计算台风路径与实测路径吻合较好，除了最大风速有所差异外，模拟与实测风速和风向总体是比较一致的。总体上看，风场模拟结果较好地反应了实际“达维”台风场，可以进一步作为波浪和潮流计算的输入条件。

2.2.2 波浪与潮流场模拟

图4 测站A风速、风向结果对比图Fig.4 Comparison of wind speedsand directions of station A

利用计算得到“达维”台风风场驱动波浪模型和潮流模型并进行波流耦合计算。图5显示了测站A的有效波高模拟结果与实测结果的对比图。图6是测站A潮位计算结果与实测过程比较情况。图7显示了测站B的模拟表层流速和流向与实测结果的对比情况。可知，波浪和流场模拟结果与实测结果吻合较好。

图5 测站A有效波高对比图Fig.5 Comparison of significant waveheightsof station A

图6 测站A潮位对比图Fig.6 Comparison of water elevations of station A

图7 测站B流速、流向表层对比图Fig.7 Comparison of current velocities and directions of station B

2.2.3 泥沙浮泥模拟

连云港徐圩航道浮泥测站位置见图2，沿航道包含E1、E2、E3、E4、E5和E6六个测点。“达维”台风风后测量采用音叉密度仪和分层取样器确定浮泥厚度[13]。

利用考虑“达维”台风影响的波浪、潮流耦合模型来驱动泥沙、浮泥模型，得到测站B的含沙量底层结果如图8所示。含沙量最大时刻的底层含沙量分布图见图9。比较悬沙与波浪、水流变化过程可知，含沙量变化趋势与波浪的变化趋势一致。最大含沙量出现在波浪最大值附近，说明波浪起着重要的掀沙作用。由于含沙量测量是在风后进行的，此时测站波浪已经明显减小，因此实测底层含沙量也比较小（见图8中8月5日8时—8月6日6时计算与实测含沙量）。

图8 测站B含沙量底层对比图Fig.8 Comparison of sediment concentration of bottom layer of station B

图9 含沙量底层分布图Fig.9 Distribution of sediment concentration of bottom layer

图10 徐圩航道测站浮泥厚度变化历时图Fig.10 Thickness change of fluid mud of Xuweichannel stations

图10 显示了测点浮泥厚度随时间变化的过程。浮泥厚度变化趋势与含沙量变化趋势趋于一致，当8月2日18时含沙量受波浪作用急剧增大时，测站开始产生浮泥，之后至8月3日12时含沙量一直很大，此时淤积起主要作用。8月3日18时至4日0时，浮泥厚度趋于稳定。8月4日6时之后，台风影响基本消散，含沙量变得很小，此时浮泥厚度开始下降，固结起主要作用。图中同时反映了8月5日8时测站浮泥厚度对比情况，可见计算值与实测值吻合较好。可以认为浮泥模型较好地反映了“达维”台风期间航道浮泥厚度先增大后减小的变化趋势。洪柔嘉[14]定义重度为1 060～1 200 kg/m3的低浓度淤泥为浮泥。本文模型中模拟计算的浮泥密度随时间推移增加，由1 080 kg/m3逐渐增加至1 180 kg/m3，始终属于浮泥重度范围内，体现了浮泥的固结作用。

3 结语

利用包含数据同化模块WRFDA的大气模型WRF，波浪模型SWAN，水动力模型FVCOM，以及悬沙运动和浮泥运动及固结模型建立起的浮泥运动集成模型，模拟了“达维”台风作用下连云港的黏性泥沙及浮泥运动情况。模拟结果与实测值吻合较好，表明该模型能够合理地反映“达维”台风作用下连云港黏性泥沙及浮泥运动过程。波浪具有重要的掀沙作用，含沙量随波高的增大而增大，当波高减小时，含沙量也随之降低。浮泥的厚度变化趋势与含沙量的变化趋势趋于一致，浮泥厚度先增加后减小，前期淤积起主导作用，后期固结起主导作用。浮泥的密度随时间推移增加，体现了浮泥的固结作用。