CN105956343B - 潮沟边壁冲刷及塌落过程的一维模拟方法 - Google Patents

Abstract

本发明公开了一种潮沟边壁冲刷及塌落过程的一维模拟方法，该方法包括：根据潮流模型参数建立一维潮流浅水方程，得到流速和潮位；根据流速和潮位模拟潮流对潮沟壁的侵蚀掏空过程，获取潮流侵蚀掏空后的潮沟壁边界线；根据边界线将潮流侵蚀掏空后的潮沟壁三角网格剖分为多个网格单元；采用有限单元法求解弹性力学平衡方程得到各网格单元的应力；采用摩尔‑库伦破坏准则判断剖分的网格单元是否发生剪切破坏和拉破坏，若是则修正应力判断发生剪切破坏和拉破坏破坏区是否贯穿，若是则模拟结束，否则返回求解流速和潮位。本发明可以模拟出潮沟壁的破坏过程、判断潮沟壁失稳的破坏类型、确定破坏区。

Description

潮沟边壁冲刷及塌落过程的一维模拟方法

技术领域

本发明涉及土木工程领域，尤其涉及一种潮沟边壁冲刷及塌落过程的一维模拟方法。

背景技术

潮沟壁的失稳涉及两个过程，首先是潮流对潮沟壁下部的侵蚀掏空，形成悬臂；随之悬臂在自身重力作用下失稳塌落。因其频繁发生在潮沟迁移摆动过程中，对潮沟演变有着特殊意义，进而影响着潮滩的稳定。

目前国内外对于潮沟边壁稳定性计算并没有成熟的方法，而对于河岸边壁的稳定性计算主要有3种方法：1)认为破坏类型为剪破坏，考虑断面安全系数；2)认为破坏过程为绕轴崩塌，考虑悬臂的临界宽度；3)计算岸壁内应力分布，从而确定破坏类型、破坏过程及破坏位置。

发明内容

发明目的：本发明针对现有技术存在的问题，提供一种潮沟边壁冲刷及塌落过程的一维模拟方法，可以模拟出潮沟壁的破坏过程、判断潮沟壁失稳的破坏类型、确定破坏区。

技术方案：本发明所述的潮沟边壁冲刷及塌落过程的一维模拟方法包括：

根据潮流模型参数建立一维潮流浅水方程，得到流速和潮位；

根据流速和潮位模拟潮流对潮沟壁的侵蚀掏空过程，获取潮流侵蚀掏空后的潮沟壁边界线；

根据边界线将潮流侵蚀掏空后的潮沟壁三角网格剖分为多个网格单元；

采用有限单元法求解弹性力学平衡方程得到各网格单元的应力；

采用摩尔-库伦破坏准则判断剖分的网格单元是否发生剪切破坏和拉破坏，若是则修正应力；

判断发生剪切破坏和拉破坏破坏区是否贯穿，若是则模拟结束，否则返回求解流速和潮位。

进一步的，所述根据潮流模型参数建立一维潮流浅水方程，得到流速和潮位，具体包括：

获取潮流模型参数，具体包括潮差A、曼宁系数n、滩面坡度i、外海水深H；

采用有限差分法求解一维浅水方程，得出流速u和潮位ξ，靠陆一侧采用干湿网格边界模拟潮间带的周期性露滩；其中，一维浅水方程为：

式中，h为总水深，为潮位和滩面高程之差，C为谢才系数，

进一步的，所述根据流速和潮位模拟潮流对潮沟壁的侵蚀掏空过程，获取潮流侵蚀掏空后的潮沟壁，具体包括：

采用Partheniades公式计算得到潮流作用下的潮沟壁侧侵蚀速率ε；其中，Partheniades公式为：

式中，Q为泥沙侵蚀通量，ρ_s为泥沙密度，M是侵蚀系数，a为经验系数，取值为1，τ_c为临界冲刷切应力，τ_b为水流的切应力，且τ_b＝ρc_Du²，ρ为海水密度，c_D为拖曳系数，u为流速；

在潮沟壁垂向设立N个监测点，每一次求解出流速u和潮位ξ时，对每一个监测点，判断其侵蚀与否；

当判定为受到侵蚀时，计算出该监测点侵蚀后退距离，并赋予该监测点侵蚀后的新坐标；

在一个潮周期后，将赋予新坐标的N个监测点重新连接作为潮流侵蚀掏空后潮沟壁的边界线。

进一步的，监测点判断侵蚀的方法为：若潮位淹没该监测点且水流提供的切应力τ_b大于临界冲刷切应力，则判定为收到侵蚀。

进一步的，所述根据边界线将潮流侵蚀掏空后的潮沟壁三角网格剖分为多个网格单元，具体包括：

将初始潮沟壁网格划分为多个网格单元，分别记录单元编号、节点编号及坐标、单元所包含的节点编号，其中，所述节点是指三角网格单元的顶点；

寻找新生成的潮沟壁边界线外的节点，并删除包含该节点的网格单元；

寻找周围网格单元数少于2个的节点，并删除包含该节点的网格单元，完成潮流侵蚀掏空后的潮沟壁网格剖分。

进一步的，所述采用有限单元法求解弹性力学平衡方程得到各网格单元的应力，具体包括：

输入土力学参数，具体包括土体重度γ、杨氏弹性模量E、泊松比μ、粘聚力c、内摩擦角φ；

计算作用于潮流侵蚀掏空后的潮沟壁的外荷载，具体包括开挖荷载P和静水压力P_w；其中，

P_W＝ρgh

式中，B_e为网格单元的应变矩阵，σ_e为网格单元的应力，N_e为网格单元的形函数矩阵，f_e为网格单元重力，Ne为剩余网格单元总数，Ve为每一个网格单元体积，ρ为海水密度，g为重力加速度，h为总水深；

根据开挖荷载P和静水压力P_w采用有限单元法求解弹性力学平衡方程得到各网格单元的应力。

进一步的，所述采用摩尔-库伦破坏准则判断剖分的网格单元是否发生剪切破坏和拉破坏，若是则修正应力，具体包括：

采用摩尔-库伦破坏准则判断各网格单元是否发生剪切破坏，若是，将网格单元类型设置为剪切破坏单元，并修正应力，以及将网格单元弹性模量设为小于第一阈值的值；其中，剪切破坏单元应力修正公式为：

式中，σ′_x为修正后的正应力，τ′_xy为修正后作用于垂直潮沟壁方向的平面的剪应力，τ_xy为作用于垂直潮沟壁方向的平面的剪应力，σ_x为垂直于潮沟壁方向的正应力，σ_y为沿重力方向的正应力，σ₁为大主应力，σ₃为小主应力，为内摩擦角；

采用摩尔-库伦破坏准则的负半轴判断各网格单元是否发生拉破坏，若是，将单元类型设置为拉破坏单元，并修正应力，以及将网格单元弹性模量设为小于第一阈值的值；其中，拉破坏单元应力修正公式为：

其中，σ′_x为修正后的正应力，τ′_xy为修正后作用于垂直潮沟壁方向的平面的剪应力。

有益效果：本发明与现有技术相比，其显著优点是：本发明提供的潮沟壁稳定计算方法，结合河口海岸动力学与土力学，综合考虑潮流对潮沟壁下部的侵蚀掏空过程和潮沟壁在自身重力作用下的坍塌过程，可以模拟出潮沟壁的破坏过程、判断潮沟壁失稳的破坏类型、确定破坏区，为潮沟壁稳定计算提供新方法，促进了潮沟蜿蜒摆动研究的发展。

附图说明

图1是本发明的一个实施例的流程示意图；

图2是本发明实施例的初始潮沟壁网格剖分及潮流侵蚀掏空后潮沟壁的网格剖分图；其中：图2(A)为初始潮沟壁网格划分图；图2(B)～(C)为潮流侵蚀掏空后潮沟壁网格剖分图，其中：虚线—潮流侵蚀掏空后的潮沟壁边界线；灰色线—优化前的网格边界线；粗线—优化后的网格边界线；

图3是本发明实施例的潮沟壁破坏过程、破坏区及现场观测图；其中：图3(A)～(E)分别为模拟得出的第1、2、4、6、8个潮周期后潮沟壁破坏方式及破坏区；图3(F)～(G)为现场实拍的潮沟壁破坏图。

具体实施方式

如图1所示，本实施例的潮沟边壁冲刷及塌落过程的一维模拟方法包括步骤S1-S5，其中，步骤S1和S2是模拟潮沟壁下部的侵蚀掏空过程，步骤S3和S4和S5模拟潮沟壁在自身重力作用下的坍塌过程。本实施例的具体步骤为：

S1、根据潮流模型参数建立一维潮流浅水方程，得到流速和潮位。

该步骤具体包括：(1)获取潮流模型参数，具体包括潮差A、曼宁系数n、滩面坡度i、外海水深H；(2)采用有限差分法求解一维浅水方程，得出流速u和潮位ξ，靠陆一侧采用干湿网格边界模拟潮间带的周期性露滩；其中，一维浅水方程为：

式中，h为总水深，为潮位和滩面高程之差，C为谢才系数，

例如，潮流模型参数可以如表1所示；

表1 潮流模型参数

S2、根据流速和潮位模拟潮流对潮沟壁的侵蚀掏空过程，获取潮流侵蚀掏空后的潮沟壁边界线。

该步骤具体包括：(1)采用Partheniades公式计算得到潮流作用下的潮沟壁侧侵蚀速率ε；(2)在潮沟壁垂向设立N个监测点，具体可以为30个，每一次求解出流速u和潮位ξ时，对每一个监测点，判断其侵蚀与否；(3)当判定为受到侵蚀时，计算出该监测点侵蚀后退距离，并赋予该监测点侵蚀后的新坐标；(4)在一个潮周期后，将赋予新坐标的N个监测点重新连接作为潮流侵蚀掏空后潮沟壁的边界线。其中，监测点判断侵蚀的方法为：若潮位淹没该监测点且水流提供的切应力τ_b大于临界冲刷切应力，则判定为收到侵蚀。

其中，Partheniades公式为：

式中，Q为泥沙侵蚀通量，ρ_s为泥沙密度，M是侵蚀系数，a为经验系数，取值为1，τ_c为临界冲刷切应力，τ_b为水流的切应力，且τ_b＝ρc_Du²，ρ为海水密度，c_D为拖曳系数，u为流速。具体在实施例，可以设置为：M为0.001m³/N·s，临界启动应力τ_c为0.1N/m²，经验系数a为1，海水密度ρ为1024kg/m³，拖曳系数c_D为0.002。

S3、根据边界线将潮流侵蚀掏空后的潮沟壁三角网格剖分为多个网格单元。

该步骤具体包括：(1)将初始潮沟壁网格划分为多个网格单元，分别记录单元编号、节点编号及坐标、单元所包含的节点编号，其中，所述节点是指三角网格单元(即三角形)的顶点；(2)寻找新生成的潮沟壁边界线外的节点，并删除包含该节点的网格单元；(3)寻找周围网格单元数少于2个的节点，并删除包含该节点的网格单元，完成潮流侵蚀掏空后的潮沟壁网格剖分。

例如，当初始潮沟壁尺寸为宽20m，高3m，初始网格划分为共2623个单元、1412个节点，如图2(A)所示；图2(B)～(C)给出了网格剖分示意图。

S4、采用有限单元法求解弹性力学平衡方程得到各网格单元的应力。

其中，该步骤具体包括：(1)输入土力学参数，具体包括土体重度γ、杨氏弹性模量E、泊松比μ、粘聚力c、内摩擦角φ；例如，土力学参数如表2所示。(2)计算作用于潮流侵蚀掏空后的潮沟壁的外荷载，具体包括开挖荷载P和静水压力P_w；其中，

P_W＝ρgh

式中，B_e为网格单元的应变矩阵，σ_e为网格单元的应力，N_e为网格单元的形函数矩阵，f_e为网格单元重力，Ne为剩余网格单元总数，Ve为每一个网格单元体积，ρ为海水密度，g为重力加速度，h为总水深；(3)根据开挖荷载P和静水压力P_w采用有限单元法求解弹性力学平衡方程得到各网格单元的应力。

表2 初始岩土力学参数

S5、采用摩尔-库伦破坏准则判断剖分的网格单元是否发生剪切破坏和拉破坏，若是则修正单元应力。

该步骤具体包括：

(1)采用摩尔-库伦破坏准则判断各网格单元是否发生剪切破坏，若是，将网格单元类型设置为剪切破坏单元，并修正单元应力，以及将网格单元弹性模量设为小于第一阈值的值；其中，剪切破坏单元应力修正公式为：

式中，σ′_x为修正后的正应力，τ′_xy为修正后作用于垂直潮沟壁方向的平面的剪应力，τ_xy为作用于垂直潮沟壁方向的平面的剪应力，σ_x为垂直于潮沟壁方向的正应力，σ_y为沿重力方向的正应力，σ₁为大主应力，σ₃为小主应力，为内摩擦角；

(2)采用摩尔-库伦破坏准则的负半轴判断各网格单元是否发生拉破坏，若是，将单元类型设置为拉破坏单元，并修正单元应力，以及将网格单元弹性模量设为小于第一阈值的值；其中，拉破坏单元应力修正公式为：

其中，σ′_x为修正后的正应力，τ′_xy为修正后作用于垂直潮沟壁方向的平面的剪应力。

S6、判断发生剪切破坏和拉破坏破坏区是否贯穿，若是则模拟结束，否则返回求解流速和潮位。

结果如图3(A)～(E)所示，分别为第1、2、4、6、8个潮周期后潮沟壁破坏方式及破坏区，从顶部开裂到形成贯穿的破坏区共经历了8个潮周期，本发明模拟出的潮沟壁破坏方式与现场观测的结果相吻合。

以上所揭露的仅为本发明一种较佳实施例而已，不能以此来限定本发明之权利范围，因此依本发明权利要求所作的等同变化，仍属本发明所涵盖的范围。

Claims (6)

1.一种潮沟边壁冲刷及塌落过程的一维模拟方法，其特征在于该方法包括：

根据潮流模型参数建立一维潮流浅水方程，得到流速和潮位；

根据流速和潮位模拟潮流对潮沟壁的侵蚀掏空过程，获取潮流侵蚀掏空后的潮沟壁边界线；

根据边界线将潮流侵蚀掏空后的潮沟壁三角网格剖分为多个网格单元；

采用有限单元法求解弹性力学平衡方程得到各网格单元的应力；

采用摩尔-库伦破坏准则判断剖分的网格单元是否发生剪切破坏和拉破坏，若是则修正应力；

判断发生剪切破坏和拉破坏破坏区是否贯穿，若是则模拟结束，否则返回求解流速和潮位；

其中，所述根据流速和潮位模拟潮流对潮沟壁的侵蚀掏空过程，获取潮流侵蚀掏空后的潮沟壁，具体包括：

采用Partheniades公式计算得到潮流作用下的潮沟壁侧侵蚀速率ε；其中，Partheniades公式为：

式中，Q为泥沙侵蚀通量，ρ_s为泥沙密度，M是侵蚀系数，a为经验系数，取值为1，τ_c为临界冲刷切应力，τ_b为水流的切应力，且τ_b＝ρc_Du²，ρ为海水密度，c_D为拖曳系数，u为流速；

在潮沟壁垂向设立N个监测点，每一次求解出流速u和潮位ξ时，对每一个监测点，判断其侵蚀与否；

当判定为受到侵蚀时，计算出该监测点侵蚀后退距离，并赋予该监测点侵蚀后的新坐标；

在一个潮周期后，将赋予新坐标的N个监测点重新连接作为潮流侵蚀掏空后潮沟壁的边界线。

2.根据权利要求1所述的潮沟边壁冲刷及塌落过程的一维模拟方法，其特征在于：所述根据潮流模型参数建立一维潮流浅水方程，得到流速和潮位，具体包括：

获取潮流模型参数，具体包括潮差A、曼宁系数n、滩面坡度i、外海水深H；

采用有限差分法求解一维浅水方程，得出流速u和潮位ξ，靠陆一侧采用干湿网格边界模拟潮间带的周期性露滩；其中，一维浅水方程为：

式中，g为重力加速度，h为总水深，为潮位和滩面高程之差，C为谢才系数，

3.根据权利要求1所述的潮沟边壁冲刷及塌落过程的一维模拟方法，其特征在于：监测点判断侵蚀的方法为：若潮位淹没该监测点且水流提供的切应力τ_b大于临界冲刷切应力，则判定为收到侵蚀。

4.根据权利要求1所述的潮沟边壁冲刷及塌落过程的一维模拟方法，其特征在于：所述根据边界线将潮流侵蚀掏空后的潮沟壁三角网格剖分为多个网格单元，具体包括：

将初始潮沟壁网格划分为多个网格单元，分别记录单元编号、节点编号及坐标、单元所包含的节点编号，其中，所述节点是指三角网格单元的顶点；

寻找新生成的潮沟壁边界线外的节点，并删除包含该节点的网格单元；

寻找周围网格单元数少于2个的节点，并删除包含该节点的网格单元，完成潮流侵蚀掏空后的潮沟壁网格剖分。

5.根据权利要求1所述的潮沟边壁冲刷及塌落过程的一维模拟方法，其特征在于：所述采用有限单元法求解弹性力学平衡方程得到各网格单元的应力，具体包括：

输入土力学参数，具体包括土体重度γ、杨氏弹性模量E、泊松比μ、粘聚力c、内摩擦角φ；

计算作用于潮流侵蚀掏空后的潮沟壁的外荷载，具体包括开挖荷载P和静水压力P_w；其中，

式中，B_e为网格单元的应变矩阵，σ_e为网格单元的应力，N_e为网格单元的形函数矩阵，f_e为网格单元重力，Ne为剩余网格单元总数，Ve为每一个网格单元体积，ρ为海水密度，g为重力加速度，h为总水深；

根据开挖荷载P和静水压力P_w采用有限单元法求解弹性力学平衡方程得到各网格单元的应力。

6.根据权利要求1所述的潮沟边壁冲刷及塌落过程的一维模拟方法，其特征在于：所述采用摩尔-库伦破坏准则判断剖分的网格单元是否发生剪切破坏和拉破坏，若是则修正应力，具体包括：

采用摩尔-库伦破坏准则判断各网格单元是否发生剪切破坏，若是，将网格单元类型设置为剪切破坏单元，并修正应力，以及将网格单元弹性模量设为小于第一阈值的值；其中，剪切破坏单元应力修正公式为：

式中，σ′_x为修正后的正应力，τ′_xy为修正后作用于垂直潮沟壁方向的平面的剪应力，τ_xy为作用于垂直潮沟壁方向的平面的剪应力，σ_x为垂直于潮沟壁方向的正应力，σ_y为沿重力方向的正应力，σ₁为大主应力，σ₃为小主应力，为内摩擦角，c为粘聚力；

采用摩尔-库伦破坏准则的负半轴判断各网格单元是否发生拉破坏，若是，将单元类型设置为拉破坏单元，并修正应力，以及将网格单元弹性模量设为小于第一阈值的值；其中，拉破坏单元应力修正公式为：

其中，σ′_x为修正后的正应力，τ′_xy为修正后作用于垂直潮沟壁方向的平面的剪应力。