杜镇瀚,钟启明,董海洲,单熠博
(1.河海大学土木与交通学院,江苏 南京 210098;2.南京水利科学研究院岩土工程研究所,江苏 南京 210024;3.水利部水库大坝安全重点实验室,江苏 南京 210029)
在一定的地质地貌条件下,由于降水、地震等原因引起的山崩、滑坡、泥石流等堵截山谷、河道,造成贮水而形成的湖泊统称为堰塞湖;阻塞山谷、河道的堆积体称为堰塞体[1]。我国位于强地震活动带的青藏高原、中南及西南山区为堰塞湖的高发区[2],如2000年4月9日发生在青藏高原的易贡堰塞湖(图1),由5 520 m高程的雪山向下高速滑动,历时约10 min,滑程8 km,堆积于约2 190 m高程的易贡藏布江,形成坝高54 m、长2 500 m、最大蓄水量可达28.8亿m3、实际蓄水量达15.34亿m3的滑坡堰塞湖[3]。2008年5月12日,汶川大地震造成了257处堰塞湖,其中规模最大、潜在威胁最严重的是唐家山堰塞湖(图2),最大蓄水量可达3.16亿m3,实际蓄水量达2.3亿m3,而其形成恰逢雨季,增大了堰塞湖的溃决风险,对下游170万人民造成了巨大的威胁[4-5]。2018年10月10日和11月3日,金沙江上游白格地区分别发生了两次大规模山体滑坡,形成白格堰塞湖(图3),最大蓄水量分别达2.49亿m3和7.56亿m3,实际蓄水量分别为2.49亿m3和5.78亿m3,溃决后对下游造成了重大的经济损失[6-7]。
图1 易贡堰塞湖Fig.1 Yigong dammed lake
图2 唐家山堰塞湖 Fig.2 Tangjiashan dammed lake
图3 “11·3”白格堰塞湖Fig.3 “11·3” Baige dammed lake
堰塞体虽然由土石材料构成,但其与人工填筑的土石坝存在较大的差别[8]:①坝体形态方面,堰塞体沿河流运动方向一般较长,且顶部凹凸不平;②结构材料方面,大部分堰塞体结构复杂,不均匀性强,颗粒级配宽泛;③水动力条件方面,由于堰塞体上没有泄洪设施,极易发生漫顶溃决。据统计[9],84.4%的堰塞湖在形成1 a内发生溃决,风险远高于人工修筑的水库大坝。堰塞湖一旦溃决,往往造成严重的洪水灾害,给下游人民生命财产安全带来严重的威胁。合理评估堰塞湖的风险有利于科学的应急响应,减少社会财产损失,保障人民生命安全并降低环境生态破坏,因而具有重要的理论意义和实际价值。
我国堰塞湖风险评估工作从21世纪初开始兴起,经历2008年“5·12”汶川地震后得到了广泛关注。经过20余年的发展,理论和实践方面都取得了丰硕的成果,但是还未形成完整的体系和严谨的理论框架。本文从堰塞湖风险评估的概念入手,围绕堰塞湖灾害链的演化过程展开论述,对堰塞湖风险评估的研究现状进行分析总结并展望了堰塞湖风险评估的未来发展趋势。
关于风险的定义,不同领域的学者有着不同的理解,但通常从两方面特征进行描述:①可能性,即发生的概率;②后果,即造成的损失[10]。2009年,ISO31000:2009《风险管理国际标准》将风险定义为“一个事件后果与其发生可能性的组合”。我国现行的国家风险管理标准GB/T 24353—2009《风险管理原则与实施指南》中也采用了这一定义。根据该定义,风险可以用后果的期望值来量化,即灾害发生的频率和损失的乘积。对于堰塞湖,风险可表示为堰塞湖的溃决概率与溃决损失的乘积,即:
R=PC
(1)
式中:R为堰塞湖风险;P为堰塞湖溃决概率;C为堰塞湖溃决损失。
堰塞湖风险评估过程可分为3个主要步骤:①对堰塞湖进行危险性评价,得到堰塞湖溃决的概率;②通过堰塞湖溃决过程与洪水演进模拟获取洪水淹没情况;③堰塞湖溃决损失评估,即对淹没区的损失情况进行定量评估,将溃决概率与溃决损失相乘得到堰塞湖的风险,并进行综合评估。
对于堰塞湖而言,其溃决危险性评价可转化为挡水构筑物堰塞体稳定性评价。堰塞体稳定性主要与其物质组成[11]、结构特征[12]、堰塞湖水动力条件[13]以及当地的降水条件和次生灾害[14-15]相关。围绕着这些影响因素,国内外学者针对堰塞体稳定性评价方法展开了一系列研究,目前根据其研究特性可分为定性评价和定量评价。
堰塞体稳定性的定性评价,主要以其形成机制、物质组成以及坝体几何形态为依据,初步判断坝体的抗冲蚀能力,最终评价堰塞体的整体稳定性[16]。定性评价方法的主要特点是不经过数学计算,通过卫星遥感、无人机航测以及地形检测等技术手段,结合历史资料得到结论。
定性评价方法主要可以分为工程类比法和历史分析法。工程类比法是指将需要研究的堰塞体根据形成条件和地质条件与其他堰塞体进行类比,从而评价该堰塞体的稳定性;历史分析法是指通过对堰塞体本身的形成历史以及发展趋势进行预测分析,从而揭示其发展规律以评定堰塞体的稳定性。
自“5·12”汶川地震后,我国一批学者以及水利行业标准根据不同的评价指标对堰塞湖危险性进行了分级。选取的评价指标主要包括堰塞湖规模、流域面积、径流量、堰塞体高度、物质组成、坝体结构等,一般将堰塞湖危险等级分为极高危险、高危险、中危险和低危险4个等级(表1)。
近年来,随着卫星遥感监控、无人机航测、地表位移实时监测等技术的发展,许多学者虽然没有将堰塞湖危险性进行等级划分,但结合多方面因素实现了对堰塞湖地貌形态、物质组成、几何结构的快速测算和风险的实时分析预警(表1)。
表1 堰塞湖危险性定性评价方法
目前,单纯开展定性评价的研究较少,而定性评价大多作为定量评价的基础,也是堰塞湖风险评估的第一步,完整且迅速地获取地形信息和遥感数据是溃决洪水分析和损失评估的基础。
2.2.1 数理统计法
数理统计法是根据大量堰塞湖的资料数据确定参数,再提出数学表达式或判别准则,通过堰塞体的稳定性评价来确定堰塞湖的危险性。国外学者对该方面研究起源较早,1999年,Casagli 等[11]提出了利用基于堰塞体体积和堰塞湖流域面积等2个参数的堆积指标法来判断堰塞体的稳定性。其后,诸多学者也根据收集的堰塞湖基础资料,采用逻辑回归方法,提出一系列的判别方法,主要考虑的因素包括堰塞湖的水动力学指标、堰塞体的地貌学指标和物质组成等,见表2。
表2 基于数理统计的堰塞体稳定性评价方法
就目前研究现状来看,逻辑回归拟合分析是数理统计法的主流,但该方法依赖样本的可靠性及拟合因子的合理性和代表性,随着堰塞湖基础资料的不断丰富,稳定性评价的准确性也逐渐提高。
2.2.2 物理模型法
堰塞湖危险性评价物理试验主要以小尺度简化水槽试验为主[23-26]。近年来,国内学者基于不同的致灾因子,尝试了还原真实应力场的离心模拟试验[27]以及考虑地震荷载的振动台试验[28]。总的来说,物理模型法是根据堰塞体的基础资料,采用相似的物理材料制成模型来还原堰塞体的失稳过程,是评价堰塞体稳定性最直观的一种方法。但该方法一般将堰塞体作为形状规则、材料均匀的坝体处理,与真实的堰塞体有较大的差异,加上模型缩尺的问题,模拟的结果能否合理反映实际的情况仍值得深入研究。
2.2.3 数值分析法
数值分析法一般通过坝体信息和材料结构特征分析堰塞体的稳定性。目前,堰塞体稳定性分析常用的数值方法包括极限平衡法、有限单元法、有限差分法、离散单元法和非连续变形分析法等[16],这些方法属于确定性分析方法。但离散单元法和非连续变形分析法在堰塞体上的应用相对较少,表3为一些确定性数值分析法。近些年来,学者们尝试将不确定分析方法,主要包括模糊层次分析法、熵值法、模糊综合分析法等应用到堰塞湖的危险性评价研究中,但不确定性分析方法起步较晚,数学理论和工程实际的结合还存在许多问题,尚需改进。
表3 基于数值分析的堰塞体稳定性评价方法
从上述分析可以发现,确定性分析方法重点在于分析外荷载作用下堰塞体自身的稳定性状态,但大多采用人工土石坝的模拟方法,未能充分考虑堰塞体与土石坝不同的坝料组成和结构特征以及水动力条件。不确定性分析方法发展较晚,但是随着信息技术和人工智能的发展,该方法对于研究考虑各种不确定因素作用下堰塞湖的危险性有着良好的发展前景。
数值模拟是堰塞湖溃决过程预测的重要手段,总的说来,数学模型一般分为3类:第一类是参数模型,第二类是基于溃决机理的简化数学模型,第三类是基于溃决机理的精细化数学模型[45]。
3.1.1 参数模型
参数模型大多基于溃坝案例数据进行统计分析,采用经验公式计算溃坝相关参数(如溃口峰值流量、溃口尺寸、溃坝历时等)。由于参数模型公式简单、计算快捷,也常用于溃坝致灾后果的快速评价。1977年,Kirkpatrick[46]提出了第一个预测土石坝溃口峰值流量的经验公式,随后一些学者也提出了一系列关于土石坝溃决的参数模型[47]。由于结构和材料的复杂性,目前可专门用于堰塞湖溃决参数预测的参数模型较为缺乏。1985年,Costa等[48-49]在10个堰塞湖溃决案例的基础上,提出了3个预测溃口峰值流量的回归方程,通过堰塞体高度和堰塞湖蓄水量两个参数或其组合进行预测;随后,又建立了堰塞湖势能与溃口峰值流量之间的关系。1997年,Walder等[50]建立了堰塞湖下泄水量和水位下降高度与溃口峰值流量的关系。早期的参数模型仅简单考虑了堰塞湖的水动力条件,且只给出了溃口峰值流量的表达式。2012年,Peng等[51]开发了一种溃决参数快速预测模型,可综合考虑堰塞体形态、堰塞湖水动力条件以及堰塞体材料的冲蚀特性。利用该模型可以预测溃口峰值流量、溃口最终尺寸(顶宽、底宽、深度)以及溃坝历时。表4为典型的堰塞湖溃口峰值流量预测参数模型。
表4 堰塞湖溃口峰值流量参数模型
虽然参数模型可简单快捷地对溃坝参数进行预测,但却无法提供溃决洪水流量过程线。
3.1.2 简化数学模型
1965年,Cristofano[52]建立了第一个均质坝漫顶溃决数学模型,其后,各国学者提出了一系列模拟土石坝溃坝的数学模型[47]。
对于堰塞湖的溃决过程模拟,目前大多借助土石坝溃决模型,但由于堰塞体与人工填筑的均质坝在材料特性和结构特点上差异巨大,其溃决机理也明显不同[47,53]。近年来,有关学者在深入研究易贡、唐家山、小岗剑等拥有堰塞湖溃决实测资料案例的基础上,提出一系列考虑堰塞体物质组成、坝料冲蚀特性和坝体结构特征的基于堰塞湖溃决机理的简化数学模型(表5),并应用于2018年10—11月白格和加拉堰塞湖的应急处置。
表5 堰塞湖溃决过程简化数学模型
总之,简化数学模型的优势在于一定程度上考虑了堰塞湖的溃决机理,且计算速度较快,但大多数简化数学模型在溃口冲蚀过程存在假设,无法真正模拟溃口的演化规律及溃坝过程中的水土耦合作用。
3.1.3 精细化数字模型
为了描述溃坝过程中的水土耦合作用,近年来,随着计算机性能的提升及泥沙科学和计算流体力学的发展,出现了一些以非平衡坝料输移理论为基础、基于浅水假设的精细化溃坝数学模型。该类模型主要基于水流连续性方程、动量方程与能量方程,耦合泥沙运动方程,目前基于水动力坝料冲蚀方程的一维、沿深度平均二维和三维数学模型的研究取得了明显的进展,可以更详细地模拟土石坝的溃决过程[60-63]。为了处理不连续混合流态组成的漫顶水流,通常采用近似黎曼解法和全变差递减法(TVD)等激波捕捉方法,并采用有限体积法、水平集法、光滑颗粒流体动力学法等数值模拟方法对控制方程进行求解。该类方法可考虑溃坝过程中的水土耦合作用,实现对溃坝过程的精细化模拟,并可模拟复杂边界条件和溃口演化规律,目前普遍用于颗粒较为均匀的土石坝溃决过程的模拟,是堰塞湖溃决过程模拟的发展方向。
溃坝洪水兼具激波和稀疏波的特征,按离散基本原理可分为特征线法(MOC)、有限差分法(FDM)、有限元法(FEM)和有限体积法(FVM)。特征线法、有限差分法和有限元法在很多缓流问题的计算中取得了很大的成功,但不完全适合求解溃坝洪水的强间断流动。
有限体积法是20世纪80年代以来发展起来的一种新型微分方程离散方法,综合了有限差分法和有限元法的优点,对由一个或多个控制体组成的任意区域以至整个计算区域都严格满足物理守恒定律。有限体积法不是直接对方程组进行数值离散,而是从积分形式的守恒性方程组出发,采用非结构网格进行离散,在控制体边界上形成间断解的Riemann问题。在求出每个控制体边界沿法向输入(出)的流量和动量通量后,对每个控制体分别进行水量和动量的平衡计算,得到计算时段末各控制体平均水深和流速。
Riemann解的Godunov格式是目前求解大梯度流动的主流格式,应用较多的有Roe格式[64]、Osher格式[65]、HLL格式[66]、HLLC格式[67]等;另外,有限差分法捕捉激波的高分辨率方法以及许多格式也可直接利用到有限体积法中。
一些通用的商业软件如DHI-MIKE、InfoWorks、HEC-RAS等也可用来模拟溃坝洪水的演进过程。
我国学者自21世纪初开始对溃决损失进行系统的研究。刘永志等[68]从灾害链的角度对洪涝灾害进行了系统研究,总结了灾害链的意义;黄国如等[69]对洪涝灾害的风险指标进行了分析,总结了风险等级的区划方法。总体而言,溃决损失研究仍然处在探索阶段,且主要针对人工修筑的水库大坝,堰塞湖的溃决损失研究更是匮乏。但水库大坝与堰塞湖在溃决洪水演进时有着共通性,对于下游受灾群众而言,溃坝产生的危害以洪涝灾害的形式呈现,因此水库大坝溃决损失评价方法具有借鉴意义。参考水库大坝溃决损失的分类,将堰塞湖溃决损失评估分为生命损失评估和经济损失评估两方面。
生命损失是指溃决洪水淹没范围内因洪水灾害导致死亡的人口数。生命损失受到定性与定量因素的耦合作用,造成致灾因子体系间作用关系的不确定性与灰色性,增加了定量评估的难度[70]。当前的生命损失评估方法大多以主要影响因素为参数建立模型,Peng等[71-72]将生命损失评估模型分为3个类型,即经验模型、物理模型和折衷模型。
4.1.1 经验模型
经验模型是指通过对历史数据的统计,结合数学方法建立生命损失与影响因素之间的关系模型。1988年,Brown等[73]通过统计世界多国的溃坝历史数据,并考虑风险人口和警报时间等因素,最早提出了生命损失经验公式(B&G法)。随后,各国学者考虑了更多的影响因素,主要包括风险人口密度、警报时间、溃坝洪水严重性、风险人口对洪水理解程度、溃坝时间、天气、坝高、库容、风险人口距坝址距离、建筑物易损性和应急预案等,提出了一系列的经验公式对生命损失进行估算(表6)。
表6 生命损失评估方法
4.1.2 物理模型
物理模型是模拟人类在洪水中的行为,以个体为研究对象,探索人类在动水中的稳定程度,研究过程十分复杂,较少应用于实际案例之中,且人类在洪水中的行为与生命损失的关联仍需进一步研究考证。国外对物理模型的研究较多,国内目前缺少关于此类模型的相关研究(表6)。
4.1.3 折衷模型
折衷模型是介于经验模型和物理模型的一种模型,在一定程度上考虑了生命损失产生的机理,并以此将洪水淹没区域分为若干个子区域。在每个子区域内以历史数据或专家判断校准,从而得到生命损失和关键因素之间的关系。
总之,经验模型一般是宏观的,并未强调机理的研究,但是使用方便且模型发展较为成熟;物理模型主要注重于机理研究,研究个体生命在洪水灾害中的不稳定性,该类方法大多较为复杂,且无法考虑危险个体的主观因素;折衷模型可结合生命损失机理来研究损失率和洪水特性之间的关系,是未来生命损失评估方法的发展方向。近些年国外研究多集中于折衷模型,国内则多侧重于经验模型。在未来生命损失评估中,可以针对我国国情,增加对堰塞湖风险机理的研究,研究出适合我国实际情况的生命损失评价体系。
自20世纪70年代起,国内外针对溃决经济损失评估展开了一系列的研究,但关于堰塞湖溃决经济损失评估的研究相对比较匮乏,鉴于溃决洪水的相通性,可借鉴人工土石坝的经验[89]。将经济损失分为直接经济损失和间接经济损失,分析方法可分为数理统计法和模糊数学法,目前的研究均侧重于直接经济损失。
4.2.1 数理统计法
数理统计法是建立在大量调查数据的基础上,依赖于社会经济资料和各行业财产损失率资料的完整性,研究洪水深度、流速以及淹没时间等因素和损失之间的关系。数理统计法工作量较大、操作复杂且花费高昂(表7)。
4.2.2 模糊数学法
通过构建多指标数学模型来评估洪水经济损失的方法是目前学术界的主流。模糊数学法根据洪灾风险的形成机制,结合区域内的洪水特征、气象条件、社会经济分布情况等,构建基于数学模型的综合评估方法(表7)。目前模糊数学法实用性较广,且随着大数据时代的到来,其具有很大的提升空间。但是在构建模型的过程中,很多参数需要主观定义或简化处理,而经济损失是动态发展的,目前的损失研究主要停留在静态评估上,模糊数学法仍存在精度的问题。
表7 经济损失评估方法
目前的经济损失评估方法大多未考虑溃决洪水对生态环境的影响,但生态损失从长远来看也是一种经济损失,生态环境和经济之间的联系是不可忽略的。生态损失评估方法的研究难点主要在于:①如何将生态价值量化,从而清晰明确地对生态损失进行定量评估;②如何对风险区生态损失的影响因子进行分析与提炼,并量化影响因子,构建完善的指标体系。
a.对于堰塞湖而言,发生溃决灾害是一种概率性事件,而一种完整的灾害链体系不可能得到精确的溃决概率,因此需要深入研究堰塞湖的孕育-致灾机制,克服致灾因子分析的不确定性,结合大量堰塞湖基础数据,通过不断改进算法获取逼近真实的堰塞湖溃决概率。
b.堰塞湖溃决洪水模拟的准确性问题。堰塞湖溃决洪水模拟涉及复杂的水土耦合作用,由于堰塞体材料颗粒分布的不确定性,准确预测不同致灾因子作用下堰塞湖的溃决过程是下游洪水演进和灾害评估合理性的关键。应根据堰塞体的形成机理,充分利用物探手段,现场获取堰塞体的材料和结构特征。
c.堰塞湖溃决损失评估的定量化方法。由于溃决损失研究起步较晚,缺乏相关影响因素的历史数据,如何利用有限的数据对损失评价模型中的参数进行量化处理是目前的一大难题。随着大数据时代的到来和多学科交叉研究的兴起,各邻域的数据、方法相结合,加上算法的精进,一定程度上弥补了历史数据缺乏的问题。
d.堰塞湖风险定量评估体系。目前堰塞湖风险尚未形成一个完整的定量评估体系,危险性评价和灾害损失评估相互独立,且定量化程度不高。可以从人工智能的角度,通过智能网络探索合理的堰塞湖风险整体评估模型。
e.堰塞湖风险动态评估体系。堰塞湖风险是一个动态的过程,但目前的危险性评价和损失评价大多是静态分析,忽略了灾害本身的不确定性发展过程。后续研究可考虑人类行为和社会组织撤离的动力学特性,结合卫星遥感、无人机航拍等技术进行堰塞湖风险的动态评估。
本文基于堰塞湖灾害链的发展演化过程,从风险评估的概念入手,进行堰塞湖风险评估的定性和定量总结。在堰塞体溃决危险性评价方面,目前国内研究较多且近些年多集中于定量评价中,基于数理统计方法的参数模型以及结合信息科技和人工智能的不确定性方法是目前发展的主流。连接堰塞体溃决危险性和洪涝灾害损失之间的堰塞湖溃决过程以及洪水演进模拟是堰塞湖风险评价的枢纽,溃决模型的研究目前主要集中于简化模型上,可以在考虑堰塞湖溃决机理的情况下快速模拟溃口形态和流量演化过程;而洪水演进部分目前多借助于商业软件开展数值分析。在堰塞湖溃决损失评估方面多借鉴洪涝灾害的损失估算方法;在生命损失估算方面,国外研究多集中于折衷模型,而国内则多侧重于经验模型;在经济损失估算方面,研究多集中于直接经济损失,目前学术界的研究主流是构建多指标数学模型。
总体而言,堰塞湖风险评估的基本框架已初步形成,但仍需在堰塞湖危险性评价中的不确定性问题、堰塞湖溃决洪水模拟的准确性问题、堰塞湖溃决损失评估的定量化方法、堰塞湖风险定量评估体系、堰塞湖风险动态评估体系等方面开展系统深入的研究。