二维潮流泥沙数学模型 在挖入式港池中的应用研究

谈宝林 王俊虎

摘 要：本文根据长江口南港水域历年水文及地形资料建立的长江口二维潮流泥沙数学模型，模拟分析计算中船长兴造船基地二期工程挖入式港池方案实施后，周边水域的水动力变化以及港池泥沙淤积强度，港池内水体的交换能力，为工程项目设计及后续运营提供技术支持，为类似项目提供参考。

关键词：长兴二期;挖入式港池;泥沙回淤;水体交换

中图分类号：U653.3          文献标识码：A            文章编号：1006—7973（2021）12-0120-04

1项目背景

中船长兴造船基地二期工程位于上海市崇明区长兴镇，西起跃进港以东、东至横沙小港以西，南临长江。岸线自上游至下游依次为挖入式港池、2#舾装码头、2#船坞、3#舾装码头、4#舾装码头。

本次拟建挖入式港池：在用地上游最西侧利用海庆河入海口现状，港池大致成EN-WS走向，港池宽度240m。挖入式港池内布置两座舾装码头，其中港池西侧为5#舾装码头，长度720m、宽度为20m;港池东侧为6#舾装码头，长度为760m、宽度为20m。港池北侧布置1座室内船坞，其余为工作船码头。

2建立长江口潮流泥沙数学模型

2.1计算方法

2.1.1二维浅水控制方程

在笛卡尔直角坐标系下，根据静压和势流假定，沿垂向平均的二维潮流基本方程可表述为如下形式：

（1）

（2）

（3）

（4）

（5）

2.1.2基本方程的离散及求解

方程（1～4）垂向积分变成二维形式，写成向量形式：

（2-6）

2.1.3边界处理

边界分为开边界和闭边界。由于本文采用的是有限体积法，水位、流速布置在网格中心点，网格边界上没有布置变量，因此不能够通过网格边界处理边界条件，需用到特殊的边界处理方法。

2.1.3.1开边界

对于边界处的网格，可求，关键是求解，开边界又分为急流开边界和缓流开边界，因本文模型为缓流开边界，故下文只讨论缓流开边界，根据相容关系处理。

2.1.3.2闭边界

采用镜像法处理。在闭边界外侧虚拟一个单元，边界上的两侧的法向流速相反，切向流速相同，即表示单元法向和切向流速。

2.1.4关键技术问题

2.1.4.1动边界的处理

本次研究采用冻结法，根据计算单元水深判断是否露滩，当水深小于某一控制水深时，单元潮位“冻结”不变，要进行下一时刻计算前，被冻结的单元水深由周边有效水深进行修正，如果水深大于控制水深则重新参与计算，为避免水量不平衡，动边界控制水深采用0.01m。

2.1.4.2糙率处理

糙率是潮流计算的主要计算参数之一，反映了潮流运动过程中的阻力特性，糙率选取正确与否对计算结果有直接影响。糙率是一个综合参数，与床面泥沙特性、水深及地形形态都有一定关系，本项研究中根据经验选用了附加糙率公式，考虑水深变化后的糙率响应。

2.1.4.3桩基模拟

由于本次挖入式港池两侧建有舾装码头，结构为高桩承台，因此本报告桩基的模拟采用概化方式。一般数学模型中对于桩基模拟常用的方法有加密网格法、附加糙率法、附加阻力法和码头区不过水法等，本次模拟中为更好地模拟桩基引起的阻力，采用码头区地形概化和加糙相结合的方法，首先将桩群看作是阻水建筑物，将码头平台阻水部分当作断面突然缩水的建筑物考虑， 进而，根据桩径计算桩基的阻水体积，将其换算至桩所在的网格单元，计算网格地形凸起高度，修改网格地形。

2.2数学模型建立

本次数学模型采用嵌套模型计算方式。大模型为中国近海模型，嵌套模型为长江口近海模型，长江口近海模型的边界条件由中国近海模型提供。

2.2.1 中国近海模型

中国近海模型，南边界至越南的Qui Nhon—马来西亚的Subic Bay连线，台湾岛，东南侧边界沿巴布延群岛—八示戈岛—石桓岛—冲绳群岛，并与日本的南九州地区相连，东侧边界由日本的福冈与韩国的釜山相连。

模型采用三角形网格作为计算单元，网格总数为271156个，最小网格边长约为244m，最大网格边长约为28598m。

2.2.2长江口近海模型

局部小模型包括整个长江口和杭州湾，数学模型上边界至南京，外海边界至长江口口外-40m等深线左右。模型长江口南京至徐六泾采用2016年实测地形资料，北支采用2010年实测地形资料，南支、南港、北港、南槽、北槽采用2019年最新实测地形资料，外海地形采用最新海图拼接。

模型采用三角形网格，共计划分单元122224个，节点总数62877个，最小网格边长约15m，最大网格边长约为14935m。

2.3模型率定

采用2017年7月10日～7月19日的长江口洪季实测同步水文资料进行验证，分别验证潮位、潮流和含沙量过程，均得到了较好的模拟，与实测值吻合较好。

2.4水流泥沙验证

采用2019年1月11日～1月30日长江口枯季实测同步水文资料进行验证。潮位各站模型计算潮位过程与实测值吻合较好。涨落潮流速偏差在10%以内。含沙量偏差在20%以内，部分测点含沙量偏差较大，主要因为现场影响含沙量因素多且复杂，数学模型无法全部考虑。

2.5航道回淤验证

以2016年北槽航道淤积资料为地形验证资料，回淤验证外海潮型采用大-中-小潮组合潮型，上游流量采用多年洪季平均流量和多年枯季平均流量組合的方式。具体计算方式为：“多年洪季平均流量+大-中-小潮潮型”计算30天，“多年枯季平均流量++大-中-小潮潮型”计算30天，然后采用时间加速计算航道年回淤量。

2016年北槽航道淤积量为5400万方，数学模型计算量为6321万方，O～V单元数学模型计算值较大，总回淤量偏差约为17%。总体来说，数学模型计算的航道淤积分布能够反映北槽的航道淤积情况。

建立的长江口大范围二维潮流泥沙数学模型，以2017年、2019年实测水文、航道回淤资料作为水动力、泥沙率定验证资料，结果满足规程要求，说明该数学模型，能够反映长江口的水沙运动特征。

3水动力变化分析

3.1计算条件

外海潮汐条件采用2017年洪季验证边界条件，采用2018年长兴站的实测潮位观测资料统计，累积频率为10%的典型大潮潮差约为3.4m，累积频率为90%的典型小潮潮差约为1.5m，数学模型计算的长兴站大潮潮差约为3.3m，略小于典型大潮，长兴站大潮潮差约为1.45m，基本接近典型小潮。

水动力模型计算外海边界采用大潮潮型，径流采用长江多年洪季平均径流;泥沙模型计算外海边界采用大中小组合潮型，径流分别采用长江多年洪季、枯季平均径流;水体交换模型计算外海边界分别采用大潮和小潮潮型，径流采用长江多年洪季平均径流。

3.2潮位变化

沿着岸线，在港池口门处取一点，港池西北侧取一点，港池的东南侧位依次取三点，作为潮位变化分析采样点。经模拟分析，从潮位过程变化看，工程前后潮位过程基本没有变化，从高低潮位变化看，工程方案实施后，长兴水道的高低潮位基本不变，因此方案的实施对周边潮位基本没有影响。

3.3流态变化

经模拟分析本底条件下，涨落潮流向受南港岸线的约束，潮流为往复流特征，涨潮方向为东南向西北，落潮方向相反，为西北向东南，由于工程区附近的长兴水道等水域水深条件较好，因此，涨落潮流较平顺，但在拟建挖入式港池附近形成局部的回流。

经模拟分析方案实施后，南港及长兴水道大范围的潮流流态变化相对较小，与本底基本一致。挖入式港池由于水域的开挖作用，涨落潮流态呈现出明显的回流特征，具体表现在：涨潮时，在东南向西北的涨潮流作用下，港池口门处出现顺时针的回流区，同时港池内部在口门回流的作用下，形成逆时针方向的二次回流;落潮时，在西北向东南向的落潮流作用下，港池内形成大范围的逆时针回流。与涨潮回流不同的是，涨潮回流为二次回流，而落潮回流未出现二次回流，主要原因为落潮动力较强，造成港池口门区的水流流速相对较大，因而回流的范围增大。

3.4流场变化

经模拟分析挖入式港式方案实施后，对南港主槽及长兴水道水域潮流基本不影响。港池内部的流速变化是由于港池开挖后，其内部自身的水动力引起，而非与本底的流速对比变化。

3.5港池及近岸水流条件分析

经模拟分析，5#、6#舾装码头前沿流速基本在0.15～0.20m/s之间，流速相对较弱，口门处流速相对较大，涨落急流速可达到0.60～0.70m/s左右。对于挖入式港池内部，自口门向港池内部，涨落急流速均呈逐渐减小趋势，对于口门附近，涨落急流速基本在0.20m/s左右，而港池内部涨落急流速在0.05～0.10m/s左右，总体较小;港池口门处，涨落急流速可达到0.60～0.70m/s左右;港池上游，涨急流速约为0.90m/s左右，落急流速在1.30m/s左右，港池下游涨急流速基本在1.00m/s左右，落急流速在1.20～1.50m/s左右，落急流速总体大于涨急。码头前沿最大横流均在0.10m/s以内，幅度较小。

4挖入式港池回淤计算

4.1含沙量分布

经模拟统计分析，港池自口门向底部，含沙量呈逐渐减小趋势，但变化幅度相对较小，口门区平均含沙量约为0.23kg/m³，港池末端含沙量约为0.18kg/m³。

4.2泥沙回淤分析

经模拟统计分析，港池口门处回淤强度较大，往港池内部，回淤强度逐渐减小，口门处最大回淤强度约为2.01m/a，港池末端回淤强度约为0.76m/a，挖入式港池内的平均回淤强度约为1.29m/a，回淤量约为22.9万方/a。5#舾装码头前沿淤积强度约为1.07m/a，6#舾装码头前沿淤积强度约为1.27m/a，港池南侧淤积强度总体大于北侧。

4.3大风天骤淤计算

采用麦德逊（Madsen）和格蘭特（Grant）的波浪推移质输沙率公式计算底沙骤淤，当海向大风（强台风）出现时，含沙量可能呈数倍增大，采用刘家驹公式作近似推算。经计算，港池内一次大风过程（50年一遇S向波浪作用12h），推移质骤淤强度约为0.29m/次，悬移质骤淤约为0.05m/次，因此港池内一次大风过程骤淤强度约为0.34m/次。

5挖入式港池水体交换分析

本文采用数学模型示踪剂方法，初始计算时，在模拟水域内充满浓度为1的示踪剂，连续计算50h，分析港池内示踪剂浓度的变化。水体交换能力评价指标主要有水体交换率、半交换周期、更替周期等，本文采用半交换周期分析水体交换能力。

经模拟计算分析，大潮条件下港池内的水体交换速度明显大于小潮，大潮港池内的水体交换半周期约为4.5h，2/3水体交换周期约为8h，小潮港池内的水体交换半周期约为7.5h，2/3水体交换周期约为12.5h。

6结语

在分析南港历年水文及地形资料的基础上，建立了长江口二维潮流泥沙水质数学模型，以2017年洪季、2019年枯季同步水文测验资料为率定和验证资料，模型相似性良好，能够复演南港水域的水沙运动。在此基础上，计算了挖入式港池方案实施后，周边水域的水动力变化以及港池泥沙淤积强度，分析了港池内水体的交换能力，主要结论如下：

（1）长江口南港水域潮汐类型为正规半日潮，M2分潮占绝对主导地位，涨潮流向一般为295º～313º，落潮流向为115º～122º。据1999～2019年南港实测水文测验统计，南港涨、落潮平均流速未见明显趋势性变化。采用2018年长兴站水位资料计算了台风“安比”引起的增水值，经计算，安比台风在长兴站引起的最大增水约为0.5m。

（2）南港河段涨潮平均含沙量大于落潮平均含沙量，三峡水库蓄水后，涨潮和落潮平均含沙量降低较明显。南港悬沙中值粒径大都在0.008～0.014mm之间，南港底质以细砂、极细沙为主，中值粒径大都大于0.062mm。

（3）南港河床地形变化与上游沙体运动密切相关，河槽断面形态呈“U”“W”型交替变化;1998～2019年，工程区附近长兴水道水深条件总体较好，深槽水深基本大于10m。

（4）工程方案实施后，工程附近水域高低潮位变化均在0.002m以内，几乎不变，因此方案的实施对周边潮位基本没有影响;涨潮时，港池口门处出现顺时针的回流，港池内部形成逆时针方向的二次回流，落潮时，港池内形成大范围的逆时针回流，南港主槽及长兴水道处等水域潮流基本不变;挖入式港池码头前沿最大横流均在0.10m/s以内，幅度较小。

（5）港池自口门向底部，含沙量呈逐渐减小趋势，口门区平均含沙量约为0.23kg/m³，港池末端含沙量约为0.18kg/m³。港池口门处回淤强度较大，往港池内部，回淤强度逐渐减小，口门处最大回淤强度约为2.01m/a，港池末端回淤强度约为0.76m/a，挖入式港池内的平均回淤强度约为1.29m/a，回淤量约为22.9万方/a。

（6）港池内一次大风过程（50年一遇S向波浪作用12h），推移质骤淤强度约为0.29m/次，悬移质骤淤约为0.05m/次，港池内一次大风过程骤淤强度约为0.34m/次。

（7）大潮期，港池内的水体交换半周期约为4.5h，2/3水体交换周期约为8h;小潮期，港池内的水体交换半周期约为7.5h，2/3水体交换周期约为12.5h。

参考文献：

[1] JTS/T231-2010，海岸与河口潮流泥沙模拟技术规程[S].

[2]罗小峰.中船长兴造船基地一期工程水工建筑物物理模型试验研究[Z].南京水利科学研究院，2005.

[3]陆永军，左利钦，王红川，李浩麟.波浪与潮流共同作用下二维泥沙数学模型[J].泥沙研究，2005（06）：3-14.

[4] 韩玉芳.中船长兴造船基地二期工程物理模型试验研究[Z].南京水利科学研究院，2007.