溃口近区二维数值模拟与溃坝洪水演进耦合

杨 志，冯民权

（西安理工大学陕西省西北旱区生态水利工程重点实验室，陕西 西安 710048）

溃口近区二维数值模拟与溃坝洪水演进耦合

杨 志，冯民权

（西安理工大学陕西省西北旱区生态水利工程重点实验室，陕西 西安 710048）

基于黏土心墙砂石坝的溃决过程，以及溃坝洪水传播和运动的特性，建立黑河金盆水库大坝溃口近区二维数值模型和下游地区溃坝洪水演进耦合数学模型。使用DAMBRK法计算逐渐溃坝，并应用其结果进行后续模拟。采用Abbott⁃Ionescu六点隐式有限差分格式求解一维模型，采用单元中心的有限体积法求解二维模型方程。采用侧向连接方式，将黑河两岸计算水位点与二维网格单元相连，实现一、二维模型的耦合。采用所建立的二维模型对溃口近区进行计算与模拟，得到计算区域某一时刻的水深及流速分布。应用所建耦合模型对黑河金盆水库万年一遇入库洪水漫顶致溃坝洪水进行数值模拟，得到一维河道内各断面的水位和流量变化过程，以及二维计算区域内不同时刻的水深分布图、流速矢量图和淹没范围变化过程。溃口的形成过程不仅包括漫顶水流的直接作用，同时包括溃口形成过程中两侧漩涡状水流的反冲刷作用。耦合模型可以同时兼顾河道内的水流变化以及河道外计算区域内的洪水演进过程，从而减少由于计算结果偏大或偏小所带来的防洪资源浪费和防洪措施不利等不良影响。

溃口近区；二维模拟；溃坝洪水演进；侧向连接；耦合模拟

杨志，冯民权.溃口近区二维数值模拟与溃坝洪水演进耦合［J］.水利水运工程学报，2015（1）：8－19.（YANG Zhi，FENG Ming⁃quan.2D numerical simulation of breach area and coupling simulation of dam⁃breach flood［J］.Hydro⁃Science and Engineering，2015（1）：8－19.）

溃口发展和溃口流量的计算和模拟是研究溃坝问题的基础和重点，是计算溃坝洪水演进、分析溃坝风险以及制定地区应急预案的前提。在物理试验方面，欧盟委员会通过现场和实验室试验进行了大量的溃坝研究工作［1］。新西兰、美国等国家的学者对土坝溃决机理进行了相关的试验研究［2－3］。G.Pickert等［4］基于详细的试验数据，研究均质非黏性大坝漫顶溃决溃口发展的两个阶段。段文刚等［5］采用“埋入式轻型冲蚀捕捉器”来动态记录溃决过程，在室内试验的基础上，分析了无黏性非连续宽级配砂砾料土坝溃决机理。在溃口计算方面，V.Nourani等［6］使用不同材料建造室内大坝模型，采用人工神经网络技术和遗传算法观测土石坝溃口流量峰值和溃口最终深度和宽度。方崇惠等［7］建立了瞬时溃坝最大流量与堰流关系，理论推导得到可以计算大坝全溃、横向局部溃、垂向局部溃及横垂向局部溃的瞬时最大流量新通式。在数值模拟方面，B. Dewals等［8］将一个二维的全动态模型和简化的集总模型相结合，模拟计算多级复杂水库逐渐溃坝。张社荣等［9］建立了一维数值模型，模拟了土心墙堆石坝漫顶逐渐溃过程，以及混凝土重力坝漫顶不溃与漫顶瞬时溃过程。袁晶等［10］利用动网格技术，建立了可变网格下的堤防溃口展宽二维数学模型。梁艳洁等［11］采用平面二维水流数学模型，模拟漫顶溃堤洪水演进过程，结果表明堤防溃决属溯源冲刷，与流速密切相关。

为分析溃坝洪水风险、评估灾害损失以及制定地区应急预案，研究溃坝洪水演进十分必要。由于溃坝洪水漫出河道，二维数值模拟是研究主流，然而粒子法［12］、SPH法［13］等无网格技术尚不成熟，因此本文应用非结构网格法进行二维溃坝洪水演进模拟。为了更好地发挥一维、二维计算模型的优势，并提高计算精度，国内外学者［14－16］分别耦合一维圣维南方程和不同的二维模型进行水动力计算和模拟。H.A.Gallegos等［17］将结构损伤模型与水动力模型双向耦合，预测坝体溃决和局部洪水及其相互影响。Chao Dang等［18］建立自然溃坝过程和溃坝洪水过程预测的耦合数学模型。耦合计算的关键在于一维、二维模型耦合点交互计算的流畅性，为保证这一点，刘薇等［19－21］分别采用“流量守恒条件”、“重叠计算区域”和“堰流公式”对耦合点进行处理，完成连接处的水利要素交互计算，从而实现一、二维模型的耦合。基于MIKE FLOOD模型，孙秀丽等［22］采用标准连接的方式，建立一种溃坝洪水模拟的一维、二维耦合模型；李传奇等［23］通过侧向连接方式，建立一维二维水动力耦合模型模拟复杂市区降雨积水过程。在水沙两项耦合方面，夏军强等［24］建立了动床条件下溃坝水流的二维水沙耦合模型。此外，F.Nzokou等［25－26］将水动力模型和覆冰模型耦合求解，模拟冰层覆盖工况下的溃坝水流模拟。

目前，尽管溃坝数值模型研究已取得了长足的进展，但仍存在有待深入探讨的问题。例如，缺乏将入库洪水、溃口模拟及下游河道洪水波传播模拟作为一个整体进行分析计算，无法分析三者之间的作用对溃坝水力特性的影响；另外，一、二维嵌套计算及动态反映各计算参数间的相互作用是溃坝洪水模拟中的重难点。针对以上问题，在参考前人研究成果的基础上，建立非结构网格二维模型，模拟黑河金盆水库大坝溃口近区淹没水深和流场变化；采用侧向连接方式建立一、二维耦合模型，模拟万年一遇入库洪水漫顶致溃坝洪水在下游地区的演进过程。旨在为黑河金盆水库的安全运行和调洪决策提供参考，并为洪水风险分析和灾害损失评估提供依据。

1 数学模型

1.1 一维模型原理

一维计算模型是基于垂向积分的物质（水量）和动量守恒方程，即应用一维非恒定流Saint⁃Venant方程组来模拟河流或河口的水流状态。

式中：A为过水断面面积（m2）；t为计算点时间坐标；Q为过流流量（m3／s）；x为计算点空间坐标；q为旁侧入流流量（m3／s）；α为动量校正系数；h为水位（m）；g为重力加速度（m／s2）；C为谢才系数；R为水力半径。

方程组采用Abbott⁃Ionescu六点隐式有限差分格式进行求解。该离散格式在每个网格点并不是同时计算流量和水位，而是按照顺序交替计算流量和水位，分别称作Q点和h点。Abbott⁃Ionescu格式具有计算精度高、稳定性好的特点。离散后的线性方程组用追赶法求解。

1.2 二维模型原理

由于溃坝洪水运动水平尺度远大于垂直尺度，水深、流速等水力参数沿垂直方向的变化较水平方向的变化小得多，因而可将三维流动控制方程沿水深积分，并取平均水深，从而得到二维浅水水流运动方程。

水流连续方程：

水流运动方程：

式中：t为时间（s）；x为横向空间坐标（m）；y为纵向空间坐标（m）；Z为x，y处的水位（m）；h为x，y处的水深（m）；u为x方向的流速分量（m／s）；v为y方向的流速分量（m／s）；g为重力加速度（m／s2）；fv，fu为科氏加速度项；为表面水位加速度项；为大气压力梯度项；为波浪辐射应力项；Fu，Fv为水平涡黏项；usS，vsS是源项入流产生的加速度项。

针对溃坝洪水数值计算面临不规则边界和复杂地形等问题，建立非结构三角形网格下的二维水流模型，将任意三角形作为计算单元，采用单元中心的有限体积法求解方程，单元网格的分布如图1所示。其中，点M代表各项结果变量位于单元中心，F为跨边界通量，分别垂直于所在三角形单元中相应的各边。

图1 非结构三角形网格及其单元的分布Fig.1 Unstructured triangular gridand its units distribution

对浅水水流运动方程进行归一化，则有：

式中：U为守恒变量；FIx－FVx为x方向的通量；FIy－FVy为y方向的通量；S为源项。

式中：Ui为控制单元i的守恒变量；Ai为控制单元i的面积；N为控制单元边的个数；F为通过第j条边的通量；n0为通过第j条边的法向量；ΔΓj为控制单元各边长度；Si为控制单元i的源项。

用Roe's近似黎曼解法求解跨边界通量F，使用了二阶TVD格式的限制器函数避免数值振荡。

1.3 耦合模型方式

一维河道与二维网格通过统一的连接高程实现计算单元的重叠。一维连接高程通过一维河道堤岸高程确定，即将断面的左、右堤岸高程作为耦合连接的河床高程，无断面数据处则通过插值得出；二维连接高程是通过二维的地形确定，即将二维网格单元的高程作为耦合连接的河床高程。重叠单元通过设置水力计算节点进行搭接，一维计算节点为河道左、右堤岸连线上的计算水位点，二维计算节点在一维河道断面的左、右端点坐标对应的二维网格单元中。耦合模型对每一个节点都进行水流动态计算，并根据节点处水量平衡，将每一个经过计算的水流重新分配到二维网格单元和一维计算水位点中，从而实现一维与二维的耦合计算。一维、二维模型均采用统一的时间步长，并在每一个时间步长时刻进行信息的相互交换，当水流从一维河道向二维网格流动时，由一维模型求解出计算节点处的物理量，以此作为二维模型的边界条件；反之则以二维模型在节点处的出流值，作为连接处一维模型的边界条件。为了使耦合模型计算稳定，需要确保连接点处的水力条件保持动态平衡，即要求一维模型和二维模型中的初始水位、流量条件应当相同。将二维河道中所有未知变量用河道首末断面水位线性表达，再将首末边界处的流量表达成边界水位形式。这里需要说明：耦合模型只是搭接一维模型与二维模型交互计算的一个平台，耦合之后的结果仍是一维河道内与二维淹没区域上分别进行查看与显示。

1.4 二维模型验证

1.4.1 试验概况及模拟条件 黑河河道上有两处大于90°的转弯，故本文选取L形河道室内二维溃坝试验验证所建模型对弯曲河道的适用性［28］。该试验中，计算域由一个2.40m×2.40m的正方形水库连接L形状的90°弯曲河道组成，河道宽0.5m，里程7.25m（见图2）。河道上布置有6个观测点，其编号名称及相应位

对方程（4）在每一个单元上进行积分，经进一步简化后，得：置坐标分别为：观测点P1（1.20 m，1.20 m），观测点P2（2.75m，0.70m），观测点P3（4.25m，0.70m），观测点P4（5.75m，0.70m），观测点P5（6.55m，1.50m），观测点P6（6.55 m，3.00 m）。水库初始水深为0.2 m，流速为0；下游为干底河床。水库边界及河道左、右岸为固壁边界条件，河道出口处为自由出流边界条件。水库及河道均为平底，曼宁系数取n＝0.009 5 s／m1／3。

计算域采用非结构网格进行剖分，共计2 936个三角形单元和1 600个节点。假设在t＝0时刻，大坝瞬时溃决，模拟40 s内溃坝水流的运动情况。

1.4.2 模型验证合理性分析 6个观测点实测值和水深

数值解［28］见图3。从图3可以看出，观测点的数值解和实测值比较吻合，计算误差为8.51%；模型可以准确估测下游各观测点洪水到达时间，并且可以得到观测点P2，P3和P4都存在两次水位迅速上涨的情况。下游观测点P2，P3和P4所在河段位于大坝和河道上90°弯角之间（见图2），第1次水位的快速上涨是大坝溃决后，水体发生急剧下泄引起的；第2次水位的快速上涨是因为河道弯角处的阻水作用。此外，因为下游观测点P5和P6并没有直接受到河道弯角处的阻水作用影响，所以水位迅速上涨一次，此后随着水流流出计算区域，观测点水位逐渐下降。

图2 计算区域和测点位置（单位：m）Fig.2 Calculation region and measuring points（unit：m）

模拟结果还直观反映出河道直角转弯处的阻水现象。

图3 测点水深数值解和实测值的对比Fig.3 Comparison between numerical solutions and measured values of themeasuring points depth

2 逐渐溃坝计算

2.1 调洪演算

水库调度方式如下：

（1）当库水位低于汛限水位591 m时，下泄流量为0。

（2）当库水位为汛限水位至防洪高水位（591～594 m）时，若来流量小于下游允许泄流量1 800 m3／s，按总下泄流量等于来流量调洪；如来流量超过下游允许泄流量，控制总下泄流量等于1 800 m3／s。

（3）当库水位高于防洪高水位594 m时，闸门全部开启。

依选定溃坝计算条件，当万年一遇洪水经入库调洪，静水位等于600 m时，大坝溃决。黑河金盆水库水位库容关系曲线如图4，此时库容为21 680万m3。根据黑河金盆水库万年一遇入库洪水的流量过程线（见图5），最大入库流量7 800 m3／s。经过调洪演算得到水库下泄流量过程线见图6。从图6可见，t＝0时刻开始调洪，经过16.4 h，水库水位达600 m，大坝溃决。

图4 金盆水库水位库容关系Fig.4 Water level vs storage capacity

图5 万年一遇入库洪水的流量过程Fig.5 Reservoir inflow hydrograph

图6 洪水调洪演算结果Fig.6 Flood regulating calculation

2.2 最终溃口宽度

溃坝时金盆水库蓄水量V≥100万m3，故溃口平均宽度b计算式如下：

式中：b为溃口平均宽度（m）；K为坝体材质系数，对土心墙坝取1.19；V为水库蓄水量（万m3）；B为坝址处水面宽，通常等于坝长（m）；H为溃坝前上游水深（m）；B／b一般不应超过17。该式计算结果为 b＝340.26 m。

土石坝逐渐溃决的溃口形状概化为梯形，溃口形状的确定主要取决于最终溃口底宽bm和溃口边坡m。最终溃口底宽bm＝b－mhd，其中：hd为溃坝有效水深，一般以溃口相对于坝顶的高程近似代替。m的大小与坝体材料有关，取值范围一般为0≤m≤2。b＝0，m＞0时溃口概化形状为三角形；b＞0，m＝0时为矩形；b＞0，m＞0时为梯形。溃口平均宽度b＝340.26 m时，最终溃口底宽bm＝190.53 m。

2.3 溃决历时

2.4 瞬时溃口形态和流量

溃决过程中任意t时刻的溃口瞬时底部高程ht和瞬时底部宽度bt计算公式如下

式中：h为大坝残留高度（m）；n表征溃口非线性程度，一般情况下1≤n≤4，本文取n＝1，即溃口呈线性变化；0≤t／Tf≤1。

瞬时溃口流量Qt计算式为：

式中：Ht为瞬时坝前水深（可经试算获得）（m）；θt为边坡与横向的瞬时夹角。本文的计算参数取值见表1。

表1 模型相关参数Tab.1 Parameters related to DAMBRK model

2.5 计算结果

以大坝溃决初始时刻（t＝16.4 h）作为计算开始时间，计算2种溃决历时Tf＝0.906 1 h和Tf＝0.897 5 h下的溃口下泄流量变化过程及溃口扩展变化过程，计算结果见图7。

图7 2种溃决历时条件下的模拟结果Fig.7 Simulation results of different dam breach durations

由图7可见，当Tf＝0.906 1 h时，溃口下泄流量峰值达到80 504.89 m3／s；溃口底部高程和溃口底宽于0.906 1 h分别达到528.7 m和190.534 m。当Tf＝0.897 5 h时，溃口下泄流量峰值达到84 334.09 m3／s；溃口底部高程和溃口底宽于0.897 5 h分别达到528.7 m和190.534 m。两种溃决历时计算结果影响下的洪峰流量和溃口扩展过程差别不大，但考虑到溃坝计算的特殊性，本文选取Tf＝0.897 5 h的计算结果作为溃口进一步模拟的初始条件。DAMBRK模型不仅可以给出溃口发展稳定后的最大值，并且可以给出溃口形态及流量的变化过程，计算结果更加全面，也更加具体。

3 溃口近区二维模拟

3.1 模拟条件

溃口近区在溃坝洪水演进区域中相当于一个很小的点，以溃口近区作为研究对象，旨在深入具体研究溃口及近区的流量、水深和流速变化过程，从而补充和完善洪水演进过程的模拟。二维模拟中，拟定了一个长1 000 m，宽433 m的矩形计算区域。将其进行网格剖分，每个网格平均面积约为68 m2，共计6 410个三角形单元和3 350个节点。根据调洪演算结果（图6）和DAMBRK模型计算中Tf＝0.895 7 h的结果（图7（c）），整理得到万年一遇洪水下泄流量过程线（见图8）作为溃口近区二维模拟的边界条件。模拟时段为大坝溃决初始时刻至溃口口门达到最大时刻，共计54 min。

3.2 模拟结果

通过二维模拟，给出大坝溃决（t＝0）后溃口扩展过程中t＝10，24，48和54 min时刻的瞬时水深及流场分布（见图9）。图9反映了溃口扩展与水流变化的相互作用规律，即当水流漫顶，溃坝开始；水流流过溃口后，会在溃口两侧形成漩涡对其冲刷（图9（a））；冲刷由轻到重（图9（b）），最终趋缓（图9（c））；在54 min时，由于溃口的限制，水流流速仍呈扩散状从溃口流出，但是已经几乎看不到涡漩回至溃口的水流继续冲刷溃口两侧，此时认为形成最终溃口（图9（d））。因此得到溃口的形成过程不仅包括漫顶水流的直接作用，同时包括溃口形成过程中两侧漩涡状水流的反冲刷作用。

图8 万年一遇洪水下泄流量过程线Fig.8 Releasing discharge of 10 000 years flood

图9 不同时刻水深及流速分布Fig.9 Water depth and velocity distribution at different times

4 溃坝洪水演进耦合模拟

4.1 一维模拟条件

在一维洪水演进模拟中，河道计算部分长30 920 m，溃口近区可近似看作1个断面。根据现有资料，在黑河流域上另选取6个断面，所选断面分别与坝址断面依次相距2 700，7 300，1 3700，19 100，25 000和30 920 m（见图10）。采用1万年一遇入库洪水下泄流量过程（见图8）作为一维计算入流的边界条件。在黑河入渭口设置出流条件，其为根据断面011＃的高程（图11）计算所得流量水位关系Q＝Q（z）（见图12）。

图10 断面位置Fig.10 Section positions

图11 011＃断面Fig.11 011＃section

图12 黑河入渭口流量水位关系Fig.12 Discharge⁃stage relationship curves

4.2 二维模拟条件

二维计算区域涉及黑河金盆水库坝址至黑河入渭口及渭河部分流域，总面积为1 066 km2，地形及采样点位置见图13。区域内采用三角形网格划分，并通过网格密度控制优化处理突变边界，每个网格面积约为0.065 8 km2，共有16 199个网格。由于黑河金盆水库地形复杂，糙率值难以准确估计，而糙率选取的恰当与否对计算成果有较大的影响，糙率取值范围一般为0.025～0.045，参照以往研究［27］，最终拟定二维计算选定河床糙率值为0.031。采用1万年一遇入库洪水下泄流量过程（见图8）作为二维计算入流的边界条件，其y方向流速假设为零，出流边界采用开边界条件。耦合模型采用侧向连接方式，一维计算部分包括整个黑河河道，二维计算部分则为整个二维网格区域（见图14）。二者通过水力计算节点进行搭接，溃坝水流以侧向堰流的形式从一维黑河河道流向二维计算区域。

图13 黑河金盆水库计算区域地形Fig.13 Calculation area topography of Jinpen reservoir located at Heihe River

图14 一维、二维耦合示意Fig.14 Schematic diagram of one⁃dimensional and two⁃dimensional coupling

4.3 耦合模拟结果

通过查看一、二维耦合模拟结果，得到4个典型断面的水位变化过程（图15）和5个观测断面的流量变化过程（图16）。由图15可知溃坝洪水波到达3＃，4＃，5＃和011＃断面的时间分别为模拟开始后的16.79，17.07，17.23和17.40 h；各断面水位依次于17.90，18.23，18.87和18.92 h时刻达到最大，峰值分别为426.67，419.34，415.71和405.22 m。图16中，溃坝洪水波到达5个观测断面的时间分别为16.57，16.73，16.95，17.07和17.23 h；洪水流量依次于17.73，17.89，17.94，18.07和18.92 h达到峰值，洪峰流量分别为78 480.02，43 375.16，24 654.59，12 143.31和15 154.57 m3／s。距坝址1 350.00，10 500.00，16 400.00和22 050.00 m的断面洪峰流量下降明显，距坝址27 960.30m的断面由于近黑河入渭口，河宽增加致使最大流量增大。一维模拟仅限于河道内部计算，虽然通过概化河岸地形以保证模型运行的稳定性，但计算精度上还是无法平衡。通过将一维、二维模型进行耦合，使河道内的水流可以漫过河堤参与二维计算，从而大大增加计算的准确性，同时减少由于计算结果较大带来的防洪资源浪费。

图15 典型断面水位变化过程线Fig.15 Water⁃level variation hydrographs of typical sections

图16 观测断面流量变化过程Fig.16 Hydrograph of variation in discharge

通过耦合模拟得到的大坝溃决（t＝0 h）后洪水演进1，2，3 h时的水深等值线图和流场图。计算结果表明，耦合模拟结果符合溃坝洪水演进的基本规律，洪水在河道内的演进过程清晰。图17为耦合后黑河金盆水库下游A（519 741.13，3 771 313.12），B（521 351.06，3 774 972.38），C（523 604.26，3 779 566.97），D（526 438.69，3 783 541.46），E（529 349.24，3 783 083.32）5个观测点的水深和流速变化过程。由图17可知，溃坝洪水到达A，B，C，D，E的时间分别为模拟初始时刻后的16.40，16.44，16.79，17.39和17.60 h；各点水深及流速分别于17.23，17.62，17.77，18.14和18.28 h达到峰值；最大淹没水深分别为7.85，5.97，5.47，7.02和6.65 m；最大流速分别为9.73，5.76，5.98，3.72和4.53 m／s。单纯的二维模型不能精确配准河道的断面信息，河岸处的网格生成也会影响河道内水流运动计算的准确性，通过一、二维耦合模拟，可以明确河道内部的水流运动环境，同时通过连接点进行交互计算，提高模型计算的准确性，从而减少计算结果偏小导致防洪不利的不良影响。

图17 各采样点的水深和流速变化过程Fig.17 Variation hydrographs of water depth and velocity of sampling points

5 结 语

本文以黑河金盆水库大坝作为研究对象，建立了溃口近区二维模型和下游洪泛区耦合模型，模拟了溃口发展和水流形态，补充完善洪水演进模拟的研究，模拟了河道内及计算区的洪水演进，为其溃坝风险和灾害损失评估提供了依据，主要结论如下：

（1）DAMBRK计算逐渐溃坝得到两种溃决历时Tf＝0.906 1 h和Tf＝0.897 5 h时的溃口发展过程和流量变化过程；当Tf＝0.897 5 h时，溃口下泄流量峰值达到84 334.09 m3／s，溃口底部高程和宽度达到528.7和190.534 m，并应用此结果作为溃口近区模拟和下游洪水演进耦合模拟的边界条件。

（2）溃口近区二维模拟得到溃口扩展过程中计算区域某一时刻的水深及流速分布；观察到初始时段溃坝水流在溃口两侧形成漩涡对其冲刷，这种漩涡冲刷作用随溃口扩展由轻到重发生变化，最终随溃口形成而缓慢下来，最后几乎看不到涡漩，此时认为形成最终溃口，但由于溃口的限制，水流流速仍呈扩散状从溃口流出；因此得到溃口的形成过程不仅包括漫顶水流的直接作用，同时包括溃口形成过程中两侧漩涡状水流的反冲刷作用，从而进一步解释溃口的发展机理。

（3）采用侧向连接方式的耦合模型模拟下游溃坝洪水演进，得到一维河道内不同断面水位和流量过程，以及二维淹没区域内瞬时水深和流速分布。耦合模型可以同时兼顾河道内的水流变化以及河道外计算区域内的洪水演进过程，减少由于计算结果偏大或偏小所带来的防洪资源浪费和防洪措施不利等不良影响。在实际工作中，对于河网较为复杂的地区，建议采用耦合模型进行计算，可以配置更准确断面信息，并且大大加快计算速度，实现根据每条河流及地形的特殊性，设置不同糙率及相关系数，达到提高模拟精度的目的。

［1］MORRISM W，HASSAN M A A M，VASKINN K A，etal.Breach formation：Field testand laboratory experiments［J］.Journal of Hydraulic Research，2007，45（Suppl1）：9⁃17.

［2］GAUCHER J，MARCHE C，MAHDIT，et al.Experimental investigation of the hydraulic erosion of noncohesive compacted soils［J］.Journal of Hydraulic Engineering，2010，136（11）：901⁃913.

［3］DUPONT E，DEWALSB J，ARCHAMBEAU P.Experimental and numerical study of the breaching of an embankment dam［C］∥Proceedings of the 32nd IAHR Biennial Congress，2007：339⁃348.

［4］PICKERTG，WEITBRECHT V，BIZEBERSTEIN A，et al.Breaching of overtopped river embankments controlled by apparent cohesion［J］.Journal of Hydraulic Research，2011，49（2）：143⁃156.

［5］段文刚，杨文俊，王思莹，等.无黏性土坝漫顶溃决过程及机理研究［J］.长江科学院院报，2012，29（10）：68⁃72.（DUANWen⁃gang，YANGWen⁃jun，WANG Si⁃ying，et al.Process and mechanism of non⁃cohesive dam overtopping outburst［J］.Journal of Yangtze River Scientific Research Institute，2012，29（10）：68⁃72.（in Chinese））

［6］NOURANIV，HZKIMZADEH H，AMINIA B，et al.Implementation of artificial neural network technique in the simulation of dam breach hydrograph［J］.Journal of Hydroinformatics，2012，14（2）：478⁃496.

［7］方崇惠，方堃.瞬时溃坝最大流量计算新通式推导及验证［J］.水科学进展，2012，23（5）：721⁃726.（FANG Chong⁃hui，FANG Kun.A new calculation formula derivation and verification of themaximum flow of instantaneous dam break［J］.Advances in Water Science，2012，23（5）：721⁃726.（in Chinese））

［8］DEWALSB，ERPICUM S，DETREMBLEUR S，et al.Failure of dams arranged in series or in complex［J］.Natural Hazards，2011，56（3）：917⁃939.

［9］张社荣，张婷.考虑调洪的一维漫顶不溃与漫顶溃坝过程模拟［J］.自然灾害学报，2012，21（4）：33⁃39.（ZHANG She⁃rong，ZHANG Ting.A study of one⁃dimensional overtopping to dam break or not of consider flood［J］.Journal of Natural Disasters，2012，21（4）：33⁃39.（in Chinese））

［10］袁晶，张为，张小峰.可变网格下的堤防溃口展宽二维数学模型［J］.水科学进展，2013，24（3）：358⁃365.（YUAN Jing，ZHANGWei，ZHANG Xiao⁃feng.A two⁃dimensionalmathematicalmodel of burstbroadening under variable grid［J］.Advances in Water Science，2013，24（3）：358⁃365.（in Chinese））

［11］梁艳洁，余明辉，王党伟，等.无粘性土堤漫顶溃决溃口区域水力参数变化过程数模研究［J］.武汉大学学报：工学版，2012，45（2）：142⁃147.（LIANG Yan⁃jie，YU Ming⁃hui，WANG Dang⁃wei，et al.Research of hydraulic parameters change processmodulus of breach area of no clay dike overtopping outburst［J］.Engineering Journal of Wuhan University，2012，45（2）：142⁃147.（in Chinese））

［12］SHAKIBAEINIA A，JIN Y C.A mesh⁃free particle model for simulation of mobile⁃bed dam break［J］.Advances in Water Resources，2011，34（6）：794⁃807.

［13］缪吉伦，赵万星，黄成林.SPH法模拟立面二维溃坝流动应用研究［J］.重庆交通大学学报：自然科学版，2012，31（1）：121⁃123.（MIAO Ji⁃lun，ZHAOWan⁃xing，HUANG Cheng⁃lin.Simulation study on application of vertical two⁃dimensional dambreak flow in SPH method［J］.Journal of Chongqing Jiaotong University（Natural Science），2012，31（1）：121⁃123.（in Chinese））

［14］SUN Xiu⁃li，FENG Min⁃quan.One and two⁃dimension coupling simulation for dam⁃break flood［C］∥ Advances in Water Resources and Hydraulic Engineering⁃Proceedings of16th IAHR⁃APD Congress and 3rd Symposium of IAHR⁃ISHS（Ⅴ），2008：1806⁃1810.

［15］邓健，黄立文，赵前，等.基于一、二维水动力耦合模拟的三峡库区溢油预测模型研究［J］.武汉理工大学学报：交通科学与工程版，2011，35（4）：793⁃796.（DENG Jian，HUANG Li⁃wen，ZHAO Qian，et al.Study on prediction model of oil spill in Three Gorges reservoir area based on 1⁃D and 2⁃D hydrodynamic coupling simulation［J］.JournalofWuhan University of Technology（Transportation Science＆Engineering），2011，35（4）：793⁃796.（in Chinese））

［16］曲红玲.河道溃堤与溃堤波的一、二维耦合计算数值模拟［J］.水利水运工程学报，2007（4）：49⁃54.（QU Hong⁃ling.A numerical simulation of 1⁃D and 2⁃D coupling calculation of river dike and break waves［J］.Hydro⁃Science and Engineering，2007（4）：49⁃54.（in Chinese））

［17］GALLEGOSH A，SCHUBERT JE，SANDERS B F，et al.Structural damage prediction in a high⁃velocity urban dam⁃break flood：field⁃scale assessment of predictive skill［J］.Journal of Engineering Mechanics，2012，138（10）：1249⁃1262.

［18］DANG Chao，DING Yu，CHENG Zun⁃Lan，etal.A simplemodel for natural dam break：I formulation［J］.Advanced Materials Research，2011，243：838⁃843.

［19］刘薇，彭新德，张小兵.一、二维耦合水流数值模拟在防洪评价中的应用［J］.中国水运，2010，10（6）：112⁃113.（LIU Wei，PENG Xin⁃de，ZHANG Xiao⁃bing.1⁃D and 2⁃D coupling numerical simulation ofwater flow of application in evaluation of flood control［J］.ChinaWater Transport，2010，10（6）：112⁃113.（in Chinese））

［20］张大伟，程晓陶，黄金池，等.基于Godunov格式的溃坝水流数学模型［J］.水科学进展，2010，21（2）：167⁃172.（ZHANG Da⁃wei，CHENG Xiao⁃tao，HUANG Jin⁃chi，etal.Dam break flowmathematicalmodel based on Godunov format［J］. Advances in Water Science，2010，21（2）：167⁃172.（in Chinese））

［21］姜小明，李丹勋，王兴奎.基于黎曼近似解的溃堤一维－二维耦合数学模型［J］.水科学进展，2012，23（3）：214⁃220.（JIANG Xiao⁃ming，LIDan⁃xun，WANG Xing⁃kui.1⁃D and 2⁃D coupled dam break mathematicalmodel based on the Riemann approximate solution［J］.Advances in Water Science，2012，23（3）：214⁃220.（in Chinese））

［22］孙秀丽，冯民权，孙先锋.溃坝洪水的一维、二维耦合模拟［J］.华北水利水电学院学报，2008，29（4）：25⁃28.（SUN Xiu⁃li，FENG Min⁃quan，SUN Xian⁃feng.Simulation of one⁃dimensional，two⁃dimensional coupling of dam break flood［J］. Journal of North China Institute ofWater Conservancy and Hydroelectric Power，2008，29（4）：25⁃28.（in Chinese））

［23］李传奇，侯贵兵.一维二维水动力模型耦合的城市洪水模拟［J］.水利水电技术，2010，41（3）：83⁃86.（LIChuan⁃qi，HOU Gui⁃bing.Simulation of city flood coupled 1D and 2D hydrodynamic model［J］.Water Resources and Hydropower Engineering，2010，41（3）：83⁃86.（in Chinese））

［24］夏军强，王光谦，谈广鸣，等.动床条件下溃坝水流的二维水沙耦合模型［J］.水利学报，2011，42（3）：296⁃306.（XIA Jun⁃qiang，WANG Guang⁃qian，TAN Guang⁃ming，etal.The two⁃dimensionalwater and sediment coupledmodel for dam break flow inmoving bed conditions［J］.Journal of Hydraulic Engineering，2011，42（3）：296⁃306.（in Chinese））

［25］NZOKAO F，MORSE B，ROBERT J L，et al.Numerical characteristics of a coupled river ice and hydrodynamic model［J］. Canadian Journal of Civil Engineering，2011，38（4）：393⁃403.

［26］NZOKAO F，MORSE B，ROBERT J L，et al.Water wave transients in an ice⁃covered channel［J］.Canadian Journal of Civil Engineering，2011，38（4）：404⁃414.

［27］周孝德，冯民权，闵涛，等.黑河金盆水库洪水淹没范围及撤离方案研究报告［R］.西安：西安理工大学，2006.（ZHOU Xiao⁃de，FENG Min⁃quan，MIN Tao，etal.Reportof Heihe Jinpen reservoir flood submerged area and evacuation plan research［R］.Xi'an：Xi'an University of Technology，2006.（in Chinese））

［28］GOTTARDI G，VENUTELLI M.Central scheme for two⁃dimensional dam⁃break flow simulation［J］.Advances in Water Resources，2004，27（3）：259⁃268.DOI：10.1016／j.advwatres.2003.12.006.

［29］樊国瑞.石闸水库溃坝洪水分析研究［J］.电力勘测，1995（7）：25⁃30.（FAN Guo⁃rui.Analysis and research on the stone gate dam break flood［J］.Power Survey，1995，（7）：25⁃30.（in Chinese））

［30］邓昌奇，杨晓燕.土石坝漫顶逐渐溃坝模型的建立及应用［J］.吉林水利，2013（2）：36⁃42.（DENG Chang⁃qi，YANG Xiao⁃yan.Establishmentand application of embankment dam overtopping gradually break mode［J］.Jilin Water Resources，2013（2）：36⁃42.（in Chinese））

2D numerical simulation of breach area and coupling simulation of dam⁃breach flood

YANG Zhi，FENG Ming⁃quan
（State Key Laboratory of Eco⁃hydraulic Engineering of Northwest Arid Area in Shaanxi，Xi'an University of Technology，Xi'an 710048，China）

Based on the breach process of the clay⁃core wall sandy dam，as well as the characteristics of dam breach flood propagation andmovement，the 2⁃D numerical simulation of breach area and the coupling simulation of the dam⁃breach flood of Heihe Jinpen reservoir dam are presented in this paper.Themethod of DAMBRK is used to calculate gradual dam breach，and the results of DAMBRK model are used in the following simulation.A one⁃dimensionalmodel is solved by the Abbott⁃Ionescu six points implicit finite difference scheme.The unit centered finite volumemethod is adopted to solve a two⁃dimensionalmodel equation.To achieve the coupling of 1⁃D and 2⁃D model，the lateral connections are used to connect the level calculation points and the two⁃dimensional grid unit of both sides of Heihe River.The 2⁃Dmodel is developed to calculate and simulate the breach area.The instantaneous water depth and the velocity distribution of the calculation region are given by the 2⁃Dmodel.Taking the flood with frequency of 0.01%for example，the simulating results，such as water level and discharge process of each cross⁃section within the one⁃dimensional river course and water depth profile，velocity vector diagram and submerged range changing process of two⁃dimensional computational domain are presented by the coupling simulation.The formation process of the breach area includes notonly the directeffects of overtopping flow，butalso the anti⁃erosion effect of the spiral flow on both sides in the process of breach.Both flows change inner river and flood evolution of the calculation area in the outer river，which could be simultaneously taken into accountby the coupling simulation. In conclusion，the adverse effects from wasting flood resources and the inappropriate flood protection measures which are brought by too large or too small calculation results have been reduced.In engineering practice，for complex river network region，it is recommended to apply the coupling model in calculation in order to improve simulation accuracy.

breach area；2⁃D simulation；dam break flood routing；lateral connection；coupling simulation

TV122＋.4

A

1009－640X（2015）01－0008－12

10.16198／j.cnki.1009－640X.2015.01.002

2014－06－26

国家重点基础研究发展计划（973计划）项目（2012CB723201）；中央财政支持地方高校发展专项资金特色重点学科资助项目（106－00X101，106－5X1205）

杨 志（1989－），女，山西永济人，博士研究生，主要从事水环境模拟与监测工作。E⁃mail：19890517yz＠sina.cn通信作者：冯民权（E⁃mail：mqfeng＠xaut.edu.cn）