矩形与曲形坝体溃坝瞬时的水位研究

任超洋，陈善群，廖 斌

(安徽工程大学， 安徽 芜湖 241000)



矩形与曲形坝体溃坝瞬时的水位研究

任超洋，陈善群，廖 斌

(安徽工程大学， 安徽 芜湖 241000)

溃坝是一种自然灾害，且时常发生，溃坝时大量的水从缺口溢出，会对周围造成很大的影响。坝体的形状是引起溃坝的主要影响因素之一，当来流作用在坝体上之后，由于坝体形状的不同，其溃坝之后所导致的流体冲击距离与力度也不一样。为此，利用VOF模型，用界面跟踪法求解相连续方程，对矩形与曲形坝体溃坝问题进行数值模拟，并对两种形状坝体溃坝后的水位图进行分析，通过数值计算得到了各个时刻溃坝运动界面的变化和相应的速度分布。最后得出结论，矩形坝体相对于曲线形坝体较好。

形状；溃坝；VOF模型

近年来，由于自然灾害的影响，导致全国范围内北方出现大规模旱灾，而南方则暴雨连连，多处受灾且时有溃坝问题传出。最近几年内世界上有很多学者在溃坝模拟方面做了相当多的计算工作。其中，刘艾明等[1]通过利用VOF法(流体体积法)模拟二维溃坝中溃坝高度不同的问题。Martin J C等[2]对直立方水柱在重力作用下塌陷过程进行了试验并得出了结果。Stans P K[3]等对前后具有不同高程比的两块水体溃坝过程进行了试验研究。Stoker[4-5]在坝体瞬间崩溃的理想条件下得出了溃坝时在下游有水条件下的解析解。陈善群等[6]通过研究基于FIC法的不可压缩N-S方程稳定分步算法，改善了经典算法的压力稳定性。但是以上学者都没有将矩形与曲形这两种状态下的水体溃坝放在一起进行比较，也没有得出一个合适的结论。

本文基于Fluent平台，通过采用VOF模型，对剪切流越坝情况进行分析，并且对由于剪切力造成的溃坝情况进行了数值模拟。当溃坝时会有大量的水从缺口冒出，在这种情况下利用VOF模型可以非常方便地模拟仿真溃坝过程。最后将结果进行对比，并且绘制两种不同形式的坝体溃坝之后的水位图。

1 VOF模型

VOF模型是一种在固定的欧拉网格下的表面跟踪方法，通过求解单独的动量方程和处理穿过区域的每一流体的体积分数来模拟2种或3种不能混合的流体。当需要得到一种或多种互不相容流体间的交界面时，可以采用这种模型。

2 总体计算流程与三大控制方程

计算流体力学(CFD)的求解流程可以概括成以下几个方面：建立控制方程→确立初始条件以及边界条件→划分计算网格以及生成计算节点→建立离散方程→离散初始条件和边界条件→给定求解控制参数→求解离散方程，最后看结果是否收敛。若收敛，则显示和输出计算结果；若不收敛，则继续回到建立离散方程那一步。

控制方程就是由于流体流动所受到的物理守恒定律的支配，其基本的守恒定律包括以下3个：质量守恒定律、动量守恒定律(牛顿第二定律：力=质量×加速度)、能量守恒定律；所以，控制方程是3种守恒定律的数学描述。

2.1 质量守恒方程

在流体力学中质量守恒定律可以描述为单位时间内流体微元体中质量的增加，等同于同一时间间隔内流入该微元体的净质量。所以，按照该定律，可以用以下数学方程表达：

(1)

2.2 动量守恒方程(Navier-Stokes)

在流体力学中动量守恒定律可以表述为微元体中流体的动量对时间的变化率等于外界作用在该微元体上的各种力之和。按照这一定律，可以导出x、y和z3个方向的动量守恒方程为：

其中，Su、Sv、Sw是动量守恒方程的广义源项，Su=Fx+sx，Sv=Fy+sy，Sw=Fz+Sz。

2.3 能量守恒方程

能量守恒定律是包含热交换的流动系统必须满足的基本定律。该定律可以表述为微元体中能量的增加率等于进入微元体的净热量加上体力与面力对微元体所做的功。

其数学表达式为：

(5)

其中，Cp是比热容，T为温度，k为流体的传热系数，ST为流体的内热源及由于黏性作用流体机械能转化为热能的部分。

3 两种不同坝体的模型介绍以及计算分析

3.1 问题介绍以及模型建立

此问题采用2D平面计算模型，取一个长8 m，宽6 m的长方形区域，此时只受重力作用，沿y轴-9.8 m/s2，由于次模拟是模拟不同形状大坝溃坝瞬时的水位情况，所以取时间长度为1 s，即研究1 s 时长的水位变化。具体模型见下图：

图1 矩形水体

图1所示的水体是一个长4 m，宽2 m的矩形区域；图2所示的水体是一个半径为3 m的半圆形区域，通过简化的半圆形水体来代替曲线形坝体。

图2 曲形水体

对于划分网格这一环节，是CFD(计算流体力学中)必不可少的一项重要步骤，划分网格的好坏直接影响数值模拟结果的精度以及计算时间的长度，一个好的网格往往意味着节省时间和较高的计算精度。因为在使用数值模拟进行计算控制方程时，都是想办法将控制方程在空间区域上进行离散，然后求解得到的离散方程组，所以必须用到网格。在一片网格区域内，有许多的网格节点，在数值模拟方法中，这些网格小结点其实就是需要求解的未知物理量的几何位置，在离散过程中，将一个控制体积上的物理量定义且存储在该节点处。

在进行数值模拟计算时，常常会遇到如下2个选择性问题，即如何划分网格，是选择划分结构型网格还是选择划分非结构型网格。结构型网格与非结构型网格在计算时会有很大的不同，其中结构型网格往往在计算时间以及计算结果的精确度上要更胜一筹。所以，该如何将非结构型网格变成结构型网格，就需要将计算区域分快处理。先将一块较大的计算区域划分为几个较小的规则计算区域，这是因为规则的计算区域会使网格也变得规则，即所谓的结构型网格。在需要计算的地方将网格加密，而不需要计算的区域或是通过该数值模拟在某些地方不要求过多计算的区域则不需要过密的网格，在这些区域只要把网格画的疏一些即可。

本次的数值模拟不同坝体形状的溃坝瞬时水位的网格是一个很简单的2D平面规则网格。如图3所示。

图3 平面网格

3.2 矩形和曲线形坝体的溃坝数值模拟后处理与结果分析

3.2.1 数值模型建立

利用VOF模型可以建立起图4所示的气液相分离的数值模拟图。图4中左边矩形部分表示液相(water)其余部分表示气相(air)，当大坝的时候相当于左边的红色液相部分瞬间塌陷下来。当溃坝发生时，会有一个速度与水位的相对关系。而图5中则将左边的半圆形部分设置成液相(water)，将其余部分设置成气相(air)，为了区别对待，将2种不同方式的溃坝显色调成了相反的颜色，以利于比较。

图4 矩形VOF

图5 曲形VOF

3.2.2 后处理及分析

本案例研究大坝瞬时溃坝的情况，通过对比分析溃坝瞬间0.2 s、0.5 s与0.9 s这3个时刻的气液相云图以及这3个时刻的速度云图得出2种不同形状坝体溃坝的流向特征，并在随后通过对比两种不同形状坝体溃坝之后的水位图得出一个较好的坝体形状。

1) 0.2 s时的云图分析

通过图6与图7分析可知：在0.2 s时，矩形水体的底部发生了明显的变化，开始往前面溢出，而速度方面则有两处发生了较明显的变化，这两处分别在矩形水体的左侧上部与底部；而曲形水体在0.2 s时则在上端与下端分别开始有水流向下溢出的现象，而速度变化则是在水流溢出处的圆切线部分有较大变化。

图6 矩形0.2 s时的相图与速度云图

图7 曲形0.2 s时的相图与速度云图

2) 0.5 s时的云图分析

通过对图8和图9对比分析可得：在0.5 s时，两种水体的变化加大。且通过速度云图的比较可以发现：在该时刻的两种形式的水体都是在中部出现较大的速度，在速度这个点上十分的相似。

图8 矩形0.5 s时的相图与速度云图

图9 曲形0.5 s时的相图与速度云图

3) 0.9 s时的云图分析

通过图10和图11的对比可以清楚地发现：曲形水体已经达到右壁面，而矩形水体还差一点才达到右边壁面。就速度云图的比较可以发现：曲形水体在靠近右边壁面的速度更大。

4) 矩形坝体与曲形坝体溃坝时在0.2、0.5、0.9 s 三个时刻的水位图对比

如图12～14所示，通过对矩形与曲线形大坝的水位图对比，虽然在开始0.2 s时的水位变化并不明显，但可以清楚地发现在0.5 s时的矩形水体呈现一个滑梯式扩散形式。这种滑梯式的扩散方式有一个特征，即从后往前高一点的水体由于重力的影响，将重力势能传递给前面的水体，也就是说这种滑梯式的扩散方式能量传递十分均匀，且水流在0.5 s的时候才达到距离原点4.25 m的距离；而曲线形水体在0.5 s的时候水流在水平线上已达4.9 m左右，且由于曲线形水体的能量传递不均匀，最重要的问题是在0.5 s时由于中部有一段水体是凸起的，所以导致后续的能量传递相当巨大，可以从0.9 s时的水位图发现一个问题，即矩形水体在水平方向上冲击的距离有7.4 m左右，而曲线形水体则已经超过了8 m。

图10 矩形0.9 s时的相图与速度云图

图11 曲形0.9 s时的相图与速度云图

图12 0.2 s时两种形状溃坝水位对比

图13 0.5 s时两种形状溃坝水位对比

图14 0.9 s时两种形状溃坝水位对比

4 结束语

本文对两种形式坝体的溃坝数值进行模拟，可以明确得出一个结论，即从坝体溃坝之后水流冲击的距离来看，矩形坝体好于曲线形坝体。这对于溃坝这一自然灾害有了一个启示，在溃坝救灾时，往往是分秒必争，矩形坝体在溃坝之后的水流来流较缓和，而曲线形坝体溃坝之后的水流则相对激烈。本文为今后研究溃坝问题给出了数值依据。

[1] 刘艾明,刘学炎,熊鳌魁.二维不同溃坝方式的数值模拟[C]//第九届全国水动力学学术会议暨第二十二届全国水动力学研讨会论文集.2009.

[2] Martin J C,Moyce W J.Part IV.An experimental study of the collapse of liquid columns on a rigid horizontal plane[J].Philosophical Transactions of the Royal Society of London A:Mathematical,Physical and Engineering Sciences,1952,244(882):312-324.

[3] STANS P K,CHEGINI A,BARNES T C.The initial stages of dam-break flow[J].Journal of Fluid Mechanics,1998,374:407-424.

[4] STOKER J J.Water wave[M].New York:Interscience Publishers,1957.

[5] MORRIS C.Concerted action on dam-break modeling,final report[R].Oxfordshire:HR Wallingford Ltd,2000.

[6] 陈善群,廖斌.基于 FIC 法的不可压缩 NS 方程稳定分步算法[J].水动力学研究与进展(A辑),2015,30(2):129-139.

(责任编辑 陈 艳)

Research on Dam-Break Water Level for Rrectangular and Curved Dam

REN Chao-yang, CHEN Shan-qun, LIAO Bin

(Anhui Polytechnic University, Wuhu 241000, China)

Dam break is a kind of frequent natural disaster, and there is a large amount of water overflows from the gap when the dam breaks, and it will have a great impact on the surrounding. The shape is one of the main influencing factors of dam break, when the incoming flow acts on the dam body, due to the different shape of the dam, the damper caused by the impact of the distance and the intensity of the fluid is not the same. This article used the VOF model in solving continuous equations by interface tracking method, numerical simulation of dam break problem for rectangular dam and curved dam. The water level maps of two types of dams were analyzed, and simulation results focused on the evolution of the moving interface and the contour of the velocity of the dam break at different moment. Finally, it is concluded that the rectangular dam is better than the curved dam.

shape; dam break; VOF model

2016-11-25

安徽省自然科学基金资助项目(1508085QE100)

任超洋(1991—)，男，安徽芜湖人，硕士研究生，主要从事计算流体力学研究,E-mail:2906510498@qq.com；通讯作者 陈善群(1981—)，安徽合肥人，博士，主要从事计算风工程研究，E-mail:chenshanqun@126.com。

任超洋，陈善群，廖斌.矩形与曲形坝体溃坝瞬时的水位研究[J].重庆理工大学学报(自然科学)，2017(4):76-81.

format：REN Chao-yang, CHEN Shan-qun, LIAO Bin.Research on Dam-Break Water Level for Rrectangular and Curved Dam[J].Journal of Chongqing University of Technology(Natural Science)，2017(4):76-81.

10.3969/j.issn.1674-8425(z).2017.04.012

TV8

A

1674-8425(2017)04-0076-06