考虑渗流及软化效应的深埋引水隧洞围岩弹塑性分析

惠 强， 姜海波， 张玉洁， 翟栋森

(石河子大学 水利建筑工程学院， 新疆 石河子 832003)

1 研究背景

水工引水隧洞是水利枢纽的重要组成部分，随着我国溪洛渡水电站、乌东德水电站、白鹤滩水电站等大型水利工程的陆续建设，水工隧洞埋深逐渐增大，而深部岩体中往往存在丰富的地下水，故水工引水隧洞往往受地应力、渗流场综合作用，极大地影响了隧洞围岩应力分布。因此考虑渗流影响下深埋水工引水隧洞围岩的受力及应力分布的研究成为学术界与工程界关注的热点[1]。

近年来，许多学者在隧洞经典弹塑性解的基础上进行研究，推导出了考虑渗流场影响下的隧洞弹塑性解，文献[2]通过将深埋圆形隧洞简化，运用Mohr-Coulomb强度准则求解得到考虑渗流场的围岩应力解析解，探讨了渗流场对围岩应力的影响；文献[3]在文献[2]研究的基础上引入应力调整系数，求得考虑地应力重分布影响，深埋圆形透水隧洞弹塑性解析解；文献[4]运用弹塑性损伤力学理论，推导出了巷道围岩稳定性的理论解答；文献[5]通过对包括Mohr-Coulomb准则和Hoek-Brown准则在内的4种岩石强度准则进行归纳整理，得到了平面应变条件下围岩屈服方程的统一形式，在此基础上推导了渗流作用下的围岩应力场、位移场的统一解；文献[6]改进了基于Hoek-Brown准则推导出来的围岩弹塑性应力解答，使其不含有积分项。可以看到，目前学者进行围岩弹塑性计算分析时，大多选择采用未考虑中间主应力的Mohr-Coulomb准则[3]和Hoek-Brown准则[6-7]导致计算结果偏于保守[8]，并且在进行渗流场计算时，以往大多采用各向同性进行计算，即各个方向上渗透系数相等，这与工程实际存在差异，只有掌握了隧洞开挖后围岩渗流场的分布规律，才能分析出渗流作用对围岩应力及塑性区范围的影响程度，为相关工程实践提供可靠的理论指导。

本文在前人研究的基础上，基于统一强度理论准则，同时考虑应变软化、渗流及中间主应力的影响，将引水隧洞围岩分为塑性残余区、塑性软化区及弹性区，推导围岩各分区的应力及半径解析表达式，同时考虑水平方向和垂直方向上渗透系数的不同，通过定义不均匀渗透系数，分析其对围岩各方向上孔隙水压力分布规律的影响。最后通过算例分析了围岩软化、渗流作用、中间主应力系数对隧洞切向应力和塑性区半径的影响规律。研究结果可为进一步研究含水岩体下隧洞开挖围岩稳定性等问题提供参考。

2 基本假设及理论模型

2.1 力学模型

为了定性研究考虑渗流及应变软化作用下圆形引水隧洞围岩应力分布情况，本文求解时对实际问题做了如下假定：(1)计算过程按轴对称下平面应变问题考虑。(2)将含水围岩视为满足达西定律的两相介质体。(3)为耦合计算方便，文章在计算时以压应力为正值，拉应力为负值。

图1为隧洞围岩力学模型，根据围岩应力应变状态，将隧洞围岩分为3个区域，即弹性区、塑性软化区和塑性残余区，弹性区的围岩处于完整状态，当围岩应力超过岩体强度极限时，围岩处于塑性软化状态，随着变形的逐渐增大，岩体强度逐渐衰减，最终达到残余强度，此时围岩进入塑性残余区，相应力学参数取值见公式(1)、(2)。隧洞开挖半径为R0(m)，隧洞圆心至塑性残余区外边界距离为Rr(m)，隧洞圆心至塑性软化区外边界距离为Rp(m)；半径R(m)以外的稳定渗流场水压与原始渗流场外水压力Pi(Pa)相同，R可由钻孔试验获得[9]；隧洞支护力为P0(Pa)，初始地应力为σ0(Pa)，侧压力系数λ为1，即水平地应力与垂直地应力相等。

图1 围岩力学模型

2.2 应变软化模型

围岩的应变软化在细观上表现为应力作用下微小裂隙的张开、贯通及围岩性质的弱化；宏观上可以认为围岩的黏聚力、内摩擦角等参数发生改变[10]。围岩力学参数黏聚力C(Pa)、内摩擦角φ(°)分别为η的分段函数[11-12]：

(1)

(2)

图2 最大主应力变化曲线

文献[13]、[14]给出了η*的计算公式如下：

(3)

(4)

式中：ψ为岩石剪胀角，(°);E为围岩弹性模量，Pa；M为软化曲线斜率，Pa，由于岩石应力-应变曲线卸载斜率和加载斜率相差不大，故此处假设两者是相同的[12]。

2.3 统一强度理论

统一强度理论具有多种表达形式[15-16]，对于岩石类材料，采用岩石内摩擦角φ、岩石黏聚力C表示则为[16]：

=2Ccosφ

(5)

(6)

(7)

(8)

式中：σ1、σ2、σ3分别为最大主应力、中间主应力以及最小主应力；b为中间主剪应力以及其作用面上的法向正应力对材料破坏的影响程度，其取值范围为[0,1][16]。

2.4 渗流场计算

在轴对称平面稳定渗流场中，考虑到岩体在水平、垂直方向上具有不同的渗透系数，两个方向上的达西定律表示为[17]：

(9)

式中:V为渗透水流流速,m/d；h为渗流位势,m；kx、ky分别为水平、垂直方向上的渗透系数,m/d。

不可压缩水体在岩体中流动的连续性方程为：

(10)

孔隙水压力Pw=h·γw，γw为水的重度，N/m3，定义水平、垂直方向上的渗透系数比值为不均匀渗透系数ν，即kx=ν·ky，当ν=1时，岩体为各向同性。将公式(9)与(10)联立，采用柱坐标系表示，得轴对称平面稳定渗流微分方程：

cos2β)=0

(11)

(12)

当ν=1时，上式变为各向同性岩体的轴对称平面稳定渗流场，与文献[17]推导公式一致。

由边界条件：

Pw(r=R0)=0,Pw(r=R)=Pi

解得隧洞围岩孔隙水压力沿半径的分布规律为：

(13)

3 考虑渗流及软化效应的围岩弹塑性分析

由于引水隧洞围岩受内水压力、外水压力及地应力的共同作用，在不同工况下，围岩的第一主应力既可能为径向应力，也可能为切向应力[18]。当引水隧洞处于施工期或隧洞处于运行期且初始地应力大于隧洞内水压力时，若初始地应力大于隧洞内水压力，或隧洞处于施工期时，即σ0>P0时，有σθ>σr成立[19]。由于文章研究对象为施工期的深埋引水隧洞，不考虑内水压力的影响，故此时第一主应力为切向应力。对于平面应变问题，中间主应力为[20]：

(14)

式中：m为中间主应力参数，对于岩石类材料，m=1。

即在此情况下，σθ=σ1、σz=σ2、σr=σ3。由于满足公式(7)条件，将公式(1)、(2)、(14)代入公式(7)得第一主应力为切向应力时的围岩统一强度理论公式：

σθ=σrA1+A2

(15)

(16)

(17)

3.1 塑性残余区

考虑渗流作用的平衡方程为[4,17]：

(18)

式中:α为渗透水压力的作用面积系数，为安全起见，在研究岩体的破坏和稳定时，一般取α=1[21]。

将公式(13)、(18)联立求解，受限于文章篇幅，只推导得出ν=1时的应力表达式如公式(19)、(20)所示，ν≠1的应力表达式可用相同步骤推导求出。

(19)

(20)

式中:C为待定常数。将塑性残余区内边界(r=R0)上的边界条件σrpr=Po代入公式(19)，可得围岩塑性残余区应力表达式：

(21)

(22)

3.2 塑性软化区

(23)

(24)

3.3 弹性区

(25)

由切向应力与径向应力之和在围岩弹塑性交界面处(r=Rp)连续可得：

(26)

解得：

(27)

将公式(27)与(24)联立可得：

(28)

将公式(3)、(28)联立可得：

(29)

同理可得无渗透作用下完全由地应力引起的峰值应力为：

(30)

(31)

上标g表示由地应力引起的应力值，将公式(28)、(29)、(30)、(31)代入公式(25)、(21)、(23)可得由软化参数及相关岩石参数表示的考虑渗流及软化效应的围岩各分区应力表达式。

4 算例分析

以新疆某引水隧洞工程某标段为例进行算例分析。根据现场实测资料，该段岩石较坚硬，该段围岩定为Ⅲ类。本段上覆围岩厚180～800 m，洞室内有渗水和滴水现象。隧洞开挖半径R0=2 m，由钻孔实测及抽水试验可知，保持原始渗流场外水压力的围岩影响半径R=10R0=20 m。为了进一步探讨围岩软化、渗流效应、中间主应力系数对隧洞围岩应力和塑性区半径的影响规律，在此以该工程实际围岩力学参数为依据进行算例分析，围岩力学参数见表1。

表1 实例工程围岩力学参数

4.1 不均匀渗透系数对孔隙水压分布的影响

从公式(13)可以看出，影响围岩孔隙水压力分布的因素很多。下面分析不均匀渗透系数ν对围岩不同方向上孔隙水压力分布的影响。

图3为不同ν值条件下不同方向上的围岩孔隙水压力分布，其中横轴为该点距隧洞圆心距离r与原始渗流场外水压力的影响半径R的比值。由图3可知，不均匀渗透系数对不同方向上的孔隙水压力分布皆有影响。图3(a)表明，当外水压力一定时，随着不均匀渗透系数的增大，沿0°方向上(水平方向)孔隙水压力的增长速度逐渐变缓，最终趋于外水压力。沿30°方向上孔隙水压力的变化趋势与0°方向上相似(图3(b))，但是不同ν值下的增长速度差异性没有后者明显；图3(c)和图3(d)显示60°方向上与90°方向(垂直方向)孔隙水压力的变化趋势相同，都是随着不均匀渗透系数的增大，孔隙水压力的增长速度逐渐加快，90°方向上不同ν值下的增长速度差异性较60°方向上更为明显。这说明随着ν值的改变，渗流沿不同方向表现出了各向异性，且越靠近水平方向或垂直方向，孔隙水压力分布越容易受不均均匀渗透系数的影响，具体来说，当ν<1时，越靠近水平方向，孔隙水压力变化速度越快；当ν>1时，越靠近垂直方向，孔隙水压力变化速度越快。这表明，当ν<1时，即水平方向上的渗透系数小于垂直方向上的渗透系数时，在远离隧洞临空面的围岩深处水平方向上的孔隙水压力降低幅度远小于垂直方向上孔隙水压力的降低幅度，由于隧洞临空面处水压力边界条件为一定值，故计算结果显示靠近临空面处水平方向上水力梯度陡然增大。究其原因，笔者认为不同渗透系数对应着不同岩性的地层，其水力影响半径也有所不同，但在本文研究过程中，为简化计算假设其影响半径相同均为R。根据该算例结果可推测，渗透系数小的方向的水力影响半径应该反而更大，由于该部分不是本文重点，故岩层渗透系数与水力影响半径的关系可由后续研究证明。

图3 不同ν值条件下不同方向上的围岩孔隙水压力分布

4.2 渗流及软化效应对围岩不同分区范围及应力分布的影响

为揭示应变软化及渗流作用对围岩应力、塑性区半径的影响，本节将本文第3节弹塑性解(既考虑围岩应变软化也考虑渗流作用)与以下3种情况的弹塑性解进行对比讨论：不考虑应变软化、不考虑渗流作用、既不考虑软化也不考虑渗流作用。为简化计算，将软化区内岩石材料参数的衰减看作线性变化。

若不考虑渗流及软化作用，塑性区满足如下平衡方程：

(32)

求得不考虑渗流作用及软化作用围岩塑性区及弹性区应力表达式为：

(33)

(34)

塑性区半径表达式为：

(35)

当b=0时，统一强度准则变为Mohr-Coulomb强度准则，此时公式(35)称为Kastner方程[22]。

(36)

由弹塑性边界上切向应力与径向应力之和连续可得塑性区半径表达式为：

(37)

若不考虑围岩渗流作用，塑性区同样满足平衡方程(32)，弹性区视为厚壁圆筒，用上文方法可解得围岩各分区的应力表达式，这里不再赘述。最后求得塑性区半径表达式为：

(38)

(39)

当中间主应力系数b=0时，则统一强度准则退化为Mohr-Coulomb强度准则，按上文公式计算求得4种情况下P0与塑性区半径Rp之间关系曲线如图4所示。由图4可知，P0对围岩塑性区半径Rp影响显著。当b一定时，随着P0的增大，塑性区半径Rp逐渐减小直至无塑性区(由于挖半径R0=2 m，故Rp<2 m时表明塑性区不存在)。当P0=4 MPa时，4种情况所对应的Rp比P0=2 MPa时对应的Rp分别减小了28.4%、16.1%、20.8%、13.3%；而当P0=10 MPa时，4种情况所对应的Rp比P0=8 MPa时对应的Rp分别仅减小了10.9%、6.7%、9.3%、5.7%，表明在一定范围内，支护力可以减小围岩塑性区范围，且支护力越小，该影响越显著。当考虑渗流及软化时，所得Rp均大于仅考虑单一因素或两者皆不考虑时的Rp，因此当隧洞施工时，应当通过注浆处理等措施加固岩体，降低围岩软化程度及渗透系数，减小围岩塑性区范围，提高围岩稳定性。且仅考虑软化情况计算所得Rp大于仅考虑渗流情况时的Rp，而这两种情况下的Rp均大于两者皆不考虑时计算的Rp，表明渗流及材料软化都会对围岩稳定造成影响，在计算时不应忽略，且材料软化较渗流影响对围岩稳定更为不利。

图4 4种情况下P0与Rp的关系曲线(b=0)

图5为当P0=2 MPa、b=0时，4种情况下围岩切向应力分布情况，限于篇幅，且第一主应力对围岩一般起控制作用，而在此工况下切向应力为第一主应力，故只分析各参数对切向应力的影响。由图5可知，靠近隧洞临空面处，考虑材料软化情况下的切向应力值小于不考虑材料软化的切向应力值，这是由于塑性区内围岩力学参数发生劣化导致围岩承载能力下降，应力向围岩深处传递，导致塑性区范围不断扩大。同时考虑渗流场影响与不考虑渗流场影响的围岩切向应力分布在临空面附近基本一致，但靠近深处，考虑渗流场影响的围岩切向应力逐渐大于不考虑渗流场作用的围岩切向应力，这是由于不考虑渗流场时，作用在洞壁的支护反力以面力考虑，其作用范围有限，而渗透水压力实际上为体积力，作用于应力场的任意一点。

图5 4种情况下围岩切向应力分布(P0=2MPa, b=0)

图6为考虑渗流及软化条件下，当b=0时，不同P0对应的围岩切向应力分布情况，其塑性区半径见图4。由图6可知，随着P0的增大，临空面附近围岩切向应力逐渐增大，P0从2 MPa增加到10 MPa，则隧洞内壁处切向应力由8.8 MPa增加至30.9 MPa,增加了2.5倍，而峰值切向应力则由35.5 MPa减小至30.9 MPa，减小了13.0 %。

图6 不同P0对应的围岩切向应力分布(考虑渗流及软化，b=0)

4.3 中间主应力系数对围岩不同分区范围及应力分布影响

上文基于统一强度理论准则，综合考虑渗流场及围岩应变软化影响，求得隧洞围岩弹塑性解析解。已知通过改变中间主应力系数b可将该系列解转化为已知解，如b分别等于0、0.25、0.50、0.75、1.00时，该解可转化为满足Mohr-Coulomb、双剪强度准则等强度准则的特解。

通过改变中间主应力系数b分析考虑渗流及软化效应下中间主应力对围岩切向应力的影响，图7为P0=2 MPa时不同b值对应的围岩切向应力分布情况，图8为不同b值对应的塑性区半径大小。由图7、8可知，b值对围岩切向应力及塑性区范围有着显著的影响，随着b值的增大，塑性区半径逐渐减小，而临空面附近围岩切向应力则逐渐增大。考虑中间主应力效应(b=0.25、0.50、0.75、1.00)时围岩峰值切向应力相较不考虑中间主应力效应(b=0)分别增加了0.2%、0.8%、0.9%、1.8%，塑性区半径大小分别减少了11.3%、17.7%、21.8%、24.7%，这表明相较于考虑中间主应力效应，不考虑中间主应力效应时计算结果相对保守，不能有效发挥出围岩的强度。

图7 不同b值对应的围岩切向应力分布(P0=2MPa) 图8 不同b值对应的塑性区半径

5 模型验证与讨论

为验证文章理论模型的准确性，采用文章第3节工程围岩参数，对文章提出的塑性区半径计算结果与文献[5]的塑性半径进行对比分析。表中中间主应力系数b=0时，统一强度准则退化为Mohr-Coulomb强度准则；当b≠0时，Mohr-Coulomb强度准则不再适用，取b=0.25讨论两种理论计算结果对比。结果见表2。

表2 不同Pi值下塑性半径本文理论计算值与文献[5]计算值对比

从表2可以看出，在b=0时,塑性半径本文理论值计算结果与文献[5]计算结果基本一致；而在b≠0时，其结果偏大，且两种方式计算得到的塑性区半径均随外水压力的增大而扩大。总体来说，两者结果较吻合，充分说明了本文理论模型的准确性。

6 结 论

(1)本文基于统一强度理论准则，考虑渗流、应变软化和中间主应力效应影响，将围岩分为塑性残余区、塑性软化区和弹性区，推导出圆形隧道围岩塑性区以及围岩各分区的应力分布表达式，对于确定含水隧洞围岩支护范围及措施等方面提供了一定理论参考。

(2)定义不均匀渗透系数ν来定量分析双向不等渗透系数对围岩各方向上孔隙水压力分布的影响：随着ν值的改变，渗流沿不同方向表现出了各向异性；当ν<1时，越靠近水平方向，孔隙水压力变化速度越快；当ν>1时，越靠近垂直方向，孔隙水压力变化速度越快。

(3)通过算例比较，渗流及材料软化对围岩应力场分布以及塑性区大小均有不同程度的影响：由于围岩力学参数发生劣化，考虑软化较不考虑软化时的切向应力更小，塑性区半径更大，同时材料软化较渗流影响对围岩稳定更为不利。

(4)中间主应力系数b对围岩切向应力及塑性区范围有着显著的影响，随着b值的增大，塑性区半径逐渐减小，而临空面附近围岩切向应力则逐渐增大。在本文算例中，b=1时围岩峰值切向应力相较b=0增加了1.8%，塑性区半径减少了24.68%，这表明相较于考虑中间主应力效应，不考虑中间主应力效应时计算结果相对保守，不能有效发挥出围岩的强度。