CN113486557B - 一种基于内外模式的二、三维数学模型耦合模拟方法 - Google Patents

一种基于内外模式的二、三维数学模型耦合模拟方法 Download PDF

Info

Publication number
CN113486557B
CN113486557B CN202110785304.3A CN202110785304A CN113486557B CN 113486557 B CN113486557 B CN 113486557B CN 202110785304 A CN202110785304 A CN 202110785304A CN 113486557 B CN113486557 B CN 113486557B
Authority
CN
China
Prior art keywords
sand
water
model
dimensional
river reach
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110785304.3A
Other languages
English (en)
Other versions
CN113486557A (zh
Inventor
王建军
王晨阳
杨云平
刘晓菲
张宏千
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Tianjin Research Institute for Water Transport Engineering MOT
Original Assignee
Tianjin Research Institute for Water Transport Engineering MOT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Tianjin Research Institute for Water Transport Engineering MOT filed Critical Tianjin Research Institute for Water Transport Engineering MOT
Priority to CN202110785304.3A priority Critical patent/CN113486557B/zh
Publication of CN113486557A publication Critical patent/CN113486557A/zh
Application granted granted Critical
Publication of CN113486557B publication Critical patent/CN113486557B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/12—Simultaneous equations, e.g. systems of linear equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A10/00—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE at coastal zones; at river basins
    • Y02A10/40—Controlling or monitoring, e.g. of flood or hurricane; Forecasting, e.g. risk assessment or mapping

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Algebra (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Operations Research (AREA)
  • Computing Systems (AREA)
  • Fluid Mechanics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开一种基于内外模式的二、三维数学模型耦合模拟方法，包括，构建二维水沙数学模型及三维水沙数学模型；获取初始和边界条件，基于初始和边界条件，通过二维水沙数学模型进行外模式计算，基于外模式计算结果，通过三维水沙数学模型进行内模式计算，基于内模式计算结果，对二维水沙数学模型进行局部更新；通过二维水沙数学模型及三维水沙数学模型的内外模式的耦合计算，得到整体河段特征的模拟预测结果。本发明实现长河段高精度的复杂边界条件下河段特征的滩槽演变规律模拟预测，具有很强的实用性。

Description

一种基于内外模式的二、三维数学模型耦合模拟方法
技术领域
本发明涉及河段模拟技术领域，特别涉及一种基于内外模式的二、三维数学模型耦合模拟方法。
背景技术
近几十年来，人类活动的强烈干扰改变了河流自然输沙过程、输沙量和径流的年内分配，滩槽原有的形态及其演变趋势必然发生自适应调整，长江中下游原有的河势演变规律将被打破，部分河段深泓摆动频繁，河势演变更趋复杂，滩槽的冲淤特征与冲淤量、床面形态等发生改变。
随着长江干线过江通道建设数量逐渐增多，近10年建设的桥梁主桥墩均涉水，并处于活动的边心滩区域。随着桥梁建成时间的增长，桥墩周围形成局部冲刷坑，改变洲滩形态，而河道形态调整又会对水沙输移特性产生影响，使得未守护的洲滩演变规律和航道条件变化趋势更为复杂。
处于冲淤平衡或接近平衡的自然河流，受到了建筑物的约束和控制，改变了水流的方向或流速的大小，或限制了原有浅滩的变形，则原有的平衡状态被破坏，将引起有限范围内的局部冲刷。如闸墩附近，枢纽下游，丁坝坝头以及堤防工程的根部，护滩带周围等，都是容易发生局部冲刷的地方。长河段二维水沙数学模型在系统治理和上下游洲滩联动规律模拟方面具有优势，但桥墩或整治建筑物周围具有明显的三维水流特性，二维模型不能真实反映桥墩周围极限冲刷深度，由于桥墩周围床面形态模拟与实际情况不一致，进而影响桥区河段泥沙输移规律的模拟精度。
近些年来随着数学模型获得了飞速的发展，三维模型在局部绕流冲刷研究中获得了应用，并逐渐成为局部冲刷研究的重要手段之一，多名学者采用三维数学模型模拟了桥墩局部冲刷情况，得到了满意的结果。
中下游河床冲淤多变，桥梁工程和航道整治工程的实施，形成了新的河势约束边界，进一步影响桥区河段浅滩形态及航道尺度。桥梁工程与航道工程耦合作用下形成了复杂的控制边界条件，江中桥墩的存在一定程度上改变了泥沙输移路径，桥墩对于河床冲淤的影响长期存在，特别是周围局部冲刷坑形成后，桥区河段边心滩周期性变化规律将受到影响，导致现有模型对于河段特征的模拟结果不准确。
发明内容
为解决上述现有技术中所存在的对于河段特征模拟结果不准确的问题，本发明提供了一种基于内外模式的二、三维数学模型耦合模拟方法，包括：
S1.获取河段特征，基于河段特征，构建整体河段水沙模型及局部河段水沙模型，其中，整体河段水沙模型为二维水沙数学模型，局部河段水沙模型为三维水沙数学模型；
S2.获取整体河段水沙模型的初始和边界条件；
S3.基于整体河段水沙模型的初始和边界条件，对整体河段水沙模型进行外模式计算，得到外模式计算结果，其中，外模式为对整体河段特征进行模拟计算；
S4.基于外模式计算结果，对局部河段水沙模型进行内模式计算，得到内模式计算结果，其中，内模式为对局部河段特征进行模拟计算；
S5.基于内模式计算结果，对整体河段水沙模型进行局部更新；
S6.重复步骤S3-步骤S5进行迭代计算直到满足迭代条件，得到整体河段特征的模拟预测结果。
优选的，步骤S1中，构建整体河段水沙模型的具体过程包括，
基于河段特征，构建整体河段水沙模型的二维数学模型，并建立整体河段水沙模型的控制方程；
其中，整体河段水沙模型的控制方程包括水流连续方程，X方向水流运动方程、Y方向水流运动方程、悬沙不平衡输运方程、推移质不平衡输移方程和河床变形方程。
优选的，局部河段水沙模型为基于紊流模型构建的三维水沙数学模型。
优选的，步骤S2中，获取所述初始和边界条件的具体步骤包括，基于河段特征，获取整体河段水沙模型的初始条件；
其中，初始条件包括但不限于x方向流速、y方向流速、悬移质含沙量、底沙量；
基于整体河段水沙模型的初始条件，获取整体河段水沙模型的边界条件；
其中，整体河段水沙模型的边界条件包括上下游控制边界条件、固壁条件、收敛控制条件。
优选的，步骤S3中，对整体河段水沙模型进行外模式计算的具体步骤包括：
基于初始和边界条件，对整体河段水沙模型进行离散求解，并基于离散求解结果通过压力校正法进行压力～流速耦合求解，得到外模式计算结果；
其中，离散求解结果包括x方向流速、y方向流速、水位、悬移质含沙量、底沙量、河床变形；
外模式计算结果包括流速场、水位值及含沙量场。
优选的，对整体河段水沙模型进行离散求解的过程中，离散求解的方法使用有限体积法、有限差分法及有限元法其中一种。
优选的，步骤S4中，通过局部河段水沙模型进行内模式计算的具体步骤包括：
基于外模型计算结果，得到局部河段水沙模型的边界条件，基于局部河段水沙模型的边界条件，通过空间离散及对流扩散进行离散求解，得到河床变形；
基于外模型计算结果，通过推移质运动计算，得到局部河段泥沙变化量；
其中，内模式计算结果包括局部河段泥沙变化量及河床变形。
优选的，在通过推移质运动计算，得到河床变形之后，还包括对河床变形进行局部修正，具体步骤包括：
在通过推移质运动计算之后，对局部水沙模型进行网格扫描，基于扫描结果，对局部水沙模型进行局部调整，通过进行重复扫描及局部调整后，获取最终的修正结果，并基于修正结果，对河床变形进行修正。
优选的，步骤S6中，外模式计算结合内模式计算的迭代计算的具体过程包括，
获取系列水文年，对系列水文年进行划分，得到若干个试验年，对若干个试验年进行排序；
从第一个试验年开始，在每个试验年末，基于内模式计算结果，对整体河段水沙模型的初始条件进行更新，并通过重复步骤S3-S6进行迭代计算；
当对最后一个试验年末整体河段水沙的初始条件进行更新后，停止迭代计算。
本发明具有如下技术效果：
(1)本发明采用二、三维数学模型耦合计算，利用内模式和外模式分别得到局部和整体河段河床冲淤变化，提高计算效率和精度，进而实现整体河段高精度的复杂边界条件下河段特征的滩槽演变规律模拟预测。
(2)本发明通过对三维数学模型进行内模式计算对二维数学模型进行局部更新，可以实现兼顾局部河床冲淤调整和上下游之间联动规律的准确模拟。
(3)本发明中的二、三维数学模型通过所测的河段特征进行构建，可以用于内河整体河段系统治理效果模拟，也可用于桥区河段“桥梁工程+航道工程+水沙调控”综合作用下滩槽演变规律预测，具有很强的实用性。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
图1为本发明实施例提供的方法流程示意图；
图2为本发明实施例提供的模型结构示意图；
图3为本发明实施例提供的泥沙坍塌调整示意图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
为了解决在现有技术中存在河段特征的模拟结果不准确等问题，本发明提供了如下方案：
如图1、图2所述，本发明提供了一种基于内外模式的二、三维数学模型耦合模拟方法，包括：
S1.获取河段特征，基于河段特征，构建整体河段水沙模型即及局部河段水沙模型，其中，整体河段水沙模型为二维水沙数学模型，局部河段水沙模型为三维水沙数学模型；
构建整体河段水沙模型的具体过程包括，基于河段特征，构建整体河段水沙模型的二维数学模型，并建立整体河段水沙模型的控制方程；其中，整体河段水沙模型的控制方程包括水流连续方程，X方向水流运动方程、Y方向水流运动方程、悬沙不平衡输运方程、推移质不平衡输移方程和河床变形方程。
二维水流泥沙数学模型(即外模式)控制方程如下：平面二维水沙数学模型控制方程包括水流连续方程、X方向水流运动方程和Y方向水流运动方程、悬沙不平衡输运方程、推移质不平衡输移方程和河床变形方程。通过这些方程组成方程组进行联立求解，分别计算得到流速U、V，水位ζ、含沙量S，底沙N，以及河床变形η。
1)水流连续方程：
2)X方向和Y方向水流运动方程：
式中：ζ—水位(m)；t—时间(s)；h—静水深(m)；u、v—流速矢量沿X、Y方向的分量(m/s)；g—重力加速度(m/s2)；νe—水流紊动黏性系数(m2/s)；c—谢才系数(m1/2/s)；f—科氏参量(s-1)，ω为地球自转角速度，/>为地理纬度；ρ—水密度(kg/m3)。
3)悬沙不平衡输运方程
由沿水深积分，并假定由流速和含沙量沿垂线分布不均匀在积分时产生的修正系数：
引入冲淤平衡时的挟沙能力S*,得：
式中：S为单位水体垂线平均含沙量，s＝ρsc，c为单位水体体积浓度；νt＝νmt；σS＝σc为Schmidt数；ωS为泥沙沉速，下标i表示非均匀泥沙分组情况。
其中，
对于水面z＝ζ泥沙扩散通量为零边界条件：对于底部z＝-h泥沙扩散通量：/>一般认为悬沙粒径很细时，不论泥沙沿水深分布是否处于平衡状态，含沙量沿水深变化不大，上式表示为：Φs＝αωs(S*-S)
式中α＝α*Pr为系数，此表达式在泥沙输移数模计算中得到广泛运用。关于表达式中的系数α为恢复饱和系数，在数学模型计算中，垂线恢复饱和系数α取值范围为0.25～1.0，淤积状态取α＝0.25；冲刷状态取α＝1.0。
分组挟沙力Sn *计算按照窦国仁、赵士清的模式，将非均匀沙按其粒径大小分成N0组，Sn表示n组粒径得含沙量，Pn表示此粒径在悬沙总含沙量S中所占的比值：
总挟沙力：
挟沙力级配：
分组挟沙力：Sn *＝Pn *S*，式中：0<α<1，ωn为第n组粒径的沉速，ωm为非均匀平均沉速。
4)推移质不平衡输移方程
根据推移质不平衡非均匀输沙原理，通过推移质水深折算推导出底沙不平衡输沙输移方程：
对于非均匀沙，推移质不平衡输移方程采用如下形式：
式中，Ni、分别为推移质输沙量和推移质输沙能力折算成相应水深的泥沙浓度，βi为推移质泥沙恢复饱和系数，下标i表示第i组粒径泥沙对应的变量。ωsi为对于非均匀沙第i组粒径泥沙输沙率。
推移质输沙率的计算公式众多，目前较常用的有：VanRijn公式、窦国仁公式、岗恰洛夫公式等，对于非均匀沙第i组粒径泥沙输沙率根据Karim,Knenedy(1981)建议考虑隐蔽系数ηi：
总的输沙能力；/>Pbi为第i组粒径泥沙所占的百分比；为第i组粒径泥沙推移质输沙率。
窦国仁推移质输沙率公式为：
式中：Vki为第i组粒径泥沙临界启动流速：
式中：di为第i组泥沙粒径，γo为稳定干容重，γo＝1650kg/m3；εk为粘结力参数(天然沙εk＝2.56cm3/s2)；σ为薄膜水厚度，σ＝0.21×10-4cm，γ′o为床面泥沙干容重，对于细沙γ′o＝γo，C为谢才系数：Δ为床面糙度：
K为系数，对于沙质推移质k取0.01，这样对于非均匀沙，推移质不平衡输移方程中可写成：
局部河段水沙模型为基于紊流模型构建的三维水沙数学模型。紊流模型采用非平衡κ-ε模型模拟桥台绕流以及桥台局部冲刷坑内的流场。为了适应冲坑及障碍物同时存在这种复杂的边界，模型网格采用笛卡尔坐标系下非正交结构网格，空间离散采用有限体积法，对流扩散采用乘方格式离散，压力与速度耦合采用SIMPLE方法，离散的线性方程组采用强隐格式求解。模型出口控制水位取长河段二维水沙模型计算各级流量时的水位ζ，模型进口控制流速U、V由长河段二维水沙模型计算得到，模型进口沙量S和推移质输沙率qb均由长河段二维水沙数学模型提供。通过局部河段三维水沙数学模型，得到每个典型年(月)末冲刷后局部河床变形Δz三维。
在建立二维水沙数学模型与三维水沙数学模型过程中，需要对上述水沙数学模型进行网格划分，并剖分出计算区域，本发明计算对于河床地形的计算都是基于计算区域进行的计算。整体河段水沙模型为长河段二维水沙数学模型。
S2.获取整体河段水沙模型的初始和边界条件；
步骤S2中，所述获取初始和边界条件的具体步骤包括，基于河段特征，获取整体河段水沙模型的初始条件；其中，初始条件包括但不限于x方向流速、y方向流速、悬移质含沙量、底沙量；基于整体河段水沙模型的初始条件，获取整体河段水沙模型的边界条件；其中，整体河段水沙模型的边界条件包括上下游控制边界条件、固壁条件、收敛控制条件。数学模型在计算开始前，要预先给定整体河段水沙模型的初始条件和边界条件，具体方法如下：
①、初始条件
给定初始条件时刻t＝0时，计算域内所有计算变量(x方向流速U、y方向流速V、水位值ζ、悬移质含沙量Si、底沙量Ni)的初值，给出悬沙级配和分段床沙级配。
②、上、下游控制边界条件
上游进口条件：给定上游来流Qin、沙量Sin和进口推移质输沙率qb，进口各点流速Vi＝0。进口各控制点流速由下式迭代算出：
式中：Uj，hj为进口计算网格点沿y方向流速和水深，dyj为离散网格间距，Vj＝0。
下游出口条件：给定水位ζout。
③、固壁条件
流速采用非滑移边界条件，其边壁流速给定为零，即U＝V＝0；对于含沙量Si，底沙Ni计算中采用法向梯度为零条件：
④、收敛控制条件
控制连续方程最大质量源bmax和通过各断面流量Qj，
流速：水位：/>
S3.基于整体河段水沙模型的初始和边界条件，通过整体河段水沙模型进行外模式计算，得到外模式计算结果，其中，外模式为对整体河段特征进行模拟计算；
步骤S3中，对整体河段水沙模型进行外模式计算的具体步骤包括：基于初始和边界条件，对整体河段水沙模型进行离散求解，并基于离散求解结果通过压力校正法进行压力～流速耦合求解，得到外模式计算结果；其中，离散求解结果包括x方向流速、y方向流速、水位、悬移质含沙量、底沙量、河床变形；外模式计算结果包括流速场、水位值及含沙量场。对整体河段水沙模型进行离散求解的过程中，离散求解的方法使用有限体积法、有限差分法及有限元法其中一种。
为了考虑河段内上下游洲滩联动变化，首先，建立长河段二维水沙数学模型即外模式，采用水流连续方程、X方向和Y方向水流运动方程、悬沙不平衡输运方程，分别离散后，求得x方向流速U、y方向流速V、悬移质含沙量S、水位ζ，离散方法可采用有限体积法、有限差分法和有限元法等，得出正确的流速场和水位值，以及含沙量场。下面以有限体积法离散为例，说明水流运动方程(外模式)求解过程。具体离散求解过程如下，在任意一个给定的三角形控制体内对水流连续方程进行积分，并利用格林公式，得
式中：Ai为单元i的面积。
写成离散形式得
或
类似地，对X方向和Y方向水流运动方程进行积分，采用时间向前差分，得到
式中：为第i个单元的第j边上水平数值通量，采用基于黎曼解的Osher格式计算；Ai为第i个单元的面积，li,j为第i个单元的第j条边的边长；DhU、DhV分别为x、y方向水平紊动扩散项，其积分形式可表示为
对以上方程进行整理可得
通过上式，可求得水位ζ、x方向流速U、y方向流速V，通过类似离散方法，离散悬沙不平衡输运方程和对于非均匀沙的推移质不平衡输移方程，可求得S、N。
S4.基于外模式计算结果，通过局部河段水沙模型进行内模式计算，得到内模式计算结果，其中，内模式为对局部河段特征进行模拟计算；步骤S4中，通过局部河段水沙模型进行内模式计算的具体步骤包括：基于外模型计算结果，得到局部河段水沙模型的边界条件，基于局部河段水沙模型的边界条件，通过空间离散及对流扩散进行离散求解，得到河床变形；基于外模型计算结果，通过推移质运动计算，得到局部河段泥沙变化量；其中，内模式计算结果包括局部河段泥沙变化量及河床变形。
针对长河段内桥区局部开展局部三维水沙数学模型(内模式)计算，其中计算边界条件，如模型进口U、V、S悬移质含沙量S，模型出口ζ，均由长河段二维水沙数学模型(外模式)计算结果给定。通过内模式计算得到典型年或系列年试验过程中每个试验年(月)末桥墩周围河床变形Δz三维。典型水文年和系列水文年的选取应根据研究问题的需要从对滩槽形态和航道条件演变趋势不利的角度选取不同的水沙年组合作为系列水文年，系列水文年宜包含丰水丰沙、中水中沙和小水少沙等特征水文年。三维水沙数学模型(即内模式)计算方法如下：
利用最新实测资料，采用三维水流泥沙数学模型计算，研究局部河段泥沙输移规律和河床变形。
①悬沙数学模型基本方程
悬沙运动控制方程为三维对流扩散方程：
河床变形方程：
其中，AH和Kh分别为水平和垂向泥沙质量扩散系数(m2/s)，u,v,w分别为x、y、z向水流速度(m/s)，Ci为第i组份悬浮泥沙浓度(kg/m3)，ωsi为第i组份泥沙沉速(m/s)，ΔZ为床面冲淤厚度(m),qx，qy分别为X、Y向的底沙单宽输沙率[kg/(m·s)]。
②边界条件
自由水面要求含沙量的净通量为零，即在z＝ζ时，
底部边界条件表示为：
其中Di和Ei分别为第i组份泥沙淤积率和冲刷率。
冲刷率表示为：
Ei＝0,τb＜τei
其中，E0i为第i组份泥沙的床面冲刷强度，Pb是床面泥沙孔隙率，Fbi为第i组份泥沙所占比例，τb是床面剪切应力，τei是第i组份泥沙的临界冲刷应力。
泥沙的沉积作用由以下方程控制：
Di＝wsiCi
③推移质运动对局部河床变形的影响
推移质运动对局部河床冲淤的影响采用VanRijn(2007)公式计算。
其中，qb为推移质输沙率，ρsd为泥沙密度，ρw为水体密度，d50为泥沙中值粒径，fs为泥沙修正系数，D*为无量纲的泥沙粒径系数，τ′b为床面剪切力，τb,cr为床面泥沙起动临界剪切力，γ为系数。相关参数具体计算公式如下：
D*＝d50[(s-1)g/υ2]1/3
s＝ρs/ρw
τb,cr＝θcr(ρs-ρw)gd50
τ′b＝0.5ρwfcu2
fc＝8g[18log(12h/d90)]-2
其中，dsand＝0.062mm，dsand为泥沙中值粒径参数，υ为水体运动粘性，fc为水流产生的底部摩阻，d90为泥沙级配中以重量计90％较之小的泥沙粒径。
模型出口控制水位取长河段二维水沙模型计算各级流量时的水位ζ，模型进口控制流速U、V由长河段二维水沙模型计算得到，模型进口沙量S和推移质输沙率qb均由长河段二维水沙数学模型提供。
如图3所示，在通过推移质运动计算，得到河床变形之后，还包括对河床变形进行局部修正，具体步骤包括：
在通过推移质运动计算之后，对局部水沙模型进行网格扫描，基于扫描结果，对局部水沙模型进行局部调整，通过进行重复扫描及局部调整后，获取最终的修正结果，并基于修正结果，对河床变形进行修正。
在长河段二维水沙数学模型计算过程中，为了进一步提高桥墩周围泥沙输移的模拟精度，在地形替换的基础上，针对局部冲刷床面坡度较大的特殊情况，分析床面坡度对推移质输沙的影响，对已有的考虑坡度影响的推移质输沙率修正方法进行改进，在求解河床变形的过程中，除了考虑冲淤造成的床面变化外，还将对床面发生坍塌调整的过程进行处理。在斜坡床面上，除了作用在泥沙颗粒上的水流剪切力外，重力的切向分力也是推移质输运的一个动力因素。
为了考虑坍塌的影响，在数值模型中，对于非均匀沙的推移质不平衡输移方程中推移质输沙率计算结果进行局部修正。具体操作方法如下：根据泥沙模型计算完河床变形后，对床面网格进行扫描。每个时间步网格节点高程根据地貌模型调整后，都需要判断每个单元平面外法向方向与Z轴正向的夹角大小，当发现相邻角点间的倾角大于泥沙休止角时，要进行局部调整，调整模式如图所示，在扫描时发现角点A(xA，yA，zA)和角点B(xB，yB，zB)之间倾角大于休止角需要将较高的B点降至B′点，A点升A′点，使AB之间的倾角降到/>ηA、ηB表示如下：
式中，ηA、ηB分别为外模式中A、B点的河床冲淤变形。
经过对床面的一次扫描调整后可能会出现新的陡坡，因此需要对整个床面进行反复的扫描检查直至床面上所有相邻角点间的倾角均在休止角范围内。
S5.基于内模式计算结果，对整体河段水沙模型进行局部更新；
通过求解局部河段三维水沙数学模型及内模式计算中河床变形方程，得到每个典型年(月)末冲刷后局部河床变形Δz三维。
长河段二维水沙数学模型采用由悬移质冲淤引起的河床变形方程、由推移质冲淤引起的河床变形方程，分别离散后求得悬移质冲淤引起的河床变形ηsi和推移质冲淤引起的河床变形ηbi，再通过河床总的冲淤厚度方程，求得河床总的冲淤厚度η。
由悬移质冲淤引起的河床变形方程为：
式中：ηsi为第i组粒径悬移质泥沙引起的冲淤厚度。γ0为床面泥沙干容重。
由推移质冲淤引起的河床变形方程为：
式中：ηbi为第i组粒径推移质泥沙引起的冲淤厚度。
河床总的冲淤厚度方程：
长河段二维水沙数学模型和局部河段三维水沙数学模型同时进行计算，为了提高计算效率，采用不同的时间积分步长来计算求解，内、外模式计算采用的典型水文年与系列水文年一致。采用由悬移质冲淤引起的河床变形方程、由推移质冲淤引起的河床变形方程，分别离散后求得悬移质冲淤引起的河床变形ηsi和推移质冲淤引起的河床变形ηbi，再通过河床总的冲淤厚度方程，求得河床总的冲淤厚度η。在内、外模式耦合计算过程中，为保证局部河床变形数据替换完成，应保证每年(月)末内模式对桥墩周围河床形态的计算结果Δz三维要早于外模式计算结果η。
S6.重复步骤S3-步骤S5进行迭代计算直到满足迭代条件，得到整体河段特征的模拟预测结果。步骤S6中，外模式计算结合内模式计算的迭代计算的具体过程包括，获取系列水文年，对系列水文年进行划分，得到若干个试验年，对若干个试验年进行排序；在获取系列水文年过程中，选取对于滩槽演变不利的特殊水沙系列和大、中、小水年，作为长系列水文年，对于内河航道整治一般选取10年长系列，并对来水来沙情况进行概化，一般情况7～10天概化为一个流量级。以时间中年或者月为跨度进行划分，划分出试验年，按照时间顺序对试验年进行排序。从第一个试验年开始，在每个试验年末，基于内模式计算结果，对整体河段水沙模型的初始条件进行更新，并通过重复步骤S3-S6进行迭代计算；当对最后一个试验年末整体河段水沙的初始条件进行更新后，停止迭代计算。
外模式中桥墩周围局部地形修正。采用内模式计算的每个试验年(月)末桥墩周围河床变形Δz三维数据替换外模式计算结果η。即通过上述步骤建立长河段二维水沙数学模型和局部河段三维水沙数学模型后，给定内、外模式初始和边界条件，开始计算，在计算过程中，每隔7～10天即每隔一个流量级对河床地形进行一次更新，即用新时刻地形(ηt)替换上一时刻地形(ηt-1)，当计算至系列年第一年(月)末时，得到第一年(月)末的河床冲淤厚度η(第一年或月末)，将步骤S4中内模式计算得到的第一年(月)末对桥墩周围河床形态的试验结果Δz三维(第一年或月末)替换步骤S3中外模式计算的桥墩周围河床冲淤厚度η(第一年或月末)，完成第一年末桥墩周围河床地形替换；
第二年开始计算时，外模式以替换后的第一年(月)末地形作为初始地形继续计算，当外模式计算至系列年第二年(月)末时，将内模式计算得到的第二年(月)末对桥墩周围河床形态的试验结果Δz三维(第二年或月末)替换外模式桥墩周围河床冲淤厚度η(第二年或月末)，完成第二年末桥墩周围河床地形替换；
为方便理解，内模式每个试验年(月)末桥墩周围河床变形Δz三维数据均作为外模式的输入条件的具体实现方法如下：当外模式计算至系列年第一年(月)末时，将内模式计算得到的第一年(月)末对桥墩周围河床形态的计算结果Δz三维(第一年或月末)替换外模式计算的桥墩周围河床冲淤厚度η(第一年或月末)，以此作为第二年计算的起始地形；当外模式计算至系列年第二年末时，将内模式计算的第二年末对桥墩周围河床形态的试验结果Δz三维(第二年或月末)替换外模式计算的桥墩周围河床冲淤厚度η(第二年或月末)，以此作为第三年计算的起始地形；以此类推，直至系列年末整个计算过程结束。以此类推，直至系列年末整个计算过程结束。最终得到系列年末，外模式计算得到整个计算区域内河床冲淤变形η(系列年末)，完成对于河段特征的模拟预测。
本发明具有如下技术效果：本发明采用二、三维数学模型耦合计算，利用内模式和外模式分别得到局部和整体河段河床冲淤变化，提高计算效率和精度，进而实现长河段高精度的复杂边界条件下河段特征的滩槽演变规律模拟预测。本发明可以实现兼顾局部河床冲淤调整和上下游之间联动规律的准确模拟。本发明可以用于内河长河段系统治理效果模拟，也可用于桥区河段“桥梁工程+航道工程+水沙调控”综合作用下滩槽演变规律预测，具有很强的实用性。
以上显示和描述了本发明的基本原理、主要特征和优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等效物界定。

Claims (7)

1.一种基于内外模式的二、三维数学模型耦合模拟方法，其特征在于，包括：
S1.获取河段特征，基于河段特征，构建整体河段水沙模型及局部河段水沙模型，其中，整体河段水沙模型为二维水沙数学模型，局部河段水沙模型为三维水沙数学模型；
S2.获取整体河段水沙模型的初始和边界条件；
S3.基于整体河段水沙模型的初始和边界条件，对整体河段水沙模型进行外模式计算，得到外模式计算结果，其中，外模式为对整体河段特征进行模拟计算；
S4.基于外模式计算结果，对局部河段水沙模型进行内模式计算，得到内模式计算结果，其中，内模式为对局部河段特征进行模拟计算；
S5.基于内模式计算结果，对整体河段水沙模型进行局部更新；
S6.重复步骤S3-步骤S5进行迭代计算直到满足迭代条件，得到整体河段特征的模拟预测结果；采用内模式计算的每个试验年末桥墩周围和河床变形Δz三维数据替换外模式计算得到的河床总的冲淤厚度η；
步骤S2中，获取所述初始和边界条件的具体步骤包括，
基于河段特征，获取整体河段水沙模型的初始条件；
其中，初始条件包括但不限于x方向流速、y方向流速、悬移质含沙量、底沙量；
基于整体河段水沙模型的初始条件，获取整体河段水沙模型的边界条件；
其中，整体河段水沙模型的边界条件包括上下游控制边界条件、固壁条件、收敛控制条件；
其中，上、下游控制边界条件：
上游进口条件：给定上游来流Qin、沙量Sin和进口推移质输沙率qb，进口各点流速进口各控制点流速由下式迭代算出：
式中：Uj，hj为进口计算网格点沿y方向流速和水深，dyj为离散网格间距，Vj＝0；
下游出口条件：给定水位ζout；
固壁条件：
流速采用非滑移边界条件，其边壁流速给定为零，即U＝V＝0；对于含沙量Si，底沙Ni计算中采用法向梯度为零条件：
收敛控制条件：
控制连续方程最大质量源bmax和通过各断面流量Qj，流速：水位：/>
步骤S4中，推移质运动对局部河床冲淤的影响计算：
其中，qb为推移质输沙率，ρsd为泥沙密度，ρw为水体密度，d50为泥沙中值粒径，fs为泥沙修正系数，D*为无量纲的泥沙粒径系数，τ′b为床面剪切力，τb,cr为床面泥沙起动临界剪切力，γ为系数；相关参数具体计算公式如下：
D*＝d50[(s-1)g/υ2]1/3
s＝ρs/ρw
τb,cr＝θcr(ρs-ρw)gd50
τ′b＝0.5ρwfcu2
fc＝8g[18log(12h/d90)]-2
其中，dsand＝0.062mm，dsand为泥沙中值粒径参数，υ为水体运动粘性，fc为水流产生的底部摩阻；
模型出口控制水位取长河段二维水沙模型计算各级流量时的水位ζ，模型进口控制流速U、V由长河段二维水沙模型计算得到，模型进口沙量S和推移质输沙率qb均由长河段二维水沙数学模型提供；
在通过推移质运动计算，得到河床变形之后，还包括对河床变形进行局部修正，具体步骤包括：
在通过推移质运动计算之后，对局部水沙模型进行网格扫描，基于扫描结果，对局部水沙模型进行局部调整，通过进行重复扫描及局部调整后，获取最终的修正结果，并基于修正结果，对河床变形进行修正；
在数值模型中，对于非均匀沙的推移质不平衡输移方程中推移质输沙率计算结果进行局部修正；具体操作方法如下：根据泥沙模型计算完河床变形后，对床面网格进行扫描；每个时间步网格节点高程根据地貌模型调整后，判断每个单元平面外法向方向与Z轴正向的夹角大小，当发现相邻角点间的倾角大于泥沙休止角时，要进行局部调整，在扫描时发现角点A(xA，yA，zA)和角点B(xB，yB，zB)之间倾角大于休止角需要将较高的B点降至B′点，A点升A′点，使AB之间的倾角降到/>ηA、ηB表示如下：
式中，ηA、ηB分别为外模式中A、B点的河床冲淤变形；
对整个床面进行反复的扫描检查直至床面上所有相邻角点间的倾角均在休止角范围内。
2.根据权利要求1所述基于内外模式的二、三维数学模型耦合模拟方法，其特征在于：
步骤S1中，构建整体河段水沙模型的具体过程包括，
基于河段特征，构建整体河段水沙模型的二维数学模型，并建立整体河段水沙模型的控制方程；
其中，整体河段水沙模型的控制方程包括水流连续方程，X方向水流运动方程、Y方向水流运动方程、悬沙不平衡输运方程、推移质不平衡输移方程和河床变形方程。
3.根据权利要求1所述基于内外模式的二、三维数学模型耦合模拟方法，其特征在于：
局部河段水沙模型为基于紊流模型构建的三维水沙数学模型。
4.根据权利要求1所述基于内外模式的二、三维数学模型耦合模拟方法，其特征在于：
步骤S3中，对整体河段水沙模型进行外模式计算的具体步骤包括：
基于初始和边界条件，对整体河段水沙模型进行离散求解，并基于离散求解结果通过压力校正法进行压力～流速耦合求解，得到外模式计算结果；
其中，离散求解结果包括x方向流速、y方向流速、水位、悬移质含沙量、底沙量、河床变形；
外模式计算结果包括流速场、水位值及含沙量场。
5.根据权利要求4所述基于内外模式的二、三维数学模型耦合模拟方法，其特征在于：
对整体河段水沙模型进行离散求解的过程中，离散求解的方法使用有限体积法、有限差分法及有限元法其中一种。
6.根据权利要求1所述基于内外模式的二、三维数学模型耦合模拟方法，其特征在于：
步骤S4中，通过局部河段水沙模型进行内模式计算的具体步骤包括：
基于外模型计算结果，得到局部河段水沙模型的边界条件，基于局部河段水沙模型的边界条件，通过空间离散及对流扩散进行离散求解，得到河床变形；
基于外模型计算结果，通过推移质运动计算，得到局部河段泥沙变化量；
其中，内模式计算结果包括局部河段泥沙变化量及河床变形。
7.根据权利要求1所述基于内外模式的二、三维数学模型耦合模拟方法，其特征在于：
步骤S6中，外模式计算结合内模式计算的迭代计算的具体过程包括，
获取系列水文年，对系列水文年进行划分，得到若干个试验年，对若干个试验年进行排序；
从第一个试验年开始，在每个试验年末，基于内模式计算结果，对整体河段水沙模型的初始条件进行更新，并通过重复步骤S3-S6进行迭代计算；
当对最后一个试验年末整体河段水沙的初始条件进行更新后，停止迭代计算。
CN202110785304.3A 2021-07-12 2021-07-12 一种基于内外模式的二、三维数学模型耦合模拟方法 Active CN113486557B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110785304.3A CN113486557B (zh) 2021-07-12 2021-07-12 一种基于内外模式的二、三维数学模型耦合模拟方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110785304.3A CN113486557B (zh) 2021-07-12 2021-07-12 一种基于内外模式的二、三维数学模型耦合模拟方法

Publications (2)

Publication Number Publication Date
CN113486557A CN113486557A (zh) 2021-10-08
CN113486557B true CN113486557B (zh) 2023-07-21

Family

ID=77938720

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110785304.3A Active CN113486557B (zh) 2021-07-12 2021-07-12 一种基于内外模式的二、三维数学模型耦合模拟方法

Country Status (1)

Country Link
CN (1) CN113486557B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116911208B (zh) * 2023-06-27 2024-05-24 中山大学 一种海洋桩基局部冲刷坑预测方法、系统、设备及介质

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105894106A (zh) * 2015-01-05 2016-08-24 国家电网公司 一种海洋模式和气象模式的一体化耦合方法
CN109994364A (zh) * 2018-01-03 2019-07-09 原子能和替代能源委员会 微机械或纳米机械颗粒检测设备

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5537641A (en) * 1993-11-24 1996-07-16 University Of Central Florida 3D realtime fluid animation by Navier-Stokes equations
CN107480384A (zh) * 2017-08-24 2017-12-15 北方民族大学 河道水流泥沙平面二维数值模拟方法和系统
CN108532532B (zh) * 2018-04-12 2020-10-16 福建省水利水电勘测设计研究院 感潮河道防洪防潮岸线制定方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105894106A (zh) * 2015-01-05 2016-08-24 国家电网公司 一种海洋模式和气象模式的一体化耦合方法
CN109994364A (zh) * 2018-01-03 2019-07-09 原子能和替代能源委员会 微机械或纳米机械颗粒检测设备

Also Published As

Publication number Publication date
CN113486557A (zh) 2021-10-08

Similar Documents

Publication Publication Date Title
Qi et al. Wave induced instantaneously-liquefied soil depth in a non-cohesive seabed
CN108629055B (zh) 一种基于饱和输沙原理的沙质内河航道回淤量预报方法
CN107480384A (zh) 河道水流泥沙平面二维数值模拟方法和系统
CN111428401A (zh) 一种堰塞湖溃决过程模拟方法
CN113486557B (zh) 一种基于内外模式的二、三维数学模型耦合模拟方法
CN110135069B (zh) 输水隧洞输水时的泥沙特征获取方法、装置、计算机设备
CN113536643B (zh) 一种基于数字孪生的长河段滩槽演变预测方法及系统
CN109271672A (zh) 一种河-湖-泵站相互影响作用下的河道水面线计算方法
CN106959261A (zh) 一种预测沉积物颗粒分布与配比的方法
Duc et al. Numerical simulation of contraction scour in an open laboratory channel
CN115034468B (zh) 一种淤泥质海岸丁坝或防波堤拆除后海床冲淤变化预测方法
Sam et al. River sedimentation and morphology modeling–The state of the art and future development
Diab Experimental investigation on scouring around piers of different shape and alignment in gravel
Badano et al. A coupled hydrosedimentologic model to assess the advance of the Parana River Delta Front
Arafat et al. Hydrodynamics and morphological changes numerical model of the Jeneberang estuary
Paquier et al. River bed deformation calculated from boundary shear stress
CN112485106B (zh) 一种控制土体状态参数的物理模型分层制备与试验方法
Xie et al. Mathematical model for flood routing in Jingjiang River and Dongting Lake network
CN114357903A (zh) 基于河床形态参数的冲积河流阻力计算、模拟方法及装置
Khorrami et al. Development of a non-coupled algorithm for simulating long-term sedimentation in the Zonouz dam reservoir, Iran
Isaac et al. Sediment management of run-of-river hydroelectric power project in the Himalayan region using hydraulic model studies
Level Calibration of CCHE2D for sediment simulation of Tarbela Reservoir
Sobhkhiz et al. Numerical modeling of the effect of pile geometry and pile cap leveling on the local scour under inclined pier group
CN113807035B (zh) 一种平面涡旋水流作用下河床冲淤的计算方法
CN116933960B (zh) 一种沙坝潟湖型潮汐汊道航槽选线方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant