某水电站引水隧洞突水数值模拟

徐栋栋，邬爱清，孙玉杰

（1.长江科学院水利部岩土力学重点实验室，武汉 430010；2.胜利油田胜利工程设计咨询有限责任公司，山东东营 257026）

某水电站引水隧洞突水数值模拟

徐栋栋1，邬爱清1，孙玉杰2

（1.长江科学院水利部岩土力学重点实验室，武汉 430010；2.胜利油田胜利工程设计咨询有限责任公司，山东东营 257026）

离散元软件UDEC可以用来模拟裂隙岩体的开挖以及进行水力全耦合分析。采用UDEC来模拟裂隙岩体开挖后在水力耦合作用下渗流流量与其对应水压力的变化过程并预测可能发生的突水灾害。结果表明开挖洞室以后，在围岩渗流与应力耦合作用下，围岩中裂隙隙宽、裂隙中水压及其渗透流量三者相互作用、相互依赖。裂隙隙宽的减小使得结构面水力梯度变大，作用在裂缝上的渗透压力增大，促进导水裂缝扩展，裂隙连通性增加。裂缝隙宽增大，渗透能力增强，渗流量增大，其渗流压力相应降低。在一定条件下，裂隙隙宽的改变可导致局部水力通道的形成，高压水头从局部涌出，从而促进突水灾害的形成。西南某水电站在深部裂隙岩体中开挖引水隧洞，该处地应力高且外水压大，容易引起突水灾害，对其进行了突水数值模拟，提示了一些可能发生突水的位置。

裂隙岩体；引水隧洞；突水；UDEC

1 概 述

西南某水电站引水隧洞具有埋深大、洞线长、洞径大的特点，是该水电站枢纽最重要的组成部分。隧洞所在区域内地质条件复杂、地下水活跃，高压地下水是隧洞开挖及运营期间将面临的关键技术问题之一。深部岩体渗透特性主要取决于岩体中的裂隙网络形态。隧洞开挖前，岩体中的地下水与围岩应力处于一种相对平衡状态。由于隧洞的开挖，一方面使地下水排泄有了新的通道，加速了水循环，破坏了原有的补给-运移-排泄系统的平衡；另一方面，造成围岩应力重新分布，部分结构面由于增压闭合，其岩体卸荷松弛或产生剪切滑移，人为破坏了原有的地下水渗流条件，使得隧洞自身成为地下水以不同形式（渗出、滴流、股流及大范围突水等）向外排泄的地下廊道，有可能形成突水灾害。

突水产生是由于开挖以后在岩体渗流场与应力场的耦合作用下，裂隙的隙宽发生了改变，在一定的条件下隙宽的改变导致局部水力通道的形成，高压水从局部涌出，从而促进了突水灾害的形成。关于岩体渗流场与应力场的耦合研究，主要是连续介质条件下的研究成果［1-8］。这里的连续介质包括两个含义，其一是，将分析对象本身概括为连续介质，采用线性与非线性数值计算方法，研究岩体的渗流与应力耦合问题；其二，鉴于岩体裂隙对岩体渗透特性的重要影响，采用基于连续介质模型的数值方法，通过嵌入岩体中裂隙的本构模型，以实现裂隙岩体的渗流与应力耦合。在这些研究成果中，关于因应力条件的改变，引起岩体中裂隙的张开、闭合与错动等裂隙几何非线性问题，未能考虑。

为在裂隙岩体渗流分析中考虑裂隙张开与闭合变形等的几何非线性问题的影响，进而研究真正意义上裂隙岩体的渗流与应力耦合机制，这里以西南某水电站深部裂隙岩体中的洞室开挖为例，采用UDEC离散单元法中关于裂隙岩体开挖模拟及水力全耦合分析模型，对裂隙岩体洞室开挖进行模拟，判定其是否突水，从而为突水的预测提供一种方法。

2 基于UDEC的渗流分析算法［9］

2.1 基本模型

利用UDEC（Universal Distinct Element Code，ITASCA，1996）水力全耦合模型来模拟裂隙中水的流动。UDEC是ITASCA公司开发的针对非连续介质的平面离散元程序，在数学求解方式上采用了与FLAC一致的有限差分法。对于水力全耦合的渗流分析，当水流主要是由裂隙网络控制时UDEC程序是非常合适的。因此本文假设岩石基质是不透水的［1］，水流主要是通过水力连通的裂隙网络涌入隧洞内。UDEC中被裂隙所包围的岩块可以被模拟为刚体或者可变形体，通过域分析流体在裂隙中的流动。图1中将域顺序标号为①～⑤，假定域内充满各向等压流体，域和域之间通过接触与临域发生作用。接触顺序标号为A-F。域①、③、④表示节理，域②表示2个节理的交点，域⑤为空洞。

图1 通过域模拟流体在节理裂隙中的流动模型Fig.1 Hydraulic simulation for a fissured rock

不计重力时，假设流体压力在流动域中的分布是均匀的。考虑重力时，流体压力则按线性分布的静水压力计算。流体的流动是由相邻流体域的压力差决定的，其中流动域中流体压力的大小由流动域中中心压力的大小决定。

2.2 渗透流速计算

按块体接触条件的不同，裂隙岩体中流体的速率有2种计算方法。

（1）点接触。点接触分为角-边接触和角-角接触。设流动域①的流体压力为p1，流动域②的流体压力为p2，则由流动域①到流动域②流体的流速为

式中kc为接触处的渗透系数；

式中：ρw为流体密度；g为重力加速度；y1，y2为2个流动域的中心坐标。

（2）边-边接触。首先定义接触长度，图1中lD和lE分别为D和E的接触长度，然后运用平行板裂隙中的立方定律计算流动速度，即

水力开度a由裂隙在无法向应力时的开度a0及在某法向应力条件下法向开度增量un组成，即：a＝a0+un。假定法向开度增量un张开为正，压缩为负，水力开度的最小值为ares，最大值为amax。裂隙宽度随节理法向应力的变化如图2所示。

2.3 裂隙水压力计算

计算过程中，每计算一个时步，重新生成系统的几何形状，而后计算出所有接触的裂隙宽度以及所有域的体积（对二维条件，取单位厚度），之后利用上面的公式计算出各接触处的流量。最后，迭加各接触点流入裂隙域的流体流量，并考虑由于周围块体的位移增量而产生的域体积的变化，按下式计算出域内的裂隙水压力：

式中：p0为前一时步的孔隙压力；Q为通过孔隙周围的所有接触点流入该孔隙的流量之和；kw为流体的体积模量；Δv＝v-v0，vm＝（v+v0）／2，其中，v和v0分别为现在时步和前一时步孔隙的体积；Δt计算时步。

计算出域内裂隙水压力后，可以计算流体作用在其周围岩块的力。将该力与诸如接触点力和外力荷载等力迭加，施加在块体的节点上。这样得到不透水岩块的总应力以及节理的有效法向应力。

图2 水力开度a与法向应力σn之间的关系Fig.2 The relationship between hydraulic opening（a）and normal stress（σn）

3 突水数值模拟

3.1 工程概况

该水电站引水隧洞贯穿属于过山隧洞，长度约16.67 km，4条引水隧洞平行布置，引水隧洞之间的中心轴线距设计为60 m，约为隧洞开挖洞径的4.6倍，隧洞之间的净岩体厚度为47 m，约为隧洞开挖洞径的3.6倍，开挖洞径13 m，全线埋深较大，一般埋深1 500～2 000 m，最大埋深约2 525 m，属于深埋长隧洞。引水隧洞所在区域地形地质条件复杂、地下水活跃，受大气降水补给，地下水极为丰富。

3.2 计算模型及参数

根据引水隧洞线的工程地质条件评价，在第3段硐深1 220 m处选取断面进行数值模拟工作。由于4条引水隧洞之间的中心间距仅为60 m，且4条引水隧洞横断面尺寸完全相同，该部分模拟只选取西端的3号和4号隧洞进行模拟。目前，4条引水隧洞只完成了进水口段以及辅助洞的开挖，交通洞正在开挖，因此在已开挖区段的岩体中，尚难获取足够的结构面统计样本。因为辅助洞与引水隧洞平行布置，二者处于同一高程线上，研究分析辅助洞的突水机理、规律，对于引水隧洞突水的预测预报和防治处理无疑具有重要意义。因此在本次模拟中，以辅助洞内的结构面量测为依据来进行建模。

3.2.1 几何模型

此处发育有Ⅱ级结构面断层F28以及Ⅲ-2级结构面断层fw18以及Ⅳ级结构面，根据地质资料中断层的产状将其投影到所选取的断面上，得到裂隙网络。其中，根据量测数据得到Ⅳ级结构面的统计参数如表1所示。

表1 某水电站辅助洞Ⅳ级结构面统计参数Tab le 1 Statistical parameters of forth-level joint group of the auxiliary tunnel of some hydropower station

模型边界的结果选定分析区域为177 m×104 m。依据断层的投影以及表1中的Ⅳ级结构面的统计参数采用Monte-Carlo方法生成图3所示的结构面网络样本。依据结构面的交切关系，生成的连通水力网络如图4所示。

图3 结构面网络样本Fig.3 The joint network modeling result

图4 连通水力网络Fig.4 The connective hydraulic network

3.2.2 模型材料参数

根据该水电站的地质资料，该处断面所采用的岩块以及结构面所采用的模型材料参数见表2和表3所示。流体的密度、动力粘滞系数和体积模量分别为1 000 kg／m3，3.5×10-4Pa·sec，2×103MPa。

表2 岩块参数Table2 The parameters of the rock blocks

表3 结构面参数Table3 The parameters of structural p lanes

3.3 初始条件及边界条件

（1）初始条件。根据引水隧洞现场地应力的测试成果，引水隧洞初始应力场以自重应力场为主，中间主应力平行于河谷方向（即近似垂直于隧洞轴线方向），其侧压力系数为1.0。由于所建模型上表面距地面1 080 m，故上表面施加29.376 MPa的压应力，上表面初始水头为230 m，底面初始水头为334 m，左右两侧面施加沿重力方向梯度变化的水头压力。

（2）边界条件。由于该处模型边界的选取是根据分析得到的，因此可以认为假想边界处距地下洞室无限远，因而可以在假想边界处采用位移边界条件限制位移发展。

3.4 突水数值模拟

考虑到模拟突水情况的实际需要，这里采用UDEC中流体不可压缩的瞬态分析算法。

突水现象最为直观的就是流量及其对应水压力的变化情况，因此选取了3号以及4号隧洞中比较有代表意义的几个点如表4所示（这几点都位于开挖面与断层交接处，更容易发生突水），用这几个点处的流量以及水压力变化情况来模拟是否有突水现象的发生。

表4 监测点坐标Table4 The coordinates of themonitoring points

3.4.1 4号隧洞突水模拟结果

监测点1和2处的流量和水压力变化曲线见图5和图6。

从图5（a）和图6（a）可以看出，位于4号隧洞左边墙中部监测点1处（隧洞左边墙与Ⅱ级结构面交切点）开挖后由于揭露了渗流通道发生直涌型突水一次，在渗流时间大约为1 s处最大流量达到630 m3／s，与该点相对应的水压力也变化到最大。但随着时间流量逐渐减小并很快稳定。在流量达到稳定前那一时刻（大约2 s处）可以明显地看出水压力发生突变，从5.2 MPa急剧地减小到2.2 MPa左右，而后也很快稳定。这是因为该点处有块体发生脱离导致水压力急剧减小，但并未形成完整的水流通道，而是比较稳定的渗流通道，水压力和流量都为恒量。见图7监测点1处块体分离。

图5 监测点1和2处流量变化曲线图Fig.5 Variations of flow rate at themonitoring point 1，2

图6 监测点1和2处水压力变化曲线图Fig.6 Variations of the pore pressure of themonitoring point1，2

图7 监测点1和2处块体位置Fig.7 The positions ofmonioring point 1 and 2 in blocks

从图5（b）和图6（b）可见位于4号隧洞顶拱右侧的监测点2处，开挖后由于揭露了渗流通道发生直涌型突水，最大流量达到185m3／s，对应点水压力也达到最大，其后流量基本上在40 m3／s上下波动，水压力也在2.5～3.7 MPa范围内波动。因为水力耦合作用一直在进行但还没达到最终的稳定，所以流量和水压力都在一个相对稳定的范围内波动。同时在2 s处同样也发生了从5.2 MPa到2.2 MPa的突变。这里突变的原因可能是块体在水力耦合作用下发生了大的变形从而引起隙宽的改变而致，而尚未产生块体的脱离，以后很有可能会分离。见图7监测点2处块体有要分离的趋势。

3.4.2 3号隧洞突水模拟结果

监测点3和4处的流量和水压力变化曲线见图8和图9。从图8（a）和9（a）可以看出，位于3号隧洞左边墙上部点3处开挖后由于揭露了渗流通道发生直涌型突水，6 s之后流量都在20 m3／s以内波动；对应的水压力在6 s以内波动比较大，随后波动范围减小。在水力耦合作用力作用下裂隙隙宽一直在调整，流量和水压力变化比较大是因为在这种作用下块体发生了大的变形导致隙宽改变较大，随着时间趋于相对稳定。

由图8（b）和9（b）可以看出，位于3号隧洞底部监测点4处，在隧洞开挖后由于揭露了渗流通道发生了突水，但流量在逐渐减小，而在6 s，15 s，19 s时发生3次突水。在水力耦合作用力作用下裂隙隙宽在调整，在渗流到6 s，15 s，19 s时水压力梯度增大导致裂隙隙宽变大，发生水力劈裂，因此对应流量加大，水压力开始减小，之后裂隙隙宽又在该作用力作用下开始减小，流量也随之减小。此后监测点3处在水力耦合作用力作用下，裂隙隙宽不断调整，有可能继续发生突水。图10可见监测点4处劈裂痕迹。

图8 监测点3和4处流量变化曲线图Fig.8 Variations of flow rate at themonitoring point 3，4

图9 监测点3和4处水压力变化曲线图Fig.9 Variations of the pore pressure at themonitoring point 3，4

图10 监测点4处块体位置Fig.10 The position of point 4 in blocks

4 结 论

3号和4号隧洞在开挖时揭露了2条岩溶断层F28以及fw18，模型中隧洞开挖后应力重新分布，临空面微裂隙在高水头地下水压力作用下张开、劈裂，因此理论上分析在3号和4号隧洞内会有纯劈裂型以及直涌型两种类型突水情况发生，数值模拟的结果也印证了该分析的正确性。

该水电站引水隧洞工程地质条件复杂，工程区位于高山峡谷岩溶区，隧洞较长，地下水埋藏深，露头少，且存在大于600～900 m的高水头压力，这使得对于突水预测预报变得十分复杂和困难。本文的研究成果对解决这个难题提供了可能性，可以根据地质资料以及现场的结构面量测选取引水隧洞揭露岩溶断层或裂隙的断面建立模型，采用本文的方法，借助于UDEC软件开展突水的数值模拟，找出可能发生突水的位置以及突水的强度。

［1］ 张有天.岩石水力学与工程［M］.北京：中国水利水电出版社，2005.（ZHANG You-tian.Rock Hydraulics and Engineering［M］.Beijing：Water Conservancy and Hydropower Press in China，2005.（in Chinese））

［2］ 忤彦卿.岩体裂隙系统渗流场与应力场耦合模型［J］.地质灾害与环境保护，1996，7（1）：31-34（WU Yanqing.The coupled stress-groundwater flow model for jointed rock mass［J］.Journal of Geological Hazards and Environment Preservation，1996，7（1）：31-34.（in Chinese））

［3］ 陈 平，张有天.裂隙岩体渗流与应力耦合分析［J］.岩石力学与工程学报，1994，13（4）：299-308.（CHEN Ping，ZHANG You-tian.Coupling analysis of seepage／stress for jointed rock［J］.Chinese Journal of Rock Mechanics and Engineering，1994，13（4）：299-308.（in Chinese））

［4］ 王 媛，徐志英，速宝玉.裂隙岩体渗流与应力耦合分析的四自由度全耦合法［J］.水利学报，1998，（7）：55-59.（WANG Yuan，XU Zhi-ying，SU Bao-yu.Fourfreedom completemethod for the seepage-stress coupled analysis in fissured rock masses［J］.Journal of Hydraulic Engineering，1998，（7）：55-59.（in Chinese））

［5］ 陈洪凯.裂隙岩体渗流研究现状（Ⅰ）（Ⅱ）［J］.重庆交通学院院报，1996，15（2）：29-34.（CHEN Hongkai.Review on seepage researches of fissured rock mass［J］.Journal of Chongqing Jiao Tong University，1996，15（2）：29-34.（in Chinese））

［6］ 王 媛，速宝玉，徐志英.等效连续裂隙岩体渗流与应力全耦合分析［J］.河海大学学报，1998，26（2）：26-30.（WANG Yuan，SU Bao-yu，XU Zhi-ying.Coupling analysis of seepage and stress inmulti-fractured rockmasses［J］.Journal of Hohai University，1998，26（2）：26 -30.（in Chinese））

［7］ Ayatollahi M S.Stress-fluid flow analysis of fractured rocks［J］.J.of Eng.Mech.，1983，109（1）：1-13.

［8］ Oda M.An equivalent continuum model for coupled stress and fluid flow analysis in jointed rock masses［J］.Water Resource Research，1986，22（12）：3235-3246.

［9］ ITSCA.UDEC-Universal Distinct Element Code［K］.Minneapdis，USA：ITSCA Consulting Group，1996.

（编辑：曾小汉）

Simulation of Sudden Blow in Diversion Tunnel of Some Hydropower Station

XU Dong-dong1，WU Ai-qing1，SUN Yu-jie2
（1.Changjiang Scientific Research Institute／Key Laboratory of Geotechnical Mechanics and Engineering of the Ministry ofWater Resources Division of Geotechnical Engineering，Wuhan 430010，China；2.Shengli Engineering Consulting Co.Ltd of the Shengli Oil Field，Dongying 257026，China）

The discrete elementmethod can be used to simulate the excavation of the fractured rock and conduct the fully-coupled analysis of hydraulics.So UDEC（universal distinct element code）can be used to simulate the changes of flow rate and the corresponding water pressure under the coupling conditions of the seepage and stress after the excavation of the fractured rock，and forecast the possible sudden blow disasters.The results show that the changes of the fracture aperture，water pressure and fracture permeability discharge of the fractured rock mass in the coupling conditions of the seepage and stress are interactive and interdependent.Because of the shortening of the fracture width，the hydraulic gradient and the permeability pressure of the fracture become large.It also promotes the expansion of cracks and fissure connectivity increases.The increase of the crack width will result in the increase of the infiltration capacity and seepage discharge and reduction of seepage pressure.Under certain conditions，changes in fracture aperture can lead to the formation of localwater-access，the emission of high pressure head from the local place，and thus promotes the formation of the sudden blow disasters.The sudden blow easily occurs after the excavation of the diversion tunnel situated at the deep fractured rock where the ground stress and water pressure are very high in some hydropower station in the southwestof China.The simulation is done in this article，and some possible positions where the sudden blow occur are revealed.

fractured rock mass；diversion tunnel；the sudden blow；UDEC

TU457

A

1001-5485（2010）08-0044-06

2009-08-27；

2009-10-23

徐栋栋（1986-），男，山东聊城人，硕士研究生，主要从事岩石力学与工程方面的研究工作，（电话）15994238327（电子信箱）xdhappy717＠163.com。