本发明涉及一种分析深埋砂土地层中盾构隧道开挖面稳定性的方法，尤其 涉及一种考虑深埋砂土地层中的盾构隧道的开挖面极限支护力计算方法，属于 盾构施工技术领域。

背景技术：

盾构施工技术因其机械化程度高、安全风险小、施工速度快、对周围环境 影响小等优点，在地下公路与铁路隧道施工中被广泛采用。隧道施工过程中， 开挖面失稳导致地表沉降过大及地表塌陷事故是隧道施工中遇到的突出问题， 确保开挖面稳定性是隧道顺利施工的首要条件。盾构隧道的开挖面稳定性都是 通过控制开挖面的支护压力来实现的，支护力压力过大或者过小都会导致开挖 面失稳，进而引起地表隆起或沉降，对施工进度和安全施工产生影响，造成巨 大经济财产损失。因此，确定合理的支护力既是实际工程中的关键问题，也一 直是盾构隧道研究工作的热点。

极限平衡法具有原理简单，计算方便等优点。Horn(1961)基于筒仓理论提出 三维楔形体计算模型，模型中通过楔形体和棱柱体的极限平衡，构建力的平衡 关系求解极限支护力，Jancsecz et al.(1994)对该理论进行了完善并给出了考虑 土拱效应的计算公式。在该模型的基础下后续学者也针对渗流(Anagnostou et al. 1994)、非匀质土(Broere et al.1998)、楔形块与土体摩擦(Kirsch et al.2005) 等条件下对该计算模型进行改进工作。

极限分析法由于考虑了岩体的应力应变关系和具有更严格的理论基础，使 得该方法的计算结果更加精确。基于塑性极限分析上限法的基本原理，Leca and Dormieux(1990)构筑了浅埋隧道的三维破坏机制，并分别推导了相应的极限支 护力的表达式。Soubra et al.(2008)提出了一种由多个刚性截锥体组成的新的三 维锥形体破坏机制，提高了极限支护力的计算精度。Mollon et al.(2010,2011)利 用空间离散技术，采用“点对点”的方法建立了开挖面的三维破坏机制，进一步 提高了计算精度。在多椎体破坏模型下，后续学者针对椎体上部荷载(Zhang et al.2015)、非匀质土(Han et al2016),纵向倾斜隧道(Zhao et al，2017)等不 同工况进行了下一步研究工作。

现有的关于开挖面稳定性的理论研究工作多是针对浅埋隧道下提出，无论 是极限平衡法还是极限分析法，这些理论在深埋隧道中都基于一定的不适用性。 如极限平衡法下的楔形体计算模型中，在深埋隧道中的主要问题表现为上部棱 柱体作用于楔形块的荷载如何计算。不考虑土拱效应的计算结果显然是与工程 实际、模型试验得出的结论不符的。而考虑土拱效应的，基于太沙基松动土压 力理论推导的三维松动土压力。虽然可以得到一定深度后不随埋深变化的结论， 但是这一理论推导是基于假定滑动区域垂直延伸至地表，且松动区域截面的截 面参数保持不变所得的，这一假定在隧道处于深埋的条件下显然是难以成立的。

极限分析法下的计算模型中，由LecaandDormieux针对浅埋隧道提出，在 Soubraet al.(2008)一文中通过改进为多椎体，可以得到埋深比C/D大于一定值 后，土体重量影响系数Nγ不变，即极限支护力不变的结论，但是一般认为，无 论是极限平衡法还是极限分析法，都认为合理的破坏模式能更加精确的确定开 挖面的极限支护力，而Soubret al.(2008)文中的破坏形状在深埋隧道中，与现有 模型试验及数值模拟研究中破坏区域对比来说，还是存在一定的出入，有待进 一步改进的。

在松散介质中开挖隧道，伴随着地层应力释放，势必在隧道上方形成一定 的松弛区，当隧道埋深较大时，该松弛区大都会在隧道上方一定范围内形成一 个稳定的平衡拱。前苏联学者普罗托吉雅柯夫也以散体介质理论为基础，提出 了著名的普氏拱理论，认为深埋隧道开挖后在隧道上方形成一个抛物线形的平 衡拱。秦建设(2005)模拟挡板下落试验得到挡板上部区域土体位移等值线呈 抛物线状。Zhanget al.(2014)在深埋隧道的塌方试验中，对塌落拱边界进行曲 线拟合也验证了深埋隧道中，塌陷失稳形成的天然拱呈抛物线形。Asgar et al. (2010)利用PIV技术对得到了隧道模型试验下的砂土位移云图，得到松动区边界 线呈现出明显的抛物线拱形状，Zhang et al.(2015)在其文章中对埋深比为1的 砂土中进行破坏区域的界定，松动区边界也近似为一个抛物线拱的形状。

目前，整体被动失稳假设对分析局部被动破坏模式具有局限性，而对局部 被动破坏的发展与失稳机理的研究较少。

技术实现要素：

本发明的目的在于提供一种考虑局部被动破坏的砂土地层浅埋隧道开挖面 极限分析方法。为实现上述目的，本发明采用的技术方案具体如下：

S1.假定土体破坏模式，提出理论模型。

假设土体均质、各向同性，满足摩尔库伦屈服准则，考虑深埋隧道的土拱 效应，认为深埋砂土盾构隧道的破坏区域轮廓如图1所示，图中隧道直径为D， 埋深为C，开挖面破坏区域由①和②两部分组成，其中①为抛物线旋转而成三 维塌陷区，任一垂直向切面为抛物线形状，满足普式拱理论，塌陷区土体以自 重作用于区域②上。区域②为一剪切滑动破坏区，该区域沿倾斜滑动面剪切破 坏。

根据图1建立如图2所示的破坏模式的几何构成，来分析深埋砂土隧道的 开挖面极限支护力。图2中，隧道直径为D，埋深为C，开挖面的支护压力σT假 设为均布荷载。开挖面前第i个刚性截椎体由第i个圆锥体被截面i+1和截面i+1 个截面所截得到(i＝1,2,3,4,5)，第6个截面为第5个圆锥体被过隧道顶部，平 行于隧道方向的平面所截。所有截面均为椭圆形。计算第6个截面上方松动区 土体重量时候，参考极限平衡法中的隧道面积与楔形块面积等效的方法，将塌 陷区下口假定为与椭圆截面6等面积的圆。

开挖面前方第i个截面和第i+1个截面的夹角为βi，第i个刚性截锥体速度 方向与水平方向的夹角为αi，第i个椭圆椎体的特征参数为θi。假定椭圆截面i 的长轴长度为ai，短轴长度为bi。

S2.根据理论计算模型，推导理论公式。

理论推导过程如下

假定角度变量为βi和α1，将α1记为α。根据几何关系有：

第i个椭圆截面的长轴ai、短轴bi、面积Ai的值有：

Ai＝πaibi (5)

如图2所示，破坏区域②由5个刚性滑块组成，刚性滑块为椭圆椎体被平 面所截形成，记为滑块i(i＝1,2,3,4,5),则滑块i的体积Vi为：

式中，hi表示椭圆椎体被截面i截取后剩余椭圆椎体的椎体高度， h5′第截面6截取后剩余椭圆椎体的椎体高度，且有

而破坏区①为抛物线旋转所得，参考根据极限平衡法的筒仓理论中，楔形 体尺寸由隧道断面尺寸确定的方法，即破坏区①下口为一标准圆，由圆面积S＝A6＝πR2可知

根据普式拱理论，塌陷区高度H＝R/f，f为普式系数，砂土中塌落拱曲线函数为f(x)，则塌陷区土体体积为：

刚性椭圆锥体的速度方向与椭圆锥体的轴线方向平行，假设第i个刚性截椭 圆锥体的绝对速度为vi，第i个刚性截椭圆锥体和第i+1个刚性截椭圆锥体的相 对速度为vi,i+1。

根据图3的速度矢量图，绝对速度与相对速度的递推关系有：

根据极限分析法，当隧道开挖面保持稳定时，根据外力功率与内能耗散功 率相等的原则，则有：

Pe≤Pv (11)

其中，Pe指的是外力虚功率，Pv指的是破坏区范围内部能量耗散率。外力 虚功率Pe包括开挖面支护力的虚功率PT、塌陷区土体重力功率Psoil以及滑动区土 体重力功率Pγ组成，即

Pe＝PT+Psoil+Pγ (12)

开挖面支护力的虚功率PT为

PT＝-A1 cosασTv1 (13)

塌陷区②土体重力功率Psoil为

Psoil＝Vsoilγsin(2β2+2β4+α)v5 (14)

破坏区②中土体重力功率Pγ为

破坏区②中能量耗散率Pv为

将式(13)～(17)代入式(12)并使外力的虚功率等于破坏区内部能量耗 散率，取α,β1,β2,β3,β4为变量，可得砂土地层中深埋盾构隧道开挖面的极限支 护力上限解为：

式中Nγ表示土体重度对开挖面极限支护力的影响系数，且

根据破坏模式的相容速度约束和几何条件约束，确定角度变量 α,β1,β2,β3,β4的变化范围，利用MATLAB软件对破坏区的形状进行优化，砂 土地层中深埋盾构隧道开挖面极限支护力的最优解。

S3.建立数值3D模型，计算数值模拟中的极限支护力，并与步骤2中的结 果进行对照。

(1)建立数值模型

采用有限差分软件FLAC3D建立数值模型如图4所示，模型的边界约束条 件为：上表面为自由边界，四周为法向约束，底部为固定约束。模型的边界约 束条件为：上表面为自由边界，四周为法向约束，底部为固定约束。本模拟计 算中材料的破坏屈服采用摩尔-库伦剪切破坏准则，土层的弹性模量25MPa，泊 松比根据v＝k0/(1+k0)，确定。盾构隧道的开挖是一个逐渐推进的 过程，考虑到本文分析的重点是不同工况下开挖面极限支护力的确定，数值计 算中，采用一次开挖到25m，同时对开挖面后方隧洞洞周各节点施加固定约束 模拟用以模拟支护，在开挖面施加与隧道中心点原始水平地应力(σ0)相同的 支护压力σT，通过不断减小支护压力σT的大小，获取开挖面中心点位移与支护 力σT的关系曲线。当开挖面支护压力减小到某一值，隧道开挖面中心点处变形 急剧增大，即认为此支护压力为极限支护压力。

(2)给出数值计算极限支护力σT大小的方法

采用数值模型方法设置30种工况，计算参数如表1所示。

表1.计算参数

数值模拟得到了30种工况下开挖面支护力与开挖面中线点水平位移的关系 曲线，如图5所示。可以看出，随着开挖面支护力的不断降低，开挖面中心点 水平位移不断增大，且当开挖面支护力降低到一定程度即极限支护力时，开挖 面中心点水平位移将急剧增大，由此可以确定30种工况下的开挖面极限支护力。

本发明提出模型可以计算砂土地层深埋隧道开挖面极限支护力，且该支护 力不随埋深变化而变化。采用同样计算参数，将本文提出模型计算得到的开挖 面极限支护力与数值模拟结果进行对比，如图6所示。可以看出，随着隧道直 径的增大，两种方法计算得到的开挖面极限支护力均不断增大，呈线性关系， 两条直线的斜率相近，本文理论模型计算结果要大于数值模拟结果，这与很多 理论模型结果是一致的，表明本文模型计算结果是偏于安全的，具有一定的安 全储备；在地层内摩擦角增大时，本文解析解和数值模拟解均随着隧道直径的 增大而减小，近似呈对数关系，且两条曲线随地层内摩擦角的变化规律较为一 致。

鉴于连续介质模型是无法直接计算得到开挖面失稳破坏区域，本文采用张 成平(2015)的方法确定开挖面前方破坏区的边界，取工况26～30中通过FLAC3D 数值模拟得到的破坏区域，与本文计算出的破坏区域进行对比。对比结果如图7 所示。图7表明了本文提出的模式计算出的破坏区与数值模拟得到的破坏区域 具有较好的一致性，且当地层内摩擦角变化时，依然吻合较好。

附图说明

图1是所提出的深埋砂土隧道的破坏区域轮廓图；

图2是所提出的三维破坏模式的几何构成图；

图3是所提出的三维破坏模式的速度矢量图；

图4是数值模拟使用的数值模型图；

图5是开挖面支护力与中心点水平位移的关系曲线；(a)工况1～5，(b) 工况6～10，(c)工况11～15，(d)工况16～20，(e)工况21～25， (f)工况26～30

图6本发明理论解与数值模拟对比；

图7内摩擦角变化时，理论计算和数值模拟得到的破坏区范围对比；

图8本发明计算方法与现有计算方法对比

图9本发明计算方法与离心实验极限支护力对比

图10本发明计算方法求得破坏区范围与离心实验对比

具体实施方式

以下结合附图和现有其他计算方法对比，进一步详细阐述本文提出计算方 法。

(1)与现有其他计算方法比较

假设无水砂土地层中，深埋隧道(C/D＝2)的土体重度为γ＝18Kn/m3， 用本文方法计算D＝6m,8m,10m,12m,14m和D＝10m，等 ＝15°,20°,25°, 30°，35°工况下的极限支护力，并将本文方法与现有研究方法进行对比。图6 给出了本文极限支护力最优解对应的极限支护力随着隧道直径D，地层内摩擦 角的变化曲线，并与已有开挖面稳定性极限支护力的成果(leca et al,1990；Anagnostou et al,1994；soubra et al,2008)进行对比。可以看出：在深埋砂 土隧道中(C/D＝2)，当以隧道直径D为变量时，本文提出的破坏机制得出的 开挖面极限支护力比leca et al.(1990)的单块体及二块体破坏模式、soubra et al (2008)的多椎体破坏模式等以往成果获得的极限支护力大；在深埋砂土隧道 中(C/D＝2)，当以地层内摩擦角为变量，本文提出的破坏机制得出的开挖面 极限支护力比leca et al.(1990)、soubra et al(2008)得出的极限支护力大；另 外，从图10可以看出，本文理论计算结果与已有成果(leca et al,1990；soubra et al,2008)类似，采用极限分析法获得的极限支护力解都远小于极限平衡法 (Anagnostou et al,1994)获得的解。且各种方法取得的极限支护力随着隧道直 接D、地层内摩擦角的整体变化一致。

(2)与离心实验对比

根据Gregor(2011)的砂土离心实验中，离心实验原型隧道 D＝5m，γ＝15Kn/m3,当C/D＝1.5时，极限支护力为8.5kPa。根据该实 验提供的参数，运用已有理论方法进行计算对比。对比结果见图9。在深埋砂土 隧道中，地层摩擦角从图9可得此时本文方法得到的极限支护力比leca et al.(1990)、soubra et al(2008)得到的大，最为接近离心实验取得的极限支 护力。验证了本文的破坏机制具有更好的精确性和优越性。

Chambon(1994)通过离心机试验得到了无粘性土地层开挖面前方土体坍塌 的主动破坏机制，并将破坏区域画了出来。如图10所示，当C/D＝1，D＝5m，将本文提出的破坏机制轮廓与上述离心机试验结果进行对比。从图中可以看出， 破坏区下端的滑动区与离心实验具有较好的一致性。上部塌陷塌陷区宽度分别 为0.35D与0.3D，较为接近。塌陷区高度差距较大，结合相关研究(Zhang et al,2014)分析，认为在离心实验进行过程中，达到极限状态时会形成一个塌落 拱，随着时间的推移，塌落拱拱脚位置变化较小，拱顶不断向上扩大发展，形 成一个高度较高的塌陷区。

参考文献

Ahmed,M.,Iskander,M.,2012.Evaluation of tunnel face stability by transparent soil models. Tunn.Undergr.Space Technol.27,101–110.

Anagnostou,G.,Kovári,K.,1994.Dieder Ortsbrust bei Erddruckschilden. Mitteilungen der Schweiz.Gesellschaft für Boden-und Felsmechanik 129,27–34.

Ansgar,K.,2010.Experimental investigation of the face stability of shallow tunnels in sand.Acta Geotechnica.5:43-46

Broere,W.,1998.Face stability calculations for a slurry shield in heterogeneous soft soils.In Negro Jr.and Ferreira，editors，Tunnels and Metropolises，Rotterdam,Balkema.

Chambon,P.,Corté,J.F.,1994.Sallow tunnels in cohesionless soil:stability of tunnel face.J. Geotech.Eng.120(7),1148–1165.

Gregor,I.,Oelin,A.,WEI,W.,Ronaldo,I.B.,2011.Centrifuge model test on the face stability of shallow tunnel.Acta Geotechnica.6:105-117.

Han,K.H.,Zhang,C.P.,Zhang,D.L.,2016.Upper-bound solutions for the face stability of a shield tunnel in multilayered cohesive–frictional soils.Comput.Geotech.79,1-9.

Horn,N.,1961.Horizontal earth pressure on the vertical surfaces of the tunnel tubes.In:National Conference of the Hungarian Civil Engineering Industry,Budapest,November,pp.7-16(in German).

Jancsecz S.Steiner W.,1994.Face Support for large Mix-Shield in Heterogeneous Ground Conditions[J].Tunnelling 94,London:531-550

Janelid,I.,Kvapil,R.,1966.Sublevel caving.Int.J.Rock Mech.Min.Sci.3(2),129–153. Kirsch,A.,Kolymbas,D.,2005.Theoretische Untersuchung zurBeautechnik 82/7,449-456

Leca,E.,Dormieux L.,1990.Upper and lower bound solutions for the face stability of shallow circular tunnels in frictional material.Géotechnique,40(4),581–606.

Li,P.F.,Wang,F.,Fan,L.F.,Wang,H.D.,Ma,G.W.,2019.Analytical scrutiny of loosening pressure on deep twin-tunnels in rock formations.Tunn.Undergr.Space Technol.83,373-380.

Messerli,J.,Pimentel,E.,Anagnostou,G.,2010.Experimental study into tunnel face collapse in sand.Physical Modelling in Geotechnics,editied by Springman,Laue& Seward,Zurich,Vol.1,pp.575-580.

Mollon,G.,Dias,D.,Soubra,A.H.,2010.Face stability analysis of circular tunnels driven by a pressurized shield.J.Geotech.Geoenviron.Eng.136(1),215-229.

Mollon,G.,Dias,D.,Soubra,A.H.,2011.Rotational failure mechanisms for the face stability analysis of tunnels driven by a pressurized shield.Int.J.Numer.Anal.Methods Geomech.35(12), 1363-1388.

Mollon,G.,Dias,D.,Soubra,A.H.,2013.Continuous velocity fields for collapse and blowout of a pressurized tunnel face in purely cohesive soil.Int.J.Numer.Anal.Methods Geomech.37, 2061–2083.

Oblozinsky,P.,KUwano,J.,2004.Centrifuge experiment on stability of tunnel face.slovak journal of civil engineering.3:23-29.

Qin,J.S.,2005.Study on Face Deformation and Collapse of Earth Pressure Shield Tunnel.Ph.D. thesis,Hohai University,Nanjing,China

Soubra,A.H.,Dias,D.,Emeriault F.,2008.Three-dimensional face stability analysis of circular tunnels by a kinematical approach.Geocongress:Characterization,Monitoring,and Modeling of Geosystems.New Orleans,(311),894-901.

Thomas,J.P.2007.Ground movements during tunneling in sand.Master thesis,Aueen’s university,Kingston,Ontario,Canada

Zhang,C.P.,Han,K.H.,Zhang,D.L.,Li,H.,Cai,Y.,2014.Test study of collapse characteristics of tunnels in soft ground in urban areas.Chinese Journal of Geotechnical Engineering. 33(12),2433-2442

Zhang,C.P.,Han,K.H.,Zhang,D.L.,2015.Face stability analysis of shallow circular tunnels in cohesive–frictional soils.Tunn.Undergr.Space Technol.50,345-357.

Zhao L.H.,Li,D.J.,Li,L.,Yang,F.,Cheng,X.,Luo,W.,2017.Three-dimensional stability analysis of a longitudinally inclined shallow tunnel face.Comput.Geotech.87,32-48.