考虑接触非线性的盾构隧道纵向分析模型及其应用

耿 萍，陈枰良，陈昌健，杨 琪，郭翔宇，何 川

（西南交通大学交通隧道工程教育部重点实验室，四川成都610031）

断层经常造成岩层之间的不整合面，可能对隧道、地基等的力学性质造成影响。活动断层发生日积月累的蠕滑则会对隧道、地面建筑等造成不可逆的影响，因此各国建筑规范对在活动断层区域进行工程建造有严格的要求。然而，对于隧道等类型的线状工程，有时不可避免需要穿越蠕滑断层，因此，在长期蠕滑作用下，结构力学性能预测尤为重要。

因为不连续的结构构造，盾构隧道力学分析为各类穿越蠕滑断层隧道分析中最为复杂的分析之一，其面临的主要困难是盾构隧道纵向接头的非线性问题。目前盾构隧道纵向接头分析模型主要有2类：管片-环缝模型和接缝界面模型。其中，管片-环缝模型中接头力学分析考虑管片和环缝参与受力，基于盾构隧道分别在轴拉力、轴压力、纯弯矩作用下进行分析，环缝分别处于全张开、全闭合、半张开状态，常采用考虑刚度折减的连续梁等效。文献［1］最早研究了该类模型，通过引入平截面假定，推导了盾构隧道纵向接头在纯轴力和纯弯矩作用下的轴向刚度和弯曲刚度解析解。接缝界面模型中接头力学分析不考虑管片参与受力，以接缝在外力作用下的力学行为为分析对象，常采用轴向、剪切、弯曲弹簧模拟。文献［2-3］对隧道纵断面方向上环间接头转动弹簧系数进行了理论推导。2类方法中管片-环缝模型应用最广泛，徐凌［4］、黄正荣［5］、张文杰等［6］、叶飞等［7-8］、李翔宇等［9］通过理论分析、模型试验等手段研究了横向变形、环缝作用范围和螺栓预紧力等因素对盾构隧道纵向等效抗弯刚度的影响，改进了该模型；黄亮等［10］采用与志波模型相同的思路，建立了类矩形盾构隧道纵向弯曲刚度的计算方法；文献［11-14］则重点关注了纵向弯曲和剪切效应；卢慈荣［15］、耿萍等［16］、张景等［17］针对弯矩和轴力共同作用下盾构隧道纵向抗弯刚度非线性变化进行了研究，建立了轴力、弯矩共同作用下的纵向接头力学分析方法。

综合来看，盾构隧道纵向结构为传递轴力、弯矩、剪力作用的构造，具有轴向、转动、切向3个自由度。现有纵向接头分析模型进行力学分析时，一般将3个自由度独立考虑，然而，在一些特殊工况下盾构隧道同时传递着轴力、弯矩和剪力，如盾构隧道在蠕滑位错荷载作用下的力学响应、在斜入射地震动作用下的地震响应，传统纵向接头分析模型不能考虑轴力、弯矩、剪力共同作用下的接头力学特性，使得传统分析方法在特殊工况下不适用。

笔者曾对轴力和弯矩共同作用下的盾构隧道纵向抗弯刚度进行了解析［16］，然而该模型亦存在如下缺点：模型由多种接触状态分别表达，需要分状态计算，且受压区计算高度复杂；仅考虑轴力和弯矩的共同作用，未考虑剪力作用。本文通过引入基于位移的接头力学分析方法，建立轴力、剪力、弯矩共同作用下考虑接触非线性的盾构隧道纵向分析模型；然后以汕头湾盾构隧道为背景，分析蠕滑作用下盾构隧道结构内力响应。与轴力和弯矩共同作用模型相比，2个模型采用了相同的基本假定，但轴力和弯矩共同作用模型以轴力和弯矩为已知量。本文模型以位移为已知量，主要有如下改进：在相同假定的条件下，能考虑轴向、切向、转动3个自由度的相互影响；采用基于位移的推导，将环间接缝多种接触状态使用统一力学方程表达。

1 盾构隧道纵向分析模型推导

1.1 盾构隧道纵向力学分析的基本问题

当考虑轴力、剪力、弯矩共同作用时，取分析区域为相邻管片环中心截面间的区域，分别为截面1和截面2，分析区域长度为ls，如图1所示，其中作用于截面1和截面2上的位移和力的大小为(ui，vi，θi)和(Ni，Vi，Mi)（i=1，2）。

图1 盾构隧道纵向力学分析Fig.1 Basic problems in longitudinal mechanical analysis of shield tunnel

分析区域控制方程可用面2与面1的相对位移表示，如式（1）所示。

式中：M、V、N为截面2上的作用力；u、v、θ为截面2与截面1的相对位移，见式（2）～（4）。

因此，盾构隧道纵向分析的基本问题可表述为：在考虑隧道结构存在环间接缝的情况下，确定如式（1）所示的3个自由度上力与位移的关系，并计算隧道结构的环间接缝张开量、螺栓拉力、管片应力的问题。

引入如下基本假设：①假定管片环中心截面在同一平面内变形；②相邻管片环中心截面范围内的平均变形服从平截面假定；③环间接缝受压侧仅管片受压，螺栓不参与受力，受拉侧管片和螺栓共同受拉；④连接件（如螺栓）在环向的离散分布转化为等效的连续均匀分布；⑤不考虑环间接缝的滑移。

1.2 基于位移的力学分析方法

通过尝试基于力的方法和基于位移的方法，发现：当采用基于力的方法时需要建立力平衡方程和变形协调方程，联立求解常得到以受压区高度为未知量的隐式超越方程，使得内力和变形的求解十分复杂；当采用基于位移的推导方法时，接头分析区域变形分布容易求得，引入考虑了接头不连续条件的物理方程后，直接积分可得内力，使得问题的分析得到简化。因此，本文采用基于位移的接头力学分析方法，如图2所示。该方法以位移作为已知量，利用平截面假定得到分析区域变形函数，然后建立盾构隧道纵向力学方程。

图2 基于位移的力学方程推导Fig.2 Derivation of mechanical equation based on displacement

1.3 盾构隧道纵向力学方程推导

1.3.1 分析区域变形函数

由平截面假定可知分析区域的挠曲变形符合梁的挠曲变形函数，轴向变形符合轴心受力构件变形函数。若将分析区域截面1固定，在截面2上施加位移U=(u，v，θ)T，则分析区域的挠曲函数和轴向变形函数分别如式（5）、式（6）［18］所示。

式中：ξ=x/ls，x为以截面1中心点为原点、指向截面2中心点坐标轴上的坐标，如图3所示；ls为管片环宽度。

由式（5）可得，各截面相对于截面1的转角为

则各截面处的曲率κ(ξ)和轴向变形率e(ξ)分别如式（8）、式（9）所示。

将ξ=x/ls代入得

将u、v、θ作用下的分析区域变形分解为常曲率变形和变曲率变形两部分，分解后常曲率变形由常曲率的转动变形和轴向变形组成，变曲率变形为线性变化曲率的转动变形，如图3所示。

图3 纵向受力分析示意Fig.3 Diagram of longitudinal mechanical analysis

对于常曲率变形部分，曲率变化率应为零，即

可得

则常曲率变形和变曲率变形对应的位移量为

式中：Up为常曲率变形部分的位移量；Us为变曲率变形部分的位移量；U为分析区域总的相对位移量。

将式（14）代入式（10），则常曲率变形的曲率函数为

将式（15）代入式（10），则变曲率变形的曲率函数为

1.3.2 物理方程及接头不连续条件

根据环间接缝受压侧仅管片受压、螺栓不参与受力的假定，可知长为ls的分析区内其抗压刚度可表示为

受拉侧由管片混凝土和螺栓共同受拉，则在长为ls的分析区内其抗拉刚度可表示为

式中：Ec为管片弹性模量；t为管片厚度；ζ为分析区拉压刚度比；kj=nKj/(2πr)，Kj为单根螺栓沿隧道纵向的抗拉刚度，n为环间接头螺栓总个数。从上述抗压和抗拉刚度定义可以发现，该刚度为沿盾构隧道横截面圆周方向的线性刚度。

由于环间接缝不连续，导致隧道沿圆周方向上的受拉刚度和抗压刚度不同。令某截面y处应变为ε(y)，则考虑接头不连续性后的截面分布力为

显然，考虑了接头不连续性后的应力函数为分段函数，这使得内力计算需要分状态进行，计算较为繁琐，因此将其改写为式（23）。在该式中，当ε(y)为正时kc项被约去，为负值时kt项被约去，因此该式与式（22）是等效的。

1.3.3 盾构隧道纵向力学方程推导

由1.3.1节分析区域变形函数的分解，分别推导常曲率变形和变曲率变形的力学方程，然后将其组合得到盾构隧道纵向力学方程。

1.3.3.1 常曲率变形力学方程

由式（11）和式（17），可得截面上各点的平均应变函数为

则考虑了接头不连续条件后的分布力函数为

式（26）采用了极坐标表达方式，其中φ的含义如图4所示。

图4 积分变量φFig.4 Integral variable φ

如图4所示，在Up作用下，面2上作用u、θ位移时的作用轴力和弯矩分别为

其中，下标p表示该力为Up作用下的力。将式（25）、（26）代入式（27）、（28）得

其中，式（31）、（32）积分有解析解，式（33）、（34）为超越积分函数。因此，式（33）、（34）引入辛普森数值积分方法，可得Up作用下的力为

令

其物理含义为，转动θ与轴向u在分析区上缘点处引起的位移之比，称为位移偏心距。

再将式（19）、（21）代入式（35）、（36），得

式中：λ1、λ2为与位移偏心距和拉压刚度比相关的系数，称轴向刚度折减系数和弯曲刚度折减系数；A为管片环横截面的面积，A=2πrt；I为管片环的近似惯性矩，I=πr3t。

由式（38）、（39）即可求得Up作用时的轴力和弯矩项。此外，盾构隧道纵向力学分析需要分析纵向接头受压区高度、环间接缝张开量、接触压力等力学响应。

若令环间接缝的上部和下部的变形量分别为d1和d2，则其表达式为

由图2可知，随着砂粒含量的增大，砂质黄土的渗透系数总体呈增大趋势，亦即随着砂粒含量的增加砂质黄土的渗透性逐渐增强。同时，渗透系数与砂粒含量关系曲线形态呈现出先增大、后平缓、再增大的折线形，砂粒含量为35%和40%为曲线形状变化的 “转折点”，即当砂粒含量大于35%时，砂质黄土的渗透系数增幅减小，达到40%时，其渗透系数急剧增大。

由d1、d2可以判断环间接缝的接触状态：①当d1＞0，d2＞0时，为无接触；②当d1≤0，d2≤0时，为全接触；③当d1≤0，d2＞0时，为上部接触；④当d1＞0，d2≤0时，为下部接触。

对于状态③上部接触和④下部接触，曲率κp(x)表示该截面处的转角变化率，e(x)为该截面的轴向变形率。图5所示为接缝下部接触示意图，所以接缝截面的受压区和受拉区的临界坐标可表示为

图5 盾构隧道环间接缝接触示意Fig.5 Schematic diagram of contact in longitudinalshield tunnel

相应的极坐标圆心角为

式中：hc以向y轴正方向为正；φc以自z轴向y轴方向为正。

环间接缝的上沿和下沿纵向应力由式（45）、（46）求得

式中：p1为上沿纵向应力；p2为下沿纵向应力。

1.3.3.2 变曲率变形力学方程及组合方程

对于变曲率变形部分，轴向变形量为零，又由环缝对应的x坐标为ls/2，代入式（18）得环缝截面曲率为

可以发现，按常曲率变形和变曲率变形分解后，常曲率部分在环缝面位置处曲率和轴向变形率不为零，由此决定了接缝面的接触状态。变曲率变形在环缝位置处的曲率函数为零，轴向变形率亦为零，环缝面不产生相对转动和轴向变形，不影响接缝面的法向接触状态，因此变曲率变形的刚度性质应与常曲率变形相同。

由此，变曲率变形部分剪力项和弯矩项为

式中，EcIλ2为常曲率变形下的盾构隧道纵向的等效弯曲刚度。

将Up与Us作用下产生的力对应求和，可得截面2施加U=(u，v，θ)位移时的力为

写为矩阵形式得

其中

式（55）、（56）即为盾构隧道纵向的力学控制方程。可以发现，该式实际上与轴向刚度为EcAλ1、弯曲刚度为EcIλ2的梁构件控制方程相同，因此可采用等效梁模拟其力学性质，其轴向和弯曲刚度为考虑折减系数λ1和λ2的刚度。

2 蠕滑作用下盾构隧道纵向力学行为分析

2.1 工程背景

汕头至汕尾铁路途经的汕头湾海底隧道位于广东省汕头市濠江区和龙湖区境内，线路呈东北—西南走向，属控制性重点隧道工程。工程场区主要有丘陵地貌、海底地貌及平原地貌，隧道全长9 771 m，其中越岭段长约5 741 m，穿汕头湾海底段约2 000 m，穿三角洲平原段约2 030 m，最大埋深约180 m，汕头湾海底隧道纵断面如图6所示。

图6 汕头湾海底隧道纵断面Fig.6 Transverse section of tunnel in Shantou Bay

盾构隧道主体结构拟采用单层装配式管片衬砌，内径12.8 m，外径14.0 m，管片环采用“9+1、1/3封顶块”型式，由1个封顶块（F）、2个邻接块（L1、L2）和7个标准块（B1～B7）组成，平均幅宽2.0 m。混凝土强度等级为C60，抗渗等级为P12，块与块之间每块设置3颗M36环向斜螺栓，螺栓机械等级10.9级。环与环之间每环布置56颗M36纵向斜螺栓，螺栓机械等级10.9级。螺栓中心线与管片接缝交点位于0.5倍厚度处，螺栓中心线与管片径向线近似呈60°，盾构隧道横断面如图7所示。

图7 汕头湾海底盾构隧道横断面Fig.7 Cross-section of tunnel in Shantou Bay

盾构段隧道横穿桑浦山活动断裂带，其平面关系如图8所示。该断层为正断层，宽约20 m。根据地勘资料及活动断裂带工程活动性专题研究报告，针对汕头湾海底隧道近场新建GNSS（全球导航卫星系统）的观测资料，通过解算得到的速度场分布结果知断层主要产生垂直运动，推断断层垂直错动速率为（0.15±0.44）mm·a-1。

图8 汕头湾海底隧道与断层关系平面Fig.8 Relationship between tunnel and fault in Shantou Bay

2.2 蠕滑作用下盾构隧道力学分析方法

采用图9所示流程分析蠕滑作用下盾构隧道的力学行为，主要流程如下：

图9 蠕滑作用下盾构隧道纵向力学分析流程Fig.9 Process of longitudinal mechanical analysis of shield tunnel under creep action

（1）建立蠕滑位错作用下的地层宏观计算模型。汕头湾隧道地层计算模型采用连续介质地层计算模型，如图10所示，地层物理力学参数见表1。通过计算得到隧道所在位置地层变形。

表1 地层物理力学参数Tab.1 Mechanical parameters of soil

图10 蠕滑作用下地层计算模型Fig.10 Soil calculation model under creep action

（2）建立盾构隧道纵向节点-接头结构计算模型。将隧道所在位置地层变形作为荷载输入，如图11所示，为考虑接头非线性，采用如图12所示流程计算。

图12 盾构隧道纵向节点-接头结构模型计算流程Fig.12 Process of longitudinal structure analysis of shield tunnel under creep action

（3）提取结构内力和节点位移，由节点位移计算纵向接头力学响应。

汕头湾盾构隧道结构计算采用如图11所示的节点-接头结构计算模型，其中接头采用轴向刚度为EcAλ1、弯曲刚度为EcIλ2的等效梁元件模拟，梁长与管片环幅宽相同。计算参数列于表2，其中r为隧道计算半径，t为衬砌厚度，ls为管片环幅宽，Kj为环间接头单根螺栓沿隧道纵向的抗拉刚度，kt为地基弹簧纵向刚度，kn为地基弹簧横向刚度。

表2 盾构隧道结构计算参数Tab.2 Structural calculation parameters of shield tunnel

图11 盾构隧道纵向节点-接头结构模型Fig.11 Calculation model of longitudinal nodejoint of shield tunnel

2.3 计算工况

根据汕头湾隧道地勘资料及活动断裂带工程活动性专题研究报告，桑浦山断裂带主要表现为垂直错动，垂直错动速率为（0.15±0.44）mm·a-1，因此计算工况选取了蠕滑速率介于（0.15～0.59）mm·a-1的3种情况，蠕滑时间取介于20～120 a的5种情况，具体工况见表3。

表3 计算工况Tab.3 Calculation condition

将蠕滑位移作为荷载，输入图10的地层计算模型，得到隧道轴向位置处的竖向（Y轴向）和水平向（X轴向）地层位移。图13为蠕滑速率为0.59 mm·a-1（工况1）时隧道位置处的地层位移，其中隧道与断裂带交叉区段的隧道纵轴向坐标为191～212m。可以看出，由于断裂带区域地层刚度相对较小，当断层上下盘发生相对位移时，断裂带区域变形最大。

图13 蠕滑速率为0.59 mm·a-1隧道位置处地层位移Fig.13 Ground displacement at tunnel location at a creep rate of 0.59 mm·a-1

2.4 盾构隧道结构力学响应

图14为计算所得的盾构隧道150～250 m段的结构内力图，可以看出：①桑浦山断裂为正断层，因此隧道结构表现为轴力为拉力，远离断层带时轴力先增大后减小，表现出上盘隧道底部受拉、下盘隧道顶部受拉；②不同蠕滑速率下隧道结构内力分布形式基本不变，仅量值大小发生变化，以0.15 mm·a-1和0.59 mm·a-1蠕滑为例，设计使用100 a时内力增大约3倍。

图14 考虑接头非线性时盾构隧道纵向内力响应Fig.14 Longitudinal internal force of shield tunnel considering nonlinear stiffness

图15比较了采用常数刚度和非线性刚度情况下结构内力响应的差异，图中λ1、λ2分别为采用常数轴向刚度和常数弯曲刚度分析时的刚度折减系数，可以看出：①轴向刚度对轴力和剪力影响大，对弯矩影响小，采用常数轴向刚度分析时，轴力分布影响大，弯矩分布基本无影响，轴力分布由非线性变为接近线性，相差最大时，非线性刚度下的轴力是常数轴向刚度下的3.3倍；②弯曲刚度对弯矩影响较大，对轴力和剪力影响较小，弯曲刚度对弯矩分布形式基本无影响，相差最大时，非线性刚度下的弯矩是常数弯曲刚度下的1.4倍；③蠕滑作用下，采用常数刚度与非线性刚度时的盾构隧道结构内力差异较大，考虑内力共同作用下的接触非线性是有必要的。

图15 采用非线性和常数刚度时盾构隧道纵向内力响应（0.30 mm·a-1，100 a）Fig.15 Longitudinal internal force of shield tunnel with nonlinear and constant stiffness(0.30 mm·a-1,100 a)

3 结论

盾构隧道因由管片环拼接而成，其纵向力学性能表现出显著的接触非线性现象，如在活动断裂带蠕滑作用下，其环间接缝变形与传统的将轴向、弯矩性能独立考虑的情况有较大差异，因此本文主要针对考虑轴力、剪力、弯矩共同作用下环间接缝接触非线性的盾构隧道的纵向力学分析方法做了研究，在此基础上分析了蠕滑作用下盾构隧道结构内力响应，主要有如下结论：

（1）针对盾构隧道纵向接头的非线性接触变形问题，提出了一种盾构隧道纵向分析模型，该模型能考虑盾构隧道纵向接头同时传递轴力、剪力、弯矩时的接触非线性，可以更精确地分析盾构隧道在复杂载荷作用下的力学响应。

（2）当考虑盾构隧道纵向同时传递轴力、剪力、弯矩时，其传力作用可等效为考虑轴向和转动刚度折减的梁单元，轴向和转动刚度折减系数与位移偏心距有关，可由位移偏心距计算梁的等效刚度。

（3）汕头湾盾构隧道在蠕滑作用下的力学分析表明，在不同蠕滑量时，隧道纵向内力分布规律相同，纵向接头非线性对隧道结构内力响应影响较大，考虑接头非线性时轴力可达到常数刚度下的3.3倍，弯矩可达到常数刚度下的1.4倍，考虑接头的非线性刚度是有必要的。

作者贡献说明：

耿萍：制定研究方案和论文撰写，全文审阅。

陈枰良：理论推导，研究框架，论文撰写。

陈昌健：非线性算法实现，论文撰写。

杨琪：数值模拟，论文审阅及修改。

郭翔宇：提供研究思路和技术指导。

何川：提供研究支持，论文审阅。