【论文分享】猴子岩水库色玉滑坡涌浪灾害链CFD–DEM耦合数值模拟发表时间：2024-11-21 16:02 论文信息 题目：猴子岩水库色玉滑坡涌浪灾害链CFD–DEM耦合数值模拟； 作者：肖华波，王泽皓，石伟明，王东坡，欧阳朝军； 期刊来源：工程科学与技术（EI）； 发表时间：2023.01
一、 研究概述 滑坡涌浪是一种常见的灾害现象，由于滑坡体与水体之间存在复杂的流固耦合作用，使得传统的单一介质模型无法进行准确求解。为此，在前期研发Massflow地表灾害动力学过程数值模拟软件的工作基础上，提出一种基于计算流体力学方法（CFD）与离散单元法（DEM）的流固耦合模型CFD–DEM，采用计算流体力学方法（CFD）求解水体流动，采用离散单元法（DEM）模拟散粒体滑坡运动，充分利用不同计算模型的优势，对滑坡及涌浪演进过程进行数值模拟分析。首先，利用该耦合模型对Robbe–Saule开展的颗粒堆积体坍塌–涌浪试验进行了相同工况下的数值计算，从颗粒坍塌运动过程、涌浪高度演化过程等方面进行对比，结果表明模拟结果与试验结果吻合良好，验证了CFD–DEM流固耦合模型的有效性。然后，将该方法应用于四川省猴子岩水库色玉滑坡–涌浪灾害的演进过程分析，重现了该事件滑坡失稳运动、涌浪产生及传播、涌浪爬升、涌浪回流的全过程。计算结果显示：计算得到的电站进水口处涌浪高度与实测数据较为接近；色玉滑坡从失稳运动至静止堆积的持续时间约为20 s，颗粒平均速度最大达到16.12 m/s；滑坡引起的涌浪约在滑坡失稳10 s后传播到对岸，之后开始沿坡面向上爬升，最大爬升高度达到27.32 m。研究表明CFD–DEM流固耦合模型能够很好地应用于模拟山区河谷大规模滑坡涌浪灾害，可为库区防灾减灾提供高效的技术支持。 二、 研究区及数据 猴子岩水电站位于四川省甘孜藏族自治州康定县境内大渡河干流上；色玉堆积体位于电站库区右岸，距下游大坝约4.3 km。色玉堆积体顺河长约650 m，横河宽约350 m，前、后缘高程分别为1820、2100 m，高差约为280 m。高程约1930～2100 m处为较宽缓的台状地形，台地前缘至坡脚地形坡度约为50～70°。堆积台地主要由冰碛堆积物（glQ3）组成，台地后缘覆盖少量崩坡积堆积物（col+dlQ4）。堆积体厚度较大，一般为83～100 m；下伏基岩主要为泥盆系中上统河心组（D2–3h）灰、灰白色中厚层白云岩、白云质灰岩夹灰岩。 2017年10月24日，色玉堆积体前缘发生滑坡，此时水库蓄水位为1832 m。滑坡失稳后的地貌特征如图1所示。由图1可知，失稳方量约80×104 m3，失稳高度最大达160 m，失稳后在下游3.7 km处的电站进水口形成1.7 m高的涌浪。截止目前，该堆积体又多次发生垮塌，并形成多处凹槽地形。由于该堆积体距离大坝相对较近，其失稳涌浪对工程和周边人类活动安全的不利影响需重点关注和研究。  图1色玉滑坡全貌 三、 数值模拟结果分 本文利用开源程序CFDEM、OpenFOAM和LIGGGHTS，基于计算流体力学方法与离散单元法构建流固耦合模型；将该耦合模型应用于真实情况下三维大尺度滑坡–涌浪灾害动力过程演进模拟；以四川省猴子岩水库色玉滑坡为例，对滑坡涌浪灾害的动力学过程进行分析，一次解决滑坡变形失稳及涌浪形成与传播的问题，为山区大区域尺度大型滑坡–涌浪复合灾害的分析研究提供新的技术手段。 通过对色玉堆积体失稳前后地形数据进行差值分析，确定模拟的物源。计算使用了2万多个半径为2 m的球形颗粒，参考色玉滑坡附近区域土工试验测试结果，设置模型中颗粒间的泊松比和杨氏模量等参数。色玉滑坡数值计算模型见图2，计算域网格大小为40 m×40 m×40 m。为减小涌浪高度计算误差，将库水位上下20 m范围内的网格进行加密，网格大小为10 m×10 m×10 m。计算参数选取见表1。  图2 色玉滑坡三维模型（红色圆点为浪高监测点） 表1 色玉滑坡计算参数 色玉滑坡运动及涌浪传播动态演进过程模拟结果如图3所示，显示了0～30 s时间内滑坡及涌浪形态。由图3可知：从滑坡体失稳运动至静止堆积的持续时间约为20 s，从滑坡入水处涌浪形成到涌浪传播至对岸的持续时间约为10 s。t=5 s左右时，在滑坡颗粒冲击作用下，涌浪开始逐渐形成；t=10 s，大部分滑坡颗粒滑下，涌浪高度进一步增加，沿滑坡发生点呈圆弧状向外扩散；t=15 s，滑坡颗粒全部进入水库，近岸水体沿滑坡颗粒运动方向下陷，涌浪在对岸已经开始爬升；t=20 s，滑坡颗粒不再运动，堆积于水库底部，此刻对岸涌浪爬升到最大高度；t=20～30 s，涌浪开始衰退，并继续向上下游传播。  图3 色玉滑坡运动及涌浪传播动态演进过程 图4展示了滑坡颗粒速度和涌浪速度在滑坡从失稳运动到静止堆积这段时间内的变化情况。图5展示了所有滑坡颗粒的最大速度及平均速度随时间的变化曲线。由图4、5可知：随着滑坡运动不断发展，重力势能转换为动能，滑坡颗粒速度开始处于快速增加阶段，平均速度在6.0 s时达到最大值，为16.12 m/s；滑坡中部运动速度明显高于两侧运动速度，滑坡前部运动速度高于滑坡后部运动速度。之后，由于滑坡前缘颗粒持续滑入水库底部，受到地形影响，颗粒在水库底部逐渐减速并开始堆积，滑坡体的平均速度也开始减小。直至滑坡发生约20 s后，平均速度减到0.04m/s，可以认为此时滑坡完全停止运动。虽然滑坡整体运动速度在6 s后减小，但仍有部分颗粒在激烈碰撞下达到较高的运动速度，单个颗粒的最大速度达到44.29 m/s，发生在10.0 s左右。 图4 滑坡颗粒及涌浪速度  图5 滑坡颗粒速度随时间变化曲线 四、结论 1在Massflow地表灾害动力学过程数值模拟软件的研发工作基础上，提出了一种CFD–DEM耦合算法，流体的控制方程为连续性方程及动量方程，颗粒的运动由牛顿第二定律控制，通过基于拖拽力模型的动量交换项实现流体与固体之间的耦合。 2）基于颗粒堆积体坍塌试验，进行了不同颗粒体积工况下的模型验证，模拟结果与试验结果吻合较好，证明了耦合模型能准确模拟颗粒与流体的相互作用过程。 3）将耦合模型应用于实际发生的大规模滑坡入水–涌浪灾害链案例——四川省猴子岩水库色玉滑坡涌浪灾害的模拟。结果表明：滑坡体从失稳运动至静止堆积的持续时间约20 s，颗粒平均速度最大达到16.12 m/s；滑坡引起的涌浪在约10 s后传播到对岸，之后开始沿坡面向上爬升，最大爬升高度达到27.32 m。CFD–DEM耦合模型不仅对室内坍塌–涌浪试验的模拟具有高效性和准确性，也可以应用于山区河谷大规模滑坡入水–涌浪灾害链这一复杂流固耦合问题的模拟，实现滑坡失稳运动及涌浪形成传播全过程分析，为山区河谷防灾减灾提供新的技术手段，为Massflow软件提供了一种新的功能扩展和有力补充。） |
