基于改进的FAST模型的千将坪滑坡涌浪研究

赵永波 ，王 健 ，黄波林 、3，谭建民 ，王世昌

ZHAO Yong-Bo1，WANG Jian2，HUANG Bo-Lin2、3,TAN Jian-Min1,WANG Shi-Chang1

（1.中国地质调查局武汉地质调查中心，武汉430205；2.防灾减灾湖北省重点实验室(三峡大学)，湖北宜昌443002；3.湖北长江三峡滑坡国家野外科学观测研究站，湖北宜昌443002）

(1.Wuhan Center of China Geological Survey,Wuhan,China,430205；2.Hubei Key Laboratory of Disaster Prevention And Mitigation，China Three Gorges University，Yichang 443002，Hubei,China；3.National field Scientific Observation Station of Landslide in The Yangtze Threee Gorges，Hubei Province，Yichang 443002，Hubei,China）

1 引言

2003年7月13日0时20分，中国三峡库区青干河支流千将坪滑坡突然下滑，造成岸坡上13人死亡，直接经济损失约8 000万元[1]。此后，千将坪滑坡作为水库蓄水后岩质顺层滑坡的典型案例被大量国内外学者进行了深入研究。研究者从滑坡的孕育环境、启动条件、高速滑动机制、变形破坏机理、稳定性评价和运动等方面进行了各种类型的物理试验和数值计算分析[2-14]。

尽管千将坪滑坡涌浪的危害在滑坡后就被注意到了，但与上述地质研究成果相比，千将坪滑坡造成的水动力和涌浪研究成果较少。千将坪滑坡涌浪数据大多来源于波浪痕迹调查和当地目击者的报告。千将坪滑坡发生时，最高涌浪达到39m，在水库传播达30 km之远，打翻青干河中22条渔船，11人落水失踪[15]。当前，尚无数值模型对千将坪滑坡造成的大面积、长距离的涌浪危害开展研究。

本文在前人FAST模型的基础上[16-24]，建立了浅水区千将坪滑坡涌浪的水波动力学模型，开展长距离大范围的涌浪危害研究，以再现千将坪滑坡涌浪过程。这一水波动力学模型将可用于水库区支流或干流浅水区滑坡产生的涌浪问题研究，为滑坡涌浪灾害预测预警提供重要技术支撑。

2 千将坪滑坡概述

千将坪滑坡位于中国三峡库区支流青干河右岸，其对岸为湖北省秭归县沙镇溪镇，该滑坡距离青干河入长江口3 km，距离三峡大坝约50 km（直线距离约40 km）（图1）。千将坪斜坡为顺向坡，坡度13°～35°，千将坪斜坡对岸为逆向陡崖，平均坡度75°，河谷呈明显不对称“V”型。未蓄水前，滑坡区青干河河谷宽度为50-80m，河床高程为90m左右。蓄水水位高程为135 m时，青干河水深约45 m，河面宽度约500 m。

图1 千将坪滑坡位置图Fig.1 Qianjiangpinglandslide location map

千将坪滑坡的剪出口在95-111 m之间[2]，滑带为区域层间剪切带演变而成[6,13]。滑带为黑色粘土和黄色粘土等泥化夹层组成[1]，有明显磨光现象。针对两层滑带粘土的环剪试验研究表明：在滑坡启动初期，滑体沿黄色粘土层发生滑动。当运动速度和剪切形变达到一定程度后，黄色粘土的阻滑力大于黑色粘土的阻滑力，滑体改为沿黑色粘土层下滑。当沿黑色粘土下滑时，由于阻滑力的大幅下降，高速滑动开始发生[1]。同时，高速运动后，千将坪滑坡滑体保持较好的完整度，这表明剪切变形大量发生在滑带中，滑坡体属于高速整体滑行而非高速流动[16]。

滑动后，根据滑坡滑动前后地貌、地物标志的遥感片对比，滑坡在主滑方向上的最大滑距约为250 m[8]。据三峡微震台网监测资料，千将坪滑坡最大滑速为16m/s。滑坡方量达大大滑速为16m/s。滑坡方量达1500万m3，滑坡滑入水体积约240万m3，滑坡体跨过青干河，滑舌直抵对岸，发生反翘，堵塞了青干河，形成近20 m高的滑坡坝（图2）。根据两岸植被被流水冲刷痕迹的最大高度量测，滑坡体处涌浪高度为39 m；下游1.2 km处的青干河大桥附近涌浪高度7.4 m；下游1.6 km处的锣鼓洞河口涌浪高度6 m，在该河口的锣鼓洞口上游1 km处仍有渔船被翻覆[4]。

3 涌浪源数值模型的建立

当前，滑坡涌浪研究的方法主要有三种，包括公式法、物理实验法、数值分析法[15,17-24]。由于所需花费时间少、经济和可视化强，数值分析方法被广泛应用。根据流体数值控制方程的不同，当前主流的数值分析方法可分为N-S方程类和水波动力方程类[25]。N-S方程能够很精细地计算滑坡涌浪产生区及邻区的涌浪危害问题；但在面对长距离大面积的计算域时，计算耗时长，计算资源要求很难满足。水波动力学方程则较好地平衡了计算资源与计算精度的关系，具有消耗计算资源少、计算速度快、计算结果精度较好等优点[26]。但水波动力学在面对各类型滑坡涌浪计算时，需要合适的滑坡涌浪源模型来计算涌浪的产生。

从入水区河谷与滑坡方量相比较而言，千将坪滑坡产生的涌浪属于狭窄河道浅水区内产生的涌浪问题。如果滑坡体体积Vs大于水体体积Vw，就可以将该滑坡定义为浅水区滑坡，反之就是深水区滑坡。与深水区产生的涌浪相比，浅水区滑坡涌浪有着比较鲜明的特点。由于相对体积较大，滑坡进入河谷后，一般会堵塞河道，甚至形成堰塞坝。河道几何特征强烈改变，滑坡体对水体的体积侵占必须考虑；因此涌浪在河道上下游的传播被隔断或部分分割[27-29]。

图2 千将坪滑坡全貌（摄于2003年7月15日）Fig.2 The total scene ofQianjiangpinglandslide(Taken in 15th Jul 2003)

因此，以往建立的深水区滑坡涌浪源模型[19]并不适合用于计算浅水区滑坡涌浪的形成。本次研究选择在以往建立的库区滑坡涌浪灾害快速评价系统程序FAST[26]中增加嵌入浅水滑坡涌浪源模型，并将FAST升级为2.0版本。

此次更新的FAST V2.0相比V1.0有了以下功能性改进：1、增加了浅水区滑坡涌浪源数值模。2、将输入地形和遥感影像步骤化，可直接生成所需的地形网格离散数据。数据格式增加国内外工程通用的ASCII格式和DXF格式。3、增加了三维图像自检环节，如出现地形孔洞，可自动插值修补，增强了操作的可持续性，极大的提高了效率。4、增加了涌浪波及范围的自动切取功能，快速勾勒河道一定影响范围区域。5、增加了一键报告功能，能够一键生成简单报告，并能通过邮箱、短信进行发送。

为了产生一个新的滑坡涌浪源，就需要计算形成一个新的波浪液面场作为初始状态提供给传播模型调用。Watts P等人[18]经过物理试验和数值试验认为滑坡在水下的作用时间t0是涌浪的产生时间。同时，在涌浪产生期间，扰动的水体波能以势能为主，扰动区外水体的动能较小。只有在涌浪产生之后的传播阶段，外围水体的动能、势能才开始变得显著。因此，在t=t0时刻，只需关注扰动区的自由液面。由于波形和波速与滑块完全淹没产生的涌浪类型一致，因此可利用FORTRAN语言，编译浅水区滑坡涌浪源模型。由波高函数和波长函数可形成初始的波浪液面势能场，波速用来建立初始涌浪波的运动场。Lamb H[30]通过大量浅水区滑坡涌浪物理模型试验，推导形成了浅水区滑坡初始涌浪源的主要控制方程（等式1-2）。等式3为经典的非线性波波速计算公式[31]。

式中，为相对速度为相对宽度，为相对厚度，lsw为滑块入水长度为滑动冲击角度正弦值，为淹没度，C为波速（m/s），λ为波长（m），a为波幅（m）。

浅水区滑坡初始涌浪源的主要控制方程可以在数学上再现二维的数值初始涌浪波场。而三维的数值初始涌浪波场可通过由滑坡的宽度、滑坡位置和河道地形控制的拓展函数进行重建而来，形成一个圆锥状的三维初始涌浪波场[18-19]。

计算时假定河面为静止状态，涌浪波的能量完全为岩土体能量传递而来，且完全由初始涌浪源控制。在计算时一般假定初始涌浪波的形成时间t0就是崩滑体停止运动的时间ts。当t=t0时，河面形成初始涌浪波，其波高为a，波长为λ，具有一定的势能和动能。此后，滑坡涌浪开始进入河道传播阶段。

FAST的浅水区滑坡涌浪源模型专门形成浅水区的涌浪波，它提供初始涌浪波液面、初始涌浪波速度场、地形文件和水文文件。FAST计算滑坡涌浪源文件后，驱动Boussinesq方程进行涌浪传播和爬高计算，同时提供驱动Boussinesq方程的无缝对接设置和程序服务。基于水波动力学Boussinesq方程,准确的涌浪源输入和合理的软件使用可以产生良好的模拟结果，大量的案例研究证实了一点[18-19,31]。因此在利用FAST平台构建浅水滑坡涌浪模型的关键点在于在涌浪源模块中建立相应的、合理的涌浪源模型，前处理、传播和淹没计算和后处理皆可由FAST的不同模块协同解决。

浅水区滑坡涌浪除了冲击产生的涌浪效应外，还有体积侵占效应不能忽视。由于河道水体浅，造成岩土体入水侵占了河道水体的大量体积。这种侵占量大，不能被忽视，浅水滑坡造成的体积侵占效应要考虑。在数值计算中，这种侵占效应可以采用滑动后的地形来完美的表现。亦即，在新建的浅水滑坡涌浪模型中采用的计算地形网格为滑动后的河道地形，利用这一方法完美解决了河道体积侵占效应。

4 滑坡涌浪分析与验证

利用千将坪作为浅水滑坡涌浪案例来验证新建浅水滑坡涌浪模型的有效性，同时分析千将坪滑坡涌浪演进全过程和衰减规律。

涌浪计算区的地形数据主要采用了1∶10000地形图和水下地形数据。图3展示了计算区地形地貌三维情形，计算区内有长江、青干河、归州河等长江干流和支流，地形复杂多变。计算域东西长约18.6 km，南北长约11.8 km。千将坪滑坡破坏时，库水位为135 m，计算复演当时的滑坡涌浪。采用26 m×26 m的栅格将计算域划分为716列、456行。千将坪滑坡失稳后的滑坡堰塞地形在涌浪计算中被当成固定河岸处理。

图3 计算区地貌三维示意图Fig.3 3D sketch map of calculation region

图4 千将坪滑坡多源计算值对比图[27,30]Fig.4 Comparison map for multisource data of Qianjiangping landslide

千将坪滑坡滑动产生初始涌浪的输入参数见表1。根据多个水库滑坡涌浪数值模拟的计算经验，设置计算8000时步、每20步记录一次的模拟方式。计算每一时步为0.25 s，模拟了2000 s的涌浪过程。

表1 浅水滑坡涌浪源主要输入参数Table 1 The main input parameters for source model for impulse wave in shallow water

千将坪滑坡涌浪的初始特征波幅为34.1 m，波长为314 m。该初始特征波形成后，随即进入传播爬高计算。图4显示了野外观测值、以往利用全N-S方程的计算值[27]和本次Boussinesq计算值有着高度相似，N-S方程计算值与本次的计算值相关性达到0.97，数值平均相差13%。浅水区滑坡涌浪模型的计算值与野外调查相差约为±2m，非常接近。浅水滑坡涌浪源模型再现了支流锣鼓洞河内也存在较大的涌浪灾害，在深入支流1.6 km的涌浪最大约6.7 m。这和当时有船只在该支流内翻沉的现象吻合。特别是，再现了长江对岸的泄滩河和归州河也受滑坡涌浪影响，河道局部有近1 m左右的浪高，这也与实际目击者反映吻合。

与Yin（2015）的N-S方程全耦合计算的涌浪产生区及邻区结果对比显示结果吻合度总体较好，说明浅水区滑坡涌浪源模型在涌浪计算的水动力方面有较好的精度[27]。从计算所需资源来看，N-S方程计算区域仅2.4km×1.68 km，计算文件有46.6GB；而本次水波动力学计算范围是18.6 km×11.8 km文件仅1.9G。因此，在目前计算资源条件下浅水滑坡涌浪源模型可以更高效的解决长距离大范围的流域涌浪灾害问题。

图5从（1）到（4）的顺序展示了千将坪滑坡失稳后涌浪的产生、发展和过程。从滑坡涌浪瞬时水面图来看，千将坪滑坡失稳后，该滑坡形成的涌浪是以滑坡体入水处为源点迅速向四周传播的（图5（1））。40 s左右进入青干河支流锣鼓洞河口，河口附近最大浪高4.6 m（图5（2））。但涌浪向上游逆流传播，至大岭斜坡附近达到锣鼓洞河内最大涌浪高度，达6.7 m。逆流至青干河上游巴木场附近最大涌浪高度仍有3.5 m。100 s左右，涌浪传递至长江河口，长江与青干河河口最大达到1.7 m。145 s左右，涌浪传递至长江对岸的泄滩河，泄滩河最大涌浪达到1 m左右，位置在沟头（图5（3））。180 s左右传递至下游的归州河河口，涌浪进入归州河内进行传递；归州河内最大涌浪高度0.7 m左右。模拟终点时间T=2000 s时，河道中最大浪高为0.24 m，显示大的涌浪波作用过程在30分钟内基本结束（图5（4））。

图5 涌浪波的河道传播过程Fig.5 Propagation process ofimpulse wave in waterway

涌浪最大浪高图是根据模拟时间段内各点的最大波高值形成的（图6）。千将坪滑坡产生的最大浪高为38.8 m，最大爬高为36.7 m。最大涌浪值空间分布显示，涌浪在河道传播中地形及水体对涌浪衰减放大有影响。在河口区域，涌浪急剧衰减。在小河向上游狭窄河道和水体较浅区域传播中有放大效应。大于1m涌浪的河道长约10.5 km，大于2 m涌浪的河道长约9.3 km，大于3 m涌浪的河道长约8 km。1 m以上涌浪大部分分布在锣鼓洞支流和青干河内，仅有零星分布在长江干流及其他支流内。

图6 千将坪滑坡最大涌浪值分布图Fig.6 The maximum impulse wave distribution map of Qianjiangping and slide

通过泓深线涌浪最大浪高图（图7）可知，在主航道上即河流横剖面方向上，涌浪分布区域可分为急剧衰减区和平缓衰减区。急剧衰减区的涌浪波高值呈断崖式下降，大致分布于涌浪源上、下游各400 m长的范围内，该区域平均100 m内涌浪下降高度约为4 m。急剧衰减区是涌浪形成危害的重点区域，一般是滑坡涌浪激发地的附近水域。平缓衰减区内涌浪波高呈缓斜线形式下降，平均100 m内涌浪下降高度只有0.1～0.2 m，平缓衰减区内的浪高表现为起伏型下降，该区域是滑坡涌浪危害的拓展区域，长度较长。

图7 干流及支流的最大波高图Fig.7 The map of maximum wave amplitude in mainstream and branch

5 结论及建议

本文构建了浅水区滑坡涌浪源数值模型，并对千将坪滑坡造成的大面积长距离的涌浪危害进行了研究，得到了以下结论及建议：

（1）千将坪滑坡产生的涌浪属于狭窄河道浅水区内产生的涌浪问题。基于浅水区滑坡涌浪物理实验的推导公式，更新了FAST涌浪计算系统，构建了相应的涌浪源数值模型。

（2）与野外涌浪爬高值、N-S全耦合计算值对比显示，多源计算值与本次的计算值相关性达到0.97，数值相差约为±2 m。浅水滑坡涌浪源模型计算结果有较好的精度。在目前计算资源条件下浅水滑坡涌浪源模型可以更高效的解决长距离大范围的流域涌浪灾害问题。

（3）模拟结果显示，千将坪滑坡产生的涌浪最大浪高为38.8 m，最大爬高为36.7 m。最大涌浪值空间分布显示，大于1m涌浪的河道长约10.5 km，大于2 m涌浪的河道长约9.3 km，大于3 m涌浪的河道长约8 km。大于1m的波浪影响范围不仅包括较长距离的青干河和锣鼓洞河，还包括部分长江干流和其他支流零星河道。

（4）急剧衰减区的涌浪波高值呈断崖式下降，该区域平均100 m内涌浪下降高度约为4 m；平缓衰减区内涌浪波高呈缓斜线形式下降，平均100 m内涌浪下降高度只有0.1～0.2 m。涌浪在河道传播中地形及水体对涌浪衰减放大有影响。在河口区域，一般为快速衰减效应；在小河向上游狭窄河道和水体较浅区域传播中有放大效应。

（5）基于水波动力学的浅水区滑坡涌浪源模型可用于水库区支流或干流浅水区滑坡产生的涌浪问题研究，为大范围长距离滑坡涌浪灾害预测预警提供重要技术支撑。建议进一步开展相关案例对比或物理实验对比研究，为推广该模型奠定基础。

[1]Wang F W,Zhang Y,Huo Z,Peng X,Wang S,Yamasaki S.Mechanism for the rapid motion of the Qianjiangping landslide during reactivation by the first impoundment of the Three Gorges Dam reservoir,China[J].Landslides,2008,5:379-386.

[2]罗先启，肖诗荣，王世梅,等.三峡库区秭归县千将坪滑坡工程地质勘查报告[R].三峡大学，2005.

[3]Jiang Q,Zhang Z,Wei W,Xie N,Zhou C.Research on triggering mechanism and kinematic process of Qianjiangping landslide[J].Disaster Advances,2012,5(4):631-636

[4]李会中，王团乐，孙立华,段伟锋.三峡库区千将坪滑坡地质特征与成因机制分析 [J].岩土力学，2006，（S2）：1239-1244.

[5]廖秋林，李晓，李守定.三峡库区千将坪滑坡的发生、地质地貌特征、成因及滑坡判据研究[J].岩石力学与工程学报，2005，（24）17：3146-3153.

[6]李守定，李晓，刘艳辉.千将坪滑坡滑带地质演化过程研究[J].水文地质工程地质，2008，2：18-23.

[7]王德咏，李科.库水位上升对千将坪滑坡的影响研究[J].灾害与防治工程,2007,2:15-19.

[8]王治华，杨日红.三峡水库区千将坪滑坡活动性质及运动特征[J].中国地质灾害与防治学报,2005,16(3)：5-11.

[9]文宝萍，申健，谭建民.水在千将坪滑坡中的作用机理[J].水文地质工程地质，2008，35(3)：12-18.

[10]邬爱清，丁秀丽，李会中,陈胜宏，石根华.非连续变形分析方法模拟千将坪滑坡启动与滑坡全过程[J].岩石力学与工程学报，2006，25（7）：1297-1303.

[11]肖诗荣，刘德富，姜福兴，姜小兰.三峡库区千将坪滑坡地质力学模型试验研究[J].岩石力学与工程学报，2010，29（5）：1023-1030.

[12]殷跃平，彭轩明.三峡库区千将坪滑坡失稳探讨[J].水文地质工程地质，2007，34(3)：51-54.

[13]张业明，刘广润，常宏,黄波林，潘伟.三峡库区千将坪滑坡构造解析及启示[J].人民长江，2004，35（9）：24-26.

[14]Zhao Y,Cui P,Hu L,Tomasz H.Multi-scale chemo-mechanical analysis of the slip surface of landslides in the Three Gorges,China[J].Science China:Technological Sciences,2011 54(7):1757-1765.

[15]殷坤龙，刘艺梁,汪洋,姜治兵.三峡水库库岸滑坡涌浪物理模型试验[J].地球科学——中国地质大学学报，37（5）：1067-1074.

[16]吴剑，张振华，刘依松.千将坪滑坡中几个高速滑坡问题的探讨 [J].三峡大学学报 (自然科学版)，2007，29（2）：120-123.

[17]Huber A,Hager W H.Forecasting impulse waves in reservoirs[C].//In:Dixneuvie`me Congress des Grands Barrages C,Florence,Italy.Commission International des Grands Barrages,Paris,1997,31:993-1005.

[18]Watts P,Grilli S T,KirbyJ T,Fryer G J,Tappin D R.Landslide tsunami case studies using a Boussinesq model and a fully nonlinear tsunami generation model[J].Natural Hazards Earth SystemScience,2005,3(5):391-402.

[19]Walder J S,Watts P,Oscar E S,Kenneth J.Water waves generated by subaerial mass flows[J].Journal of Geophysical Research,2003,108:2236-2255.

[20]FritzH M,Hager W H,Minor H E.Near field characteristics of landslide generated impulse waves[J].Journal of Waterway Port Coastal and Ocean Engineering,2004,130:287-302.

[21]Ataie-Ashtiani B,Najafi-Jilani A.Laboratory investigations on impulsive waves caused by underwater landslide[J].Coastal Engineering,2008,55(12):989-1004.

[22]Ataie-Ashtiani B,Shobeyri G.Numerical simulation of landslide impulsive waves by incompressible smoothed particle hydrodynamics[J].International Journal for Numerical Methods in Fluids,2008,56:209-232.

[23]Di Risio M,Sammarco P.Analytical modeling of landslide-generated impulse waves[J].Journal of Waterway Port Coastal and Ocean Engineering,2008,134(1):53-60.

[24]汪洋,殷坤龙,杜鹃,刘艺梁.水库滑坡运动与涌浪灾害研究[M].武汉:中国地质大学出版社，2014.

[25]Yavari-Ramshe S,Ataie-Ashtiani B.Numerical modeling of subaerialand submarine landslide-generated tsunami waves-recent advances and future challenges[J].Landslides,2016,13:1325-1368.

[26]黄波林,殷跃平,王世昌,刘广宁,陈小婷.GIS技术支持下的滑坡涌浪灾害分析研究 [J].岩石力学与工程学报,2013,32(S2):3844-3851.

[27]Yin Y P,Huang B L,Chen X T,Liu G N,Wang S C.Numerical analysis on wave generated by the Qianjiangping landslide in Three Gorges Reservoir,China[J].Landslides,2015,12:355-364.

[28]Silvia B,Marco P.Shallow water numerical model of the wave generated by the Vajont landslide[J].Environmental Modelling&Software,2011,26:406-418.

[29]Zhao T,Utili S,Crosta G B.Rockslide and impulse wave modelling in the Vajont Reservoir by DEM-CFD Analyses[J].Rock Mechanics&Rock Engineering,2015,49(6):2437-2456.

[30]Lamb H.Hydrodynamics[M].NewYork:Dover Publications,1945.

[31]Huang B L,Wang S C,Zhao Y B.Impulse waves in reservoirs generated by landslides into shallowwater[J].Coastal Engineering,2017,123:52-61.