CN108256137A - 一种强潮河口湾人工岛作业区港池航道回淤模拟方法 - Google Patents

Abstract

本发明公开了一种模拟强潮河口湾人工岛作业区港池航道回淤的方法，采用水流泥沙互馈计算模式，采用随地貌变化和潮位变化的变糙率建立水流数学模型；计算不同潮位下代表波辐射应力；优化人工岛平面形态；进行波浪、潮流、径流输沙耦合模拟；最后统计回淤量。本发明的方法波流泥沙数学模型构建过程针对半封闭强潮河口湾离岸人工岛海岸岸线曲折、潮汐强进出潮量巨大、潮滩宽广地貌复杂、旋转流作用下的潮流强、受波浪影响大、泥沙运动受波流共同作用的独特动力地貌特点，考虑随地貌变化和潮位变化的变糙率以及不同季节径流输沙的影响，能够适应与满足不同地貌形态的强潮河口湾人工岛作业区港池航道回淤模拟的要求。

Description

一种强潮河口湾人工岛作业区港池航道回淤模拟方法

技术领域

本发明涉及人工岛回淤模拟方法，特别是涉及一种强潮河口湾人工岛作业区港池航道回淤模拟方法。

背景技术

河口湾是河口海岸地区发育的众多形态各异海湾中重要的一类，其形成与地质地貌条件、水动力条件等多因素密切相关，这些控制因素相互制约同时又相互转化，在人类开发干预下往往严重影响海湾的演变趋势与滩槽稳定性。利用海湾潮滩围垦是历史上强潮河口湾开发的重要方式，过量的围垦往往不仅使得海湾面积急剧减少，还因减少纳潮量致使海湾淤积明显。为科学合理地利用海岸资源，国内外强潮海湾常在充分考虑水动力环境特征的基础上因势利导并实施相应的整治工程。为保护天然岸线，将海岸利用过程中的生态损失降到最小，离岸人工岛围填目前已成为世界主流。由于动力地貌特征复杂，强潮河口湾人工岛围填势必会使湾内水动力环境发生改变，进而可能影响其滩槽稳定。因此，有必要系统研究强潮河口湾水沙运动过程、滩槽形成与稳定机制，通过数学模型模拟和动力地貌分析评估人工岛作业区港池航道的回淤。

半封闭强潮河口湾波浪作用一般较强，波浪潮流共同作用是湾内泥沙运动的主要动力。尤其是面对海峡的强潮河口湾还因旋转流作用而潮流较强，强潮流作用下深槽发育，涨落潮流路分离，发育有拦门沙。受入湾河流输沙以及强潮汐作用影响，这类海湾通常发育面积宽广的潮滩，因而海湾地貌复杂多变，其动力地貌特征要远复杂于一般的大陆海岸。对于拟人工岛围填的半封闭强潮河口湾海区，首先要在动力地貌特征和海床演变等专业分析的基础上明确拟围填人工岛位置，通过潮流模型优化围填方案，最后从对周边滩槽的影响和进出潮量影响等角度判断人工岛围填适宜性。判断人工岛围填适宜性一条重要原则是对周边潮量和滩槽稳定影响最小。现有关于人工岛围填水动力地貌研究一般多为针对粉沙淤泥质平原海岸、沙质海岸，缺少模拟半封闭强潮河口湾人工岛作业区港池航道回淤的方法。

发明内容

本发明的目的在于解决现有技术中缺乏模拟半封闭强潮河口湾人工岛作业区回淤的方法的问题，提供一种强潮河口湾人工岛作业区港池航道回淤模拟方法。本发明针对了该类海岸岸线曲折、潮汐强进出潮量巨大、潮滩宽广地貌复杂、旋转流作用下的潮流强、受波浪影响大、泥沙运动受波流共同作用的独特动力地貌特点，考虑随地貌变化和潮位变化的变糙率以及不同季节径流输沙的影响，能够适应与满足不同地貌形态的强潮河口湾人工岛作业区港池航道回淤模拟的要求。

为解决上述技术问题，本发明采用以下技术方案。

步骤一、根据作业区水下地形实测资料及相关水文资料，选取坐标系和基本控制方程，建立二维平面潮流数学模型并进行数值求解；

所述基本控制方程包括水流连续方程、动量方程；

步骤二、建立波浪数学模型，计算不同潮位下代表波的辐射应力；根据实测资料，换算得出不同浪向所占比重，根据波浪数学模型计算所考虑的不同浪向在高、中、低潮位下的波要素，并通过线性插值得出其他潮位下的波要素，根据波浪要素计算波浪应力；

步骤三、根据数学模型模拟结果，优化人工岛平面形态；

步骤四、根据入湾径流实测输沙量，分洪水、中水和枯水统计各自平均含沙量；径流输沙与潮流耦合计算采用洪、中、枯季加权所得；径流、潮流耦合后的结果再叠加波浪作用，分径流、潮流作用以及径流潮流叠加波浪作用分别模拟，最终的模拟结果进行加权求和；

步骤五、利用GIS工具计算拟关心区面积，统计回淤量，并进行机理分析。

本发明的方法，所述步骤一中，所述二维平面潮流数学模型采用随地貌变化和潮位变化的变糙率建立，所述变糙率采用基本糙率和水深调节糙率相加的方法确定，基本糙率采用0.015～0.025；水下地形变化较大处根据水深调节糙率，即计算糙率a₁为基础糙率，a₂为调节糙率，h为水深。

进一步的，二维平面潮流数学模型的基本控制方程如下：

水流连续方程：

ξ方向动量方程：

η方向动量方程：

式中，t表示时间；g为重力加速度；ρ表示海水密度；ξ、η分别表示正交曲线坐标系中二个正交曲线坐标；u、v分别表示沿ξ、η方向的流速；h表示水深；H表示水位；C_ξ、C_η表示正交曲线坐标系中的拉梅系数；

σ_ξξ、σ_ξη、σ_ηξ、σ_ηη表示紊动应力；

ν_t表示紊动粘性系数，ν_t＝C_μk²/ε，采用k-ε模型计算ν_t；v_t＝αu*h，α＝0.5～1.0，u*表示摩阻流速；

波浪底质点速度u_w,v_w采用式(6)计算：

L_w表示波长，T_w表示波浪周期，H_w表示波高，u_w,v_w分别表示波浪底质点速度，B表示波浪与潮流相互影响系数，当波、流同向时B＝0.917，当两者互相垂直时B＝-0.198，当方向不定时，B＝0.359；

S_ξξ、S_ξη、S_ηξ和S_ηη表示波浪辐射应力张量的四个分量，表达式为：

其中，H_w、θ分别表示波高和波向，C、C_g分别表示波速和波群速。

本发明的方法，所述步骤三中，含沙量数学模型的基本控制方程如下：

悬沙不平衡输移方程：

非均匀悬移质按其粒径大小分成n₀组，S_L表示第L组粒径含沙量，P_SL表示此粒径悬沙含沙量所占的比值，则：

针对非均匀悬移质中第L组粒径的含沙量，二维悬移质不平衡输沙基本方程为：

式中，表示第L组泥沙的挟沙能力，ω_L为第L组泥沙的沉速；K₀为挟沙能力系数，α_L为第L组泥沙的含沙量恢复饱和系数；

海域泥沙沉速受含氯度影响，含沙量验证计算表明，用式(12)计算沉速，当ω_L小于0.015cm/s时，ω_L取絮凝后沉速0.015cm/s，当ω_L超过0.015cm/s时，则采用式(12)计算值；

床沙级配方程：

其中，E_m表示混合层厚度；P_mL0表示原始床沙级配，P_mL表示床沙级配；当混合层在冲刷过程中波及到原始底床时，ε₁＝0，否则ε₁＝1；

床底变形方程：

床底总冲淤厚度

由此，完成判断半封闭强潮河口湾离岸人工岛围填适宜性全过程。

本发明针对了强潮河口湾人工岛作业区港池航道海岸岸线曲折、潮汐强进出潮量巨大、潮滩宽广地貌复杂、旋转流作用下的潮流强、受波浪影响大、泥沙运动受波流共同作用的独特动力地貌特点，考虑随地貌变化和潮位变化的变糙率以及不同季节径流输沙的影响，能够适应与满足不同地貌形态的强潮河口湾人工岛作业区港池航道回淤模拟的要求。

现有技术中，多为针对人工岛海床的稳定性、地质环境稳定性、对表层沉积物粒度及黏土矿物组成特征的影响、局部冲刷或水动力环境影响或在波浪载荷作用下的稳定性以及对波浪传播、波生流、泥沙输运与岸滩演变的影响。所采用的方法有整体潮流泥沙物理模型、波流泥沙数学模型和遥感影像资料结合地貌演变分析。如现有技术通过纳潮差，即纳潮量与平均高、低潮潮位的水域平均面积之比评价人工岛的水动力环境影响，但这种方法主要适用于人工岛相对海湾面积比较大的区域。研究区域类型有岬角海湾的沙质海岸、粉沙淤泥质海岸、辐射沙脊群、河口口外海滨、沙质海岸、海峡以及基岩海湾等。

现有技术主要通过经验公式计算泥沙冲淤变化；与本发明同采用水流泥沙互馈计算模式的区别在于，现有技术研究对象为淤泥质平原海岸等；与本发明同采用贴体坐标系、波浪辐射应力、泥沙粒径分组的区别在于，现有技术研究对象为半封闭海湾、河口等。

与现有技术相比，本发明克服了以往技术中缺乏模拟半封闭强潮河口湾人工岛作业区回淤的方法的问题。波流泥沙数学模型构建过程中坐标系选取、对泥沙和波浪的处理完全针对该类海岸岸线曲折、潮汐强进出潮量巨大、潮滩宽广地貌复杂、旋转流作用下的潮流强、受波浪影响大、泥沙运动受波流共同作用的独特动力地貌特点，考虑随地貌变化和潮位变化的变糙率以及不同季节径流输沙的影响，能够适应与满足不同地貌形态的强潮河口湾人工岛作业区港池航道回淤模拟的要求。

附图说明

图1是本发明所述方法的流程图。

图2是本发明实施例泉州湾所在位置、地貌特征、人工岛及潮量断面(高程单位：m)。

图3是本发明实施例泉州湾网格间距分布图。

图4是本发明的一个实施例泉州湾二维波流泥沙数学模型范围图。

图5是人工岛方案实施年后港池航道内泥沙回淤分布(m)

具体实施方式

下面结合附图和具体实施方式对本发明作进一步的详细说明。

本发明设计方法的原理为：对于河口湾拦门沙浅滩区围垦(如人工岛)，进出潮量变化由两部分引起，一是围填造成过水断面减少，造成整个海湾纳潮水域减少，进而将导致海湾进出潮量减少；二是围填区域本身为浅滩，低潮时出露，人工岛占据部分潮棱体，围填后这部分潮棱体消失改变流场分布，引起两侧深槽流速增加，补偿进出潮量减少值，这是海湾反调节的效果。工程所在位置断面的进出潮量在围填小于等于某一临界规模时变化很小，如果围填超过该临界规模，进出潮量的变化将明显增大。

本实施例以泉州湾秀涂人工岛为例，对本发明的技术方案作进一步描述。

泉州湾位于晋江、洛阳江入海口，向东面向台湾海峡开敞，口门附近中部有大、小坠岛。泉州湾最大潮差超过6m，为强潮海湾，湾口多年平均波高0.9m、最大波高2～4m，年均波浪作用频率超过85％，波浪潮流共同作用是泉州湾泥沙运动的主要动力。泉州湾属构造成因海湾，湾内西侧发育有大面积的潮滩，占泉州湾总面积近一半。泉州湾东侧为鞋沙浅滩，北侧为秀涂角外伸滩涂，南侧为石湖角，大坠岛以东为湾外开敞海区。泉州湾岸线曲折，湾内水深多在0～-20m之间，主要等深线与海岸线走向基本平行，局部地形变化剧烈。鞋沙浅滩南北各有一条水道，南水道深窄，最深-21.3m。现，说明本发明的具体实施步骤。

本发明方法流程如图1所示，具体包括如下步骤：

步骤一、根据作业区水下地形实测资料及相关水文资料，选取坐标系和基本控制方程，建立二维平面潮流数学模型并进行数值求解；

所述基本控制方程包括水流连续方程、动量方程；

其中，所述二维平面潮流数学模型采用随地貌变化和潮位变化的变糙率建立，所述变糙率采用基本糙率和水深调节糙率相加的方法确定，基本糙率采用0.015～0.025；水下地形变化较大处根据水深调节糙率，即计算糙率a₁为基础糙率，a₂为调节糙率，h为水深。

二维平面潮流数学模型的基本控制方程如下：

水流连续方程：

ξ方向动量方程：

η方向动量方程：

式中，t表示时间；g为重力加速度；ρ表示海水密度；ξ、η分别表示正交曲线坐标系中二个正交曲线坐标；u、v分别表示沿ξ、η方向的流速；h表示水深；H表示水位；C_ξ、C_η表示正交曲线坐标系中的拉梅系数；

σ_ξξ、σ_ξη、σ_ηξ、σ_ηη表示紊动应力；

ν_t表示紊动粘性系数，ν_t＝C_μk²/ε，采用k-ε模型计算ν_t；v_t＝αu*h，α＝0.5～1.0，u*表示摩阻流速；

波浪底质点速度u_w,v_w采用式(6)计算：

L_w表示波长，T_w表示波浪周期，H_w表示波高，u_w,v_w分别表示波浪底质点速度，B表示波浪与潮流相互影响系数，当波、流同向时B＝0.917，当两者互相垂直时B＝-0.198，当方向不定时，B＝0.359；

S_ξξ、S_ξη、S_ηξ和S_ηη表示波浪辐射应力张量的四个分量，表达式为：

其中，H_w、θ分别表示波高和波向，C、C_g分别表示波速和波群速。

通过海域地质资料、动力地貌特征和海床演变等分析，人工岛应选取地质基础稳定、略呈淤积态势的海床；为最大限度地因势利导，促淤方案应避开深槽，走向应顺应涨落潮主流向。针对该河口湾岸线曲折、潮汐强进出潮量巨大、潮滩宽广地貌复杂、旋转流作用下的潮流强、受波浪影响大、泥沙运动受波流共同作用的独特动力地貌特点，如图2-4所示，拟围填人工岛3.23km²，建立大范围的台湾海峡数学模型和小范围的泉州湾数学模型，台湾海峡水动力数学模型水域面积约86960km²，共布置49752个单元、42475个结点，网格间距20～4500m。泉州湾波流泥沙数学模型口内至晋江泉州大桥，包括洛阳江，东至崇武，南至深沪湾，外海为-50m等深线，水域面积约1714km²。针对水下地形滩槽交替、复杂多变的特点，采用贴体坐标系统。利用现场实测资料对数学模型进行验证。选取2016年实测全潮作为代表潮型，利用潮流数学模型计算促淤工程区附近水动力场；

潮位变化导致泉州湾内水边界形态变化很大，特别对于潮流流速和流向影响较大。由于泉州湾低潮时出露大面积的潮滩，因此低潮位时采用较大的基本糙率，高潮位时采用较小的基本糙率。经调试，高低潮位时的基本糙率取0.013和0.015，高低潮位之间的糙率模型将通过插值给出。

步骤二、利用波浪数学模型(王红川，2013)计算不同潮位下代表波的辐射应力；根据实测资料，分析得出代表波(H_1/10＝1.45m，T＝3.8s)，根据波浪数学模型计算代表波在高、中、低潮位下的波要素，并通过线性插值得出其他潮位下的波要素。根据波浪要素计算波浪应力。

步骤三、优化人工岛平面形态；根据水动力数学模型模拟结果，优化人工岛平面形态；

以尽可能不影响天然流态作为人工岛平面形态设计的重要原则，利用建立的水动力数学模型优化人工岛平面形态，形成推荐的人工岛平面形态，以潮流矢量在潮量断面的法向分量计算某一时刻的潮量，一个潮周期中潮量的累加为进出潮量。根据潮位升降，分涨落潮量分别统计。根据工程前后主要断面进出潮量变化，分析人工岛的水动力环境影响。秀涂人工岛因围填占据过水面积，造成泉州湾口断面进出潮量减少2％左右；围填面积约3.33km²，以潮差4.27m计算，这部分潮棱体约0.14亿m3，如果不考虑港区开挖，工程区断面进出潮量应减少0.14亿m3，而工程断面计算涨潮平均进潮量仅减少0.015亿m3，落潮平均出潮量仅减少0.013亿m3，小于围填潮棱体值。这也导致工程海域内侧断面CS2涨落潮平均进出潮量变化很小。从这个角度看，该人工岛是适宜的。

步骤四、潮流、波浪、径流输沙耦合模拟；根据入湾径流实测输沙量，分洪水、中水和枯水统计各自平均含沙量。根据多年实测资料统计，晋江汛期(4～9月)径流量约占全年径流量75％，枯期(10月～翌年3月)径流约占全年径流量25％。汛期含沙量占全年含沙量70％，枯期含沙量占全年含沙量30％。2007～2014年平均汛期径流量和含沙量分别为214m³/s、0.25kg/m³，平均枯期径流量和含沙量分别为72m³/s、0.11kg/m³。根据实测资料统计，波浪年作用频率约85％。采用加权法求总冲淤量。径流输沙与潮流耦合计算采用洪、中、枯季加权所得。径流、潮流耦合后的结果再叠加波浪作用，分径流、潮流作用以及径流潮流叠加波浪作用分别模拟，最终的模拟结果进行加权求和，得出总的冲淤量。潮流计算条件采用2016年实测大中小潮过程循环计算，经试算，地貌加速因子因子(Roelvink,2006)为48，模拟过程中采用泥沙粒径分组。通过对比逐年地形变化，当地形变化量小于一定值(如10^{- 3}m)，可认为工程区附近地形冲淤达到平衡。

含沙量数学模型的基本控制方程如下：

悬沙不平衡输移方程：

非均匀悬移质按其粒径大小分成n₀组，S_L表示第L组粒径含沙量，P_SL表示此粒径悬沙含沙量所占的比值，则：

针对非均匀悬移质中第L组粒径的含沙量，二维悬移质不平衡输沙基本方程为：

式中，表示第L组泥沙的挟沙能力，ω_L为第L组泥沙的沉速；K₀为挟沙能力系数，α_L为第L组泥沙的含沙量恢复饱和系数；

海域泥沙沉速受含氯度影响，含沙量验证计算表明，用式(12)计算沉速，当ω_L小于0.015cm/s时，ω_L取絮凝后沉速0.015cm/s，当ω_L超过0.015cm/s时，则采用式(12)计算值；

床沙级配方程：

其中，E_m表示混合层厚度；P_mL0表示原始床沙级配，P_mL表示床沙级配；当混合层在冲刷过程中波及到原始底床时，ε₁＝0，否则ε₁＝1；

床底变形方程：

床底总冲淤厚度

步骤五、统计回淤量；采用GIS方法计算拟关心区面积，统计回淤量，并进行滩槽稳定评价及机理分析。

方案实施后秀涂作业区港池航道内平均回淤0.19m，回淤量为66.7万m³；泥沙回淤主要是由代表波作用下的泥沙运动所引起，仅潮流作用下的回淤很少。从位置来看，人工岛前沿回淤强度最大，往南侧回淤相对较少。根据水下地形高程参数、潮流流速流向参数、工程地质参数，确定敏感区为水下浅滩、深水航道以及已建的跨海大桥。人工岛实施后，深水航道一线的跨海大桥附近相对冲刷0.1～0.4m，人工岛西南侧头部出现冲刷，冲刷幅度约0.1～0.9m。淤积主要出现在工作船泊位区、人工岛西端北侧、人工岛与大坠岛之间、大小坠岛之间、小坠岛附近以及石湖西部。工作船泊位区、人工岛东部与大坠岛之间普遍淤积在1m以上，原因在于受人工岛建设影响而形成缓流区。其他区域淤积幅度总体小于0.1m，原因在于人工岛引起的流速变化主要在工程附近，对大范围滩槽稳定影响较小。

本发明克服了以往技术中缺乏对强潮河口湾人工岛作业区回淤模拟的问题。本发明波流泥沙数学模型构建过程中坐标系选取、对泥沙和波浪的处理完全针对岸线曲折、潮汐强进出潮量巨大、潮滩宽广地貌复杂、旋转流作用下的潮流强、受波浪影响大、泥沙运动受波流共同作用的独特动力地貌特点，考虑随地貌变化和潮位变化的变糙率以及不同季节径流输沙的影响，能够适应与满足不同地貌形态的半封闭强潮河口湾离岸人工岛围填适宜性判断要求。

以上所述仅为本发明的较佳实施例而已，凡在本发明精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (4)

1.一种强潮河口湾人工岛作业区港池航道回淤模拟方法，其特征在于，包括如下步骤：

步骤一、根据作业区水下地形实测资料及相关水文资料，选取坐标系和基本控制方程，建立二维平面潮流数学模型并进行数值求解；

所述基本控制方程包括水流连续方程、动量方程；

步骤二、建立波浪数学模型，计算不同潮位下代表波的辐射应力；根据实测资料，换算得出不同浪向所占比重，根据波浪数学模型计算所考虑的不同浪向在高、中、低潮位下的波要素，并通过线性插值得出其他潮位下的波要素，根据波浪要素计算波浪应力；

步骤三、根据数学模型模拟结果，优化人工岛平面形态；

步骤四、根据入湾径流实测输沙量，分洪水、中水和枯水统计各自平均含沙量；径流输沙与潮流耦合计算采用洪、中、枯季加权所得；径流、潮流耦合后的结果再叠加波浪作用，分径流、潮流作用以及径流潮流叠加波浪作用分别模拟，最终的模拟结果进行加权求和；

步骤五、利用GIS工具计算拟关心区面积，统计回淤量，并进行滩槽稳定评价及机理分析。

2.根据权利要求1所述的方法，其特征在于，所述步骤一中，所述二维平面潮流数学模型采用随地貌变化和潮位变化的变糙率建立，所述变糙率采用基本糙率和水深调节糙率相加的方法确定，基本糙率采用0.015～0.025；水下地形变化较大处根据水深调节糙率，即计算糙率a₁为基础糙率，a₂为调节糙率，h为水深。

3.根据权利要求1或2所述的方法，其特征在于，所述步骤一中，二维平面潮流数学模型的基本控制方程如下：

水流连续方程：

ξ方向动量方程：

η方向动量方程：

式中，t表示时间；g为重力加速度；ρ表示海水密度；ξ、η分别表示正交曲线坐标系中二个正交曲线坐标；u、v分别表示沿ξ、η方向的流速；h表示水深；H表示水位；C_ξ、C_η表示正交曲线坐标系中的拉梅系数；

σ_ξξ、σ_ξη、σ_ηξ、σ_ηη表示紊动应力；

v_t表示紊动粘性系数，v_t＝C_μk²/ε，采用k-ε模型计算v_t；v_t＝αu*h，α＝0.5～1.0，u*表示摩阻流速；

波浪底质点速度u_w,v_w采用式(6)计算：

L_w表示波长，T_w表示波浪周期，H_w表示波高，u_w,v_w分别表示波浪底质点速度，B表示波浪与潮流相互影响系数，当波、流同向时B＝0.917，当两者互相垂直时B＝-0.198，当方向不定时，B＝0.359；

S_ξξ、S_ξη、S_ηξ和S_ηη表示波浪辐射应力张量的四个分量，表达式为：

其中，H_w、θ分别表示波高和波向，C、C_g分别表示波速和波群速。

4.根据权利要求1所述的方法，其特征在于，所述步骤四中，含沙量数学模型的基本控制方程如下：

悬沙不平衡输移方程：

非均匀悬移质按其粒径大小分成n₀组，S_L表示第L组粒径含沙量，P_SL表示此粒径悬沙含沙量所占的比值，则：

S_L＝P_SLS，

针对非均匀悬移质中第L组粒径的含沙量，二维悬移质不平衡输沙基本方程为：

式中，表示第L组泥沙的挟沙能力，ω_L为第L组泥沙的沉速；K₀为挟沙能力系数，α_L为第L组泥沙的含沙量恢复饱和系数；

海域泥沙沉速受含氯度影响，含沙量验证计算表明，用式(12)计算沉速，当ω_L小于0.015cm/s时，ω_L取絮凝后沉速0.015cm/s，当ω_L超过0.015cm/s时，则采用式(12)计算值；

床沙级配方程：

其中，E_m表示混合层厚度；P_mL0表示原始床沙级配，P_mL表示床沙级配；当混合层在冲刷过程中波及到原始底床时，ε₁＝0，否则ε₁＝1；

床底变形方程：

床底总冲淤厚度