考虑流固耦合效应的浅埋矩形顶管隧道开挖稳定性分析

屈克军

(中铁十四局集团大盾构工程有限公司， 南京 211899)

城市地下空间开发中，顶管法因其对周围扰动小、施工过程安全等特点而被广泛应用。然而，矩形顶管隧道在高水压高渗透复合地层、易液化砂层等复杂地质条件施工时，因支护力设置不当而造成的开挖面失稳甚至地表坍塌事故时有发生[1]。因此，地下水渗流条件下复合地层开挖面稳定性分析和极限支护力的获取，具有重要的工程应用价值。

基于数值模拟、模型试验与理论分析等方法，中外学者对地下水渗流条件下掌子面稳定性进行了研究。李廷春等[2]应用流固耦合计算原理，验算了厦门某海底隧道的稳定性。Lee等[3]基于极限分析上限定理求解掌子面极限支护力时，建立了考虑渗流力的计算模型，并与传统极限平衡法结果相对比验证了其合理性。Perazzelli 等[4]基于条分法，改进了Anagnostou等[5]提出的楔形体破坏模型，并推导出渗流影响下的支护力闭合解。黄阜等[6]基于数值模拟获取了盾构开挖面前方的渗流场，并利用上限定理得到了安全系数上限解。高健等[7]基于有限元法，提出了考虑盾构掘进速度的二维渗流方程，与现场实验数据有较好的一致性。曹利强等[8]基于砂土地层渗流条件下孔压的分布特征，在经典楔形体模型基础上引入水压分布函数，推导了极限支护力。黄正荣等[9]利用数值模拟研究了不同水位对开挖面稳定性的影响，发现地下水的存在对开挖面支护力影响较大。此外，吕玺琳等[10]基于自行研发的离心模型试验机，研究了渗流条件下饱和粉砂中盾构开挖面失稳特征，发现地下水渗流对极限支护力有一定的影响。

综上可知，目前研究地下水渗流对掌子面稳定性影响时，大多是将稳态渗流力作为附加力施加到破坏模式的计算中，而较少考虑渗流场和应力场的耦合效应。因此，现依托南通市地铁世纪大道站2号线工程，采用FLAC3D建立考虑流固耦合效应的三维地质模型，将实际地表沉降监测值与数值模拟结果进行对比验证。基于数值模型和理论分析，对渗流影响下的浅埋顶管隧道开挖面稳定性开展研究。

1 工程概况

1.1 工程概述

世纪大道站设有4个出入口。其中，1、2号出入口过工农南路采取顶管法施工。顶管外轮廓尺寸为7.4 m(宽)× 4.9 m (高)，壁厚0.45 m，单根管节长度为1.5 m，覆土深度约5.5 m。附属2号口过街通道长度为54 m，共设置36环管节。顶管周边主要存在Φ300 mm给水管、Φ400 mm污水管以及众多燃气管道。附属2号口过街通道平面图如图1所示。

1.2 工程地质及水文条件

世纪大道站位于南通市崇川区，为三角洲平原地貌单元。该地区地表水资源丰富，主要以饱和砂质粉土为主，且顶管通道主要处于粉土地层，如图2所示。地层的物理力学参数如表1所示。

图2 地层关系Fig.2 Stratigraphic relationship

2 流固耦合数值模型

2.1 有限差分数值模型

在使用FLAC3D建立有限差分数值模型时，根据已有经验[11]，模型的上边界一般取至地表，模型的下边界到顶管底部应取2D(D为管道外径较大值)以上的距离，水平方向的长度应为7D以上；实际操作中可将模型边界适当取大。因此，模型长×宽×高取为54 m ×52 m ×26 m。数值模型如图3所示。

已开挖的土体采用空模型模拟，土体采用摩尔-库伦本构模型，此外用等代层模拟地层损失，等代层和管片均采用弹性模型。模型侧、底面为位移边界，上面为自由边界。

图1 附属2号口过街通道平面图Fig.1 Plane figure of rectangular pipe jacking machine under passing Jingshi Road

表1 地层物理力学参数Table 1 Physical and mechanical parameters of stratum

图3 模型网格划分图Fig.3 Finite difference mesh model

设置管片密度ρ=2 450 kg/m3， 弹性模量E=34.5 GPa，泊松比υ=0.3。参考已有研究[12]，等代层厚度为盾尾间隙与土体系数的乘积，等代层厚度取d=0.2 m，密度ρ=1 920 kg/m3，弹性模量E=3 MPa，泊松比υ=0.2。

结合水文地质条件，设置地下水位为2.5 m。流体均质、各向同性，同时将模型侧、底面设置为不透水边界，开挖面为透水边界。此时渗流便可发生。

2.2 数值模型适用性验证

顶管顶进工作完成后，选取3个典型断面(顶进至15、21、27 m处)的地表沉降数值模拟结果，将其正态拟合，并与现场监测数据进行对比，如图4所示。

从图4可以观察到，顶管顶进工作完成后，距离接收井越近的监测断面地表沉降范围越大，这是因为此时距离接收井越近，受到施工扰动的影响就越大，从而沉降范围越大；另外，距离始发井越近，沉降

图4 沉降对比曲线Fig.4 Settlement curves

槽中心的沉降值也近似呈现增大的趋势，这是因为随着施工完成，前方土体的超孔隙水压力发生了一定程度的消散，从而沉降值也较大。选取的3个典型断面沉降曲线均符合正态分布，且数值模拟地表沉降槽与实测数据沉降曲线大致接近，因此验证了该数值模型的适用性和可靠性，为后续极限支护力分析和破坏模式提取奠定了基础。

3 数值计算结果与分析

3.1 极限支护力分析

σ0=γwh+K0γ′h+K0γh1

(1)

γ′=γsat-γw

(2)

式中：γ为土体重度；γw为水重度；γ′为浮重度；γsat为土体饱和重度；K0为侧压力系数；h为水位至开挖面中心距离；h1为地表至开挖面中心距离。

参考极限支护力数值模拟方法[14]，将顶管隧道开挖到距始发井27 m，停止开挖，采取应力控制法得到极限支护力。得出流固耦合和未考虑渗流影响两种条件下开挖面位移变化曲线，如图5所示。

由图5可知，两种情况下的关系曲线均可划分为3个阶段：以未考虑渗流计算为例，当支护应力比在0.7～1.0范围内，开挖面位移随支护应力比的减小近似呈线性缓慢增长，可认为开挖面处于弹性阶段；支护应力比在0.5～0.7范围内，位移曲线的斜率逐渐增大，此时开挖面处于弹塑性阶段；支护应力比为0.5时，开挖面位移突增，此时发生失稳破坏。两种情况下的极限支护力分别为0.25(流固耦合计算)和0.5(未考虑渗流计算)。考虑流固耦合的情况下，当支护应力比大于0.6时，开挖面几乎不发生位移，且该情况下极限支护应力比较小。这可能是因为在渗流场的影响下，开挖面前方的土拱效应发展程度较大。

图5 开挖面位移与支护应力比关系Fig.5 Relationship between excavation face displacement and support stress ratio

3.2 开挖面破坏模式提取

提取支护应力比分别为0.5、0.4、0.3和0.2(土体失稳破坏)时土体的位移云图，如图6所示。

图6 土体位移云图Fig.6 Soil displacement nephogram

由图6可见，随λ减小，开挖面前方的滑动区向上扩展，位移不断增大；当λ减小至0.2时，土拱效应被破坏，发生整体失稳破坏。此时的失稳破坏模式与吕玺琳等[10]饱和粉砂离心试验结果以及康志军等[13]数值模拟结果相似：开挖面前方破坏区整体呈现楔形体形状，形成延伸至地表的整体失稳区。

为进一步了解土体整体失稳时的破坏形式以及土体状态，通过FLAC3D提取λ=0.2时的塑性区分布，如图7所示。由图7可知，当λ=0.2时，塑性区范围延伸至地表，土体发生整体破坏。

图7 塑性区分布Fig.7 Distribution of the plastic zone

3.3 深径比影响分析

在同一顶管高度(B=4.9 m)，不同埋深(H=5.88、5.5、3.92 m)情况下，研究渗流作用下不同深径比对开挖面稳定性的影响。为分析此变量的影响，在建模时保持水位线位置不变。

提取深径比H/B分别为0.8、1.12(本文工程实例深径比)、1.2时，开挖面前方整体失稳的破坏模式，如图8～图10所示。

从图8～图10可以看出，当顶管隧道埋深较浅时，其破坏模式大致相似，均符合顶部破坏模型。随着深径比的增大，破坏区扩展范围略微减小，土拱效应发展程度变大，但土体最终都发生了整体失稳破坏。这与Chen等[15]的研究结果一致。深径比增大到1.2时，破坏模式顶部近似为拱形。因此可推测当深径比增大到一定值时，由于土拱效应仍然发挥作用，破坏区未延伸至地表。

3.4 土体孔隙水压力分析

提取λ分别为1、0.4和0.2时的土体孔隙水压力云图，如图11所示。

从图11可以看出，在考虑渗流情况下，随λ减小，孔隙水压力减小的范围有所增大。当λ=1时，开挖面变形很小，土体孔隙水压力沿高度方向呈层状；随λ的继续减小，开挖面前方一定距离的孔隙水压力值发生显著变小，产生形似“陀螺状”的降低区域。即随着λ的不断减小，孔隙水压力影响区范围变大，因而破坏区扩展范围变大。

图8 H/B=0.8时的破坏模式图Fig.8 Failure mode diagram when H/B=0.8

图9 H/B=1.12时的破坏模式图(本文工程实例)Fig.9 Failure mode diagram when H/B=1.12 (engineering examples in this paper)

图10 H/B=1.2时的破坏模式图Fig.10 Failure mode diagram when H/B=1.2

图11 孔隙水压力云图Fig.11 Pore pressure nephogram

4 上限研究

4.1 极限分析上限法

为了对流固耦合条件下的边坡失稳破坏进行理论分析，Viratjandr等[16]认为渗透力所做的功等于孔隙水压力在破裂机制上做的功，即

(3)

(4)

将FLAC3D提取的开挖面前方土体不同节点孔隙水压力的值代入式(4)，即可得到孔隙水压力所做的功率。

4.2 破坏模式的构建及支护力的计算

基于数值模拟提取的浅埋顶管主动破坏模式，构建流固耦合条件下开挖面失稳刚性块体模型，如图12所示。

在利用极限分析上限定理求解支护力时，做出以下假定：①隧道问题的计算可看成平面应变问题；②作用在隧道顶部、底部以及侧部的荷载为均布荷载；③体积应变增量为0。

A1为破坏区域面积；v1为破坏滑落速度；L为顶管宽度；σT为极 限支护力；l10为顶管高度图12 破坏模式示意图Fig.12 Diagram of failure mode

重力做功的总功率Pγ为

(5)

支护力做功的总功率PT为

(6)

孔隙水压力做功的总功率Pu为

Pu=uHv1

(7)

内能耗散功率Pc为

Pc=cHv1

(8)

由虚功原理，可得

Pγ+Pu-PT=Pc

(9)

式中：u为加权平均孔隙水压力，其值的计算参考数值模拟结果；c和γ分别为加权平均黏聚力和加权平均重度。

可得极限支护力为

(10)

将式(10)代入Broms等[17]提出的稳定性系数公式，可得考虑渗流条件下浅埋顶管隧道稳定性系数，公式为

(11)

4.3 计算结果对比与分析

基于4.2节相关公式，分别进行极限支护力和稳定性系数随深径比变化分析，如图13和图14所示。

从图13中可以看出，在浅埋饱和砂层中，随着深径比的增加，两种情况下极限支护力均呈线性增长，但当不考虑渗流所做功率时，极限支护力有较大程度的减小。如图14所示，不考虑渗流力的影响时得到的稳定性系数相较于流固耦合计算结果大，

图13 深径比与极限支护力关系Fig.13 Relationship between depth diameter ratio and ultimate support force

图14 深径比与稳定性系数关系Fig.14 Relationship between depth diameter ratio and stability coefficient

其数值皆大于6，围岩偏向于不安全[18]。因此在浅埋顶管隧道施工过程中，应考虑流固耦合效应对开挖面稳定性的影响。

4.4 顶管配筋计算

为更好地了解管片在渗流影响下的受力特征，从而针对管片进行合理配筋，将作用于顶管的荷载进行简化为均布荷载，并绘制其弯矩图，如图15所示。

根据弯矩分配法，可知该结构弯矩最大值位于结构边角处，最大弯矩Mmax为

(12)

式(12)中：q为作用在顶管的均布线荷载，q=1.5σT，kN/m。根据式(12)，易得到在不同深径比H/B和宽高比ξ下的弯矩值，如表2所示。

考虑到工程顶管截面B×L=4 900 mm × 7 400 mm，即宽高比ξ=1.5；深径比H/B=1.12，在该截面面积下若为单筋截面配筋时，结合相关参

图15 管片弯矩图Fig.15 Bending moment diagram

表2 不同深径比和宽高比下的弯矩值Table 2 Bending moments with different H/B and ξ

数，此时截面能够抵抗最大的弯矩值为626.425kN·m，当弯矩值大于该截面单筋截面配筋所能承受的最大弯矩值时，则按双筋截面原则进行配筋。同时，根据表2可知，相比于矩形截面，方形截面在受力上更加合理。故在顶管截面选择时，可在满足构造要求的情况下优先考虑方形截面形式。

5 结论

根据南通矩形顶管过街通道项目，建立了考虑流固耦合效应影响的三维数值模型，对浅埋矩形顶管隧道开挖面的稳定性进行研究，结论如下。

(1)随支护应力比的减小，开挖面位移的发展可划分为3个阶段，且相较于无渗流计算，考虑流固耦合影响时开挖面前方土体土拱效应发挥程度较大，极限支护比远小于1，但最终破坏区延伸到地表，形成整体失稳区。

(2)当支护应力比较小时，开挖面前方一定距的土体孔隙水压力显著降低，产生形似“陀螺”状的降低区域。

(3)深径比对浅埋顶管隧道破坏模式影响不大，其均符合顶部破坏模式；但在浅埋顶管隧道施工中，随深径比的增大，所需极限支护力有所增大，相较于未考虑渗流影响时计算结果更加安全。

(4)基于弯矩分配法，给出了不同深径比和宽高比下的弯矩计算表，结果表明在深径比一定时，相比于矩形截面，方形截面在受力上更加合理。