2017 
2. 华南理工大学 亚热带建筑科学国家重点实验室，广东 广州 510640
水库大坝是人类历史上重要的水工建筑物，在兴利和调蓄洪水等方面发挥着不可替代的作用，同时对区域内国民经济发展也承担着重要角色[1]。大坝一旦发生溃决，将会给城市带来灾难性破坏，造成不可估量的损失[2-3]。近些年随着计算机和数值计算方法的迅速发展，溃坝洪水演进在复杂计算域和地形条件下的高性能数值模拟研究成为国内外学术界和工程界关注的前沿研究领域之一[4]。如何使用计算机仿真技术科学准确地预测、模拟和显示洪水淹没范围，对于防洪救灾和损失评估具有十分重要的意义。不少学者在这方面进行了深入研究并取得了一定成果[5-9]。一维模型对溃坝洪水演进进行模拟研究时，虽计算效率高、适应性好，但在模拟洪水漫过河堤进入洪泛区演进时具有很大难度。二维模型则在处理复杂水流在无固定路径的地表演进时，能够描述洪水淹没的时间、范围、流速、水深分布等信息，且具有良好精度；但在模拟洪水在河道内演进时，计算效率及适应性远不如一维模型且需要详尽的河道地形资料。鉴于此，本文以深圳市龙华新区上游民治水库为研究对象，利用MIKE FLOOD将一、二维模型进行动态耦合，发挥一、二维模型模拟洪水在河道中及地表演进时各自的优势，同时回避其计算效率、网格精度及适用性等问题，更为合理真实地模拟上游民治水库大坝溃决后溃坝洪水的演进过程。
1 模型原理MIKE FLOOD是一个耦合的水力模型，能够完整模拟一维地下排水管网系统水流过程和二维地表漫流过程。MIKE FLOOD集成了MIKE Urban CS，MIKE11及MIKE21三个独立的软件模块。根据不同的应用情境可将其中的MIKE Urban CS或者MIKE11与MIKE21进行动态耦合，耦合后的MIKE FLOOD具有以下优点：能够在暴雨洪水过程中，模拟一维河道洪水超出河堤后，从决堤口向淹没区域的二维洪水演进过程；能够模拟一维河网与二维模拟区域连接处的动量传输；能够与外部软件如ArcGIS等进行复杂的数据交互，通过外部专业软件对其输入数据进行预处理等[10]。
MIKE11用于模拟一维河道水体的流态，集成了水动力模块、降雨径流模块、构筑物模块以及溃坝模块等，几乎涵盖了河流模拟的各个方面，在模拟过程中采用六点Abbott-Ionescu有限差分格式对圣维南方程组求解[11]。MIKE21属于平面二维自由表面流模型，广泛运用于海洋、湖泊、河道及蓄滞洪区的流场、流速、水位等方面的模拟，能够获得不同水文要素的时空分布及洪水淹没信息[10]。
本文利用MIKE11与MIKE 21模型在MIKE FLOOD平台上进行动态耦合，采用侧向连接方式[12](见图 1)。侧向连接即河道通过两岸与二维模型计算区域进行水流交换，它允许MIKE21网格单元从侧面连接到MIKE11的部分河段甚至整个河段。
|
图 1 侧向连接示意 Figure 1 Lateral connection diagram |
侧向连接中一维河道和二维网格单元间的水量交换采用堰流公式近似计算[13]，计算方法如下：
| $ \mathit{q}{\rm{ = }}\mathit{WC}{\left( {{\mathit{H}_{{\rm{us}}}}{\rm{ - }}{\mathit{H}_{{\rm{ds}}}}} \right)^\mathit{k}}{\left[{{\rm{1-}}\left( {\frac{{{\mathit{H}_{{\rm{ds}}}}{\rm{-}}{\mathit{H}_{\rm{w}}}}}{{{\mathit{H}_{{\rm{us}}}}{\rm{-}}{\mathit{H}_{\rm{w}}}}}} \right)} \right]^{{\rm{0}}{\rm{.385}}}} $ | (1) |
式中：q为交换水量(m3/s)；W为宽度(m)，一般取单元格和河道相连的边长；C为堰流系数；k为堰指数；Hus和Hds分别为堰上游和下游水位(m)；Hw为堰顶高程(m)。本次研究选择MIKE11中的河岸标记作为堰的位置。
2 模型构建 2.1 研究区域民治水库位于深圳市龙华新区民治街道民治社区东南方约2.5 km，在2个小型水库民乐水库及雅宝水库(原名油柑水库)的下游，包含2个小水库在内的总集雨面积为4.5 km2。民治水库正常库容271万m3，总库容400.2万m3，正常蓄水位79.58 m，校核水位82.36 m，死水位69 m，坝顶长320 m，其设计标准为50年一遇设计，100年一遇校核。民治街道位于民治水库下游，为龙华新区与南山、福田、罗湖、龙岗五区交汇处。下辖9个社区工作站，21个社区居委会，总面积26 km2，人口42.7万。根据地形条件以及可能的淹没区域，确定研究范围扩展至民治水库下游4.25 km2。为分析受灾情况选取4个特征点进行分析比较(图 2)，分别为民治水库和民治河上游间重灾居民区横岭村(特征点1)，民治河道中游西侧向南村(特征点2)，民治河下游水尾新村附近(特征点3)，以及民治河道中游段居民区和商业区附近(特征点4)。
|
图 2 研究范围示意 Figure 2 Schematic diagram of study area |
利用MIKE11构建一维模型，从ArcGIS中提取河道的相关地理信息，按MIKE11河网文件的格式导入MIKE11中，并设定河道名称、里程数、地理标识等信息，制成河网文件。从河道横断面CAD图中提取断面资料，并根据里程数标志在河道上，对河道赋予断面信息。由于民治水库通过一条暗涵与下游民治河相连，其过流能力相对于洪水在地表演进时的流量很小，因此在进行本次溃坝洪水模拟时，根据其过流能力，在河道上游边界给定3 m3/s的恒定入流，河道水位作为下游开边界。初始水深设为1 m，初始流量为0.1 m3/s，民治河断面较为规整，水流较为通畅，两侧岸壁为土砂或石质，由于没有实测资料加以验证，本文根据天然河道糙率取值经验，取糙率为0.025。将上述所有文件导入模拟文件编辑器，设定模拟时间步长、结果输出文件名等，运行MIKE11，至此一维模型建立完毕。
2.3 二维模型利用MIKE21构建二维模型，根据计算范围提取计算边界以及河道两岸边界，制作成Mesh文件所需要的格式，导入非结构化网格生成器中。本文采用非结构三角形网格，网格划分完成后，节点数为9 314个，网格总数为17 911个，最大网格为400 m2，最小网格为79.2 m2，平均网格大小237 m2。导入研究区域的高程散点，根据网格划分进行插值，得到模拟区域的地形文件。研究区域边界取无滑移闭边界，闭边界不与外界进行水量交换，即法向流速为0，上边界以溃口流量过程作为流量边界条件。
2.4 MIKE FLOOD耦合模型将独立的MIKE11模型文件和MIKE21模型文件建立完成并可正常运行之后，将两者在MIKE FLOOD平台上进行耦合连接。首先导入需要连接的两个模型文件，然后选择侧向连接方式。非结构化网格与MIKE11进行侧向连接时，需要设定连接处的耦合线，以确定MIKE21与MIKE11的连接位置。本文将二维地形中一维河道的边界坐标导出，制作成耦合线所需要的格式，保证水量交换点即为河道两侧堤岸，同时，沿程河道的所有计算点都将与MIKE21相连。通过耦合线进行连接后，一维河道模型左侧岸线与72个MIKE21网格单元进行连接，右侧岸线与69个MIKE21网格单元进行连接。
2.5 溃口流量计算鉴于本次模拟的民治水库主坝为土石坝，发生逐渐溃的概率较大，同时由于潜在的地震、滑坡等因素引发的瞬间溃的危害更大，因此，本文采用瞬间溃和逐渐溃两种溃决方式进行模拟。对于土石坝瞬间溃的溃口流量计算，根据相关资料[14]，本文采用下式进行溃坝洪峰流量计算：
| $ {{\mathit{Q}}_{\rm{max}}}\rm{=}\frac{8}{27}\sqrt{\mathit{g}}{{\left( \frac{\mathit{B}}{\mathit{b}\rm{m}} \right)}^{\rm{0}\rm{.25}}}{{\mathit{b}}_{\rm{m}}}\mathit{H}_{0}^{1.5} $ | (2) |
式中：Qmax为溃决最大流量(m3/s)；g为重力加速度(m/s2)；B为坝长(m)；bm为最终溃口宽度(m)；H0为溃决时水深(m)。
对于逐渐溃的溃口流量计算，假设大坝溃决时，坝体先渗透变形破坏发展到一定程度后再瞬间局部破坏。坝体渗透变形先从小孔开始，逐渐扩展，此时溃坝洪水流量根据文献按下式计算[15]：
| $ {{\mathit{Q}}_{\rm{max}}}\rm{=}\mathit{A}{{\left[\rm{2}\mathit{g}\left( \mathit{H-}{{\mathit{H}}_{\rm{P}}} \right)\rm{/}\left( \rm{1+}\mathit{fL}\rm{/}\mathit{D} \right) \right]}^{\rm{0}\rm{.5}}} $ | (3) |
式中：H为水库水位(m)；A为水流流过小孔的断面面积(m2)；Hp为溃口小孔处的高程(m)；f为达西摩擦系数，根据Moody曲线由d50粒径计算；L为小孔沿水流方向的长度(m)；D为小孔宽度(m)。假设小孔初始形状为矩形，宽度逐渐线性扩展到30 m后再瞬间部分溃决，溃决后的流量再根据式(3)进行计算。
上述两种溃决方式的溃口流量过程均采用下式计算：
| $ \mathit{Q}\rm{=}{{\mathit{Q}}_{\rm{max}}}{{\left( \frac{{{\mathit{Q}}_{\rm{max}}}}{\rm{5}\mathit{V}}\mathit{t}\rm{-1} \right)}^{\rm{4}}} $ | (4) |
式中：Q为溃口流量(m3)；V为库容(m3)；t为时间(s)。
分别模拟计算4种工况下溃口流量过程线以及下游洪水演进过程，具体工况设定见表 1。
|
表 1 工况设定 Table 1 Working conditions |
根据表 1中的4种不同工况，采用式(2)~(4)，分别推算出各种工况下溃口洪水流量过程线，见图 3。由图 3可知，对于瞬间溃，当溃口宽度为100和320 m时，最大下泄流量分别为4 220, 10 098和13 912 m3/s，库容水量排空历时约52, 22和23 min。可见在溃坝初始时刻流量即达到最大，随着溃坝洪水下泄，水库水位不断下降，洪水流量也快速下降，并且溃口宽度越大，溃坝初始时刻的最大洪水流量也越大，其下泄时长越短。而对于逐渐溃，当溃坝开始后，随着渗透变形管道小孔不断扩大，下泄流量也不断增加，溃口发展到一定程度后管道上部坝体开始垮塌，类似于瞬间溃，溃口流量迅速增大达到最大值，接着流量变小，直至库区水量排空，其中渗透变形破坏阶段最大流量为1 206 m3/s，之后迅速垮塌，最大流量达3 632 m3/s，总历时约68 min。各工况下洪水流量参数见表 2。
|
图 3 各工况下溃口洪水流量过程线 Figure 3 Line break flood process under different conditions |
|
|
表 2 各工况下洪水流量参数 Table 2 Flood flow parameters under different conditions |
为了更好地体现洪水演进过程对下游的影响规律，在本次洪水演进结果分析中，选取了洪水演进过程中受灾情况比较严重以及居民集中的4个特征点进行分析比较，特征点位置详见图 2。选择淹没情况最为严重的工况3进行分析，其淹没范围见图 4。
|
图 4 瞬间全溃淹没范围 Figure 4 Instantaneous full collapse flooded range |
由图 4可知，由于地形原因，河道左侧淹没情况较为严重，溃坝开始的1 min内，下泄洪水量极大且流速快，上游居民区迅速被淹没，10 min内，淹没范围已经到达民治河中游段向南村、居民区及商业区，并继续向下游演进。30 min后，洪水开始逐渐消退，由于地势较高，洪水在下游部分地区的河道内演进时，没有漫出。各特征点在各工况下水深变化见图 5。
|
图 5 不同工况下各特征点水深变化 Figure 5 Water depth change of four feature points under different conditions |
由图 5可知，瞬间溃各特征点水深在上游洪水到达时迅速增加，并且随着水库计算水位增加，各特征点处有1~2 m的水深增加。在工况3情况下，洪水历时11 s到达离坝址较近的横岭村附近，3 min左右水深达到最高6.15 m，之后水深迅速降低，最终有少量积水。洪水历时4~5 min到达民治河中游段居民区和商业区，水深最高分别可达5.3和4.1 m。洪水历时约10 min到达民治河下游水尾新村附近，水深最高达6.1 m。民治河道中游西侧向南村由于地势较低，在此次溃坝模拟中为重灾区，其水深可达9.6 m，并且在洪水消退后，仍有3.1 m左右积水。
相对于瞬间溃，逐渐溃刚开始流量很小且增长缓慢，溃坝洪水7 min左右到达距离坝址处较近的横岭村，37 min左右水深达到最大值4.07 m，之后洪水慢慢消退，最终存在0.14 m积水。其余特征点变化较为一致，由于溃决初期洪水流量较少，远离坝址处的特征点附近几乎不受洪水影响，当发展成瞬间部分溃时，才会出现水深变化，但此时水库水位已经降低，因此各特征点的水深相较于其余工况下有明显下降，其到达各特征点时间较晚，平均有30 min左右的推迟。
3.2 淹没区洪水流速分析溃坝洪水的流速大小一定程度上反映了洪水的破坏能力，对淹没较为严重以及财产生命较为集中的4个特征点的流速进行分析，结果见图 6。
|
图 6 不同工况下各特征点流速变化 Figure 6 Velocity change of four feature points under different conditions |
由图 6可知，随着洪水向下游演进，沿程特征点的最高流速逐渐下降。其中，在工况3情况下，靠近坝址的横岭村附近流速最高可达10.5 m/s，而远离坝址，在民治河下游水尾新村附近的最高流速小于3 m/s。根据相关文献[16]，当水流流速大于5 m/s时，有较强破坏力，会对建筑物产生破坏，因此横岭村附近上游位置会因流速过大，对建筑物造成破坏。而相对于瞬间溃，逐渐溃在各特征点处产生的流速均小于5 m/s，造成的危害比瞬间溃小。
4 结语(1) 基于MIKE FLOOD平台将MIKE11和MIKE21模型进行动态耦合，可以较好地模拟溃坝洪水演进过程，且一、二维模型在连接处地形的契合度、准确性及网格划分的合理性是其模拟精确的重要保证。
(2) 采用瞬间溃及逐渐溃两种溃决方式，设定4种不同工况，得出各工况下民治水库溃坝流量过程线，瞬间溃的洪峰流量较大，出现在溃坝开始时刻，之后随着库区水位逐渐降低，下泄流量逐渐减小，直至库区水体排空；逐渐溃的洪峰流量相对较小，出现在渗透破坏变形发展至上部坝体坍塌时刻，随后下泄流量变小，直至库区水体排空。
(3) 根据淹没水深以及生命财产集中程度选取4个特征点进行分析比较，分析表明溃坝洪水对上游地区横岭村附近破坏力较大，淹没水深较深，同时其受灾速度较快，应做好相应预防措施。民治河道中游段居民和商业区附近洪水流速接近5 m/s，对建筑物有一定破坏力，左侧向南村地势较低，淹没水深情况最为严重，并且在洪水消退后仍有3 m左右积水。民治河下游地区在洪水消退后也有少量积水。
| [1] |
林金波, 金生, 丁伟业. 基于HydroInfo软件的溃坝水流模拟[J]. 水利与建筑工程学报, 2015, 13(6): 113-117. ( LIN Jinbo, JIN Sheng, DING Weiye. Dam break water flow simulation based on HydroInfo software[J]. Journal of Water Resources and Architectural Engineering, 2015, 13(6): 113-117. (in Chinese)) |
| [2] |
曹伟. 水库溃坝数值模拟及风险分析[D]. 太原: 太原理工大学, 2015. (CAO Wei. Study on dam-break numerical simulation and risk analysis[D]. Taiyuan: Taiyuan University of Technology, 2015.(in Chinese)) http://cdmd.cnki.com.cn/Article/CDMD-10112-1015604123.htm
|
| [3] |
李云, 李君. 溃坝模型试验研究综述[J]. 水科学进展, 2009, 20(2): 304-310. ( LI Yun, LI Jun. Review of experimental study on dam-break[J]. Advances in Water Science, 2009, 20(2): 304-310. (in Chinese)) |
| [4] |
任海. 溃坝洪水演进数值模拟分析研究[D]. 天津: 天津大学, 2012. (REN Hai. Numerical simulation on the flood evolution process due to dam break[D]. Tianjin: Tianjin University, 2012.(in Chinese)) http://cdmd.cnki.com.cn/Article/CDMD-10056-1013041374.htm
|
| [5] |
吴天蛟, 杨汉波, 李哲, 等. 基于MIKE11的三峡库区洪水演进模拟[J]. 水力发电学报, 2014, 33(2): 51-57. ( WU Tianjiao, YANG Hanbo, LI Zhe, et al. Modeling of flood routing for Three Gorges reservoir area based on MIKE11[J]. Journal of Hydroelectric Engineering, 2014, 33(2): 51-57. (in Chinese)) |
| [6] |
丁志雄, 李纪人, 李琳. 基于GIS格网模型的洪水淹没分析方法[J]. 水利学报, 2004(6): 56-60, 67. ( DING Zhixiong, LI Jiren, LI Lin. Method for flood submergence analysis based on GIS grid model[J]. Journal of Hydraulic Engineering, 2004(6): 56-60, 67. (in Chinese)) |
| [7] |
可友国, 雷宏军, 王永高, 等. 尖岗水库溃坝洪水计算与风险评估[J]. 人民黄河, 2008, 30(7): 36-37. ( KE Youguo, LEI Hongjun, WANG Yonggao, et al. Study on dam-break flood calculation and risk assessment for Jiangang reservoir[J]. Yellow River, 2008, 30(7): 36-37. (in Chinese)) |
| [8] |
郭凤清, 屈寒飞, 曾辉, 等. 基于MIKE21 FM模型的蓄洪区洪水演进数值模拟[J]. 水电能源科学, 2013, 31(5): 34-37. ( GUO Fengqing, QU Hanfei, ZENG Hui, et al. Flood routing numerical simulation of flood storage area based on MIKE21 FM model[J]. Water Resources and Power, 2013, 31(5): 34-37. (in Chinese)) |
| [9] |
SHIRVAN M K M, ALAMATIAN E, BAFTI F G, et al. Flow simulation of dam break and determining flooding zones using MIKE 21[J]. Ecology Environment and Conservation, 2013, 19(3): 679-684. |
| [10] |
王世旭. 基于MIKE FLOOD的济南市雨洪模拟及其应用研究[D]. 济南: 山东师范大学, 2015. (WANG Shixu. The study of simulation and application of rainwater and flood based on MIKE FLOOD model in Jinan City[D]. Jinan: Shandong Normal University, 2015.(in Chinese)) http://cdmd.cnki.com.cn/Article/CDMD-10445-1015600159.htm
|
| [11] |
魏凯, 梁忠民, 王军. 基于MIKE21的濛洼蓄滞洪区洪水演算模拟[J]. 南水北调与水利科技, 2013, 11(6): 16-19. ( WEI Kai, LIANG Zhongmin, WANG Jun. Flood routing simulation of MengWa detention basin based on MIKE21[J]. South-to-North Water Transfers and Water Science and Technology, 2013, 11(6): 16-19. (in Chinese)) |
| [12] |
喻海军. 城市洪涝数值模拟技术研究[D]. 广州: 华南理工大学, 2015. (YU Haijun. Research on numerical simulation technology of urban floods[D]. Guangzhou: South China University of Technology, 2015.(in Chinese)) http://cdmd.cnki.com.cn/Article/CDMD-10561-1015986910.htm
|
| [13] |
衣秀勇. DHI MIKE FLOOD洪水模拟技术应用与研究[M]. 北京: 中国水利水电出版社, 2014: 59-61. ( YI Xiuyong. Application and research of flood simulation technology through DHI MIKE FLOOD[M]. Beijing: China Water and Power Press, 2014: 59-61. (in Chinese))
|
| [14] |
落全富, 安莉娜. 青山水库溃坝洪水模拟计算[J]. 浙江水利科技, 2010(2): 17-19. ( LUO Quanfu, AN Lina. Simulated calculation of dam break flood for Qingshan reservoir[J]. Zhejiang Hydrotechnics, 2010(2): 17-19. (in Chinese)) |
| [15] |
周远方. 大南川水库溃坝的数值模拟研究[D]. 长沙: 长沙理工大学, 2010. (ZHOU Yuanfang. Numerical simulation study on dam-break flow Dananchuan reservoir[D]. Changsha: Changsha University of Science and Technology, 2010.(in Chinese)) http://cdmd.cnki.com.cn/Article/CDMD-10536-1011025696.htm
|
| [16] |
赖成光. 城市地区水库溃坝洪水演进数值模拟研究[D]. 广州: 华南理工大学, 2013. (LAI Chengguang. Study on numerical simulation of dam-break flood evolution in urban area[D]. Guangzhou: South China University of Technology, 2013.(in Chinese)) http://cdmd.cnki.com.cn/Article/CDMD-10561-1013319909.htm
|
2. State Key Laboratory of Subtropical Building Science, South China University of Technology, Guangzhou 510640, China