Prediction of subsurface settlement induced by shield tunnelling in sandy cobble stratum

Fan Wang,Xiuli Du,Pengfei Li

Key Laboratory of Urban Security and Disaster Engineering,Ministry of Education,Beijing University of Technology,Beijing,100124,China

Keywords: Shield tunnelling Sandy cobble stratum Subsurface settlement Volumetric deformation mode Stochastic medium theory

ABSTRACT This study focuses on the analytical prediction of subsurface settlement induced by shield tunnelling in sandy cobble stratum considering the volumetric deformation modes of the soil above the tunnel crown.A series of numerical analyses is performed to examine the effects of cover depth ratio (C/D),tunnel volume loss rate (ηt) and volumetric block proportion (VBP) on the characteristics of subsurface settlement trough and soil volume loss.Considering the ground loss variation with depth,three modes are deduced from the volumetric deformation responses of the soil above the tunnel crown.Then,analytical solutions to predict subsurface settlement for each mode are presented using stochastic medium theory.The influences of C/D,ηt and VBP on the key parameters (i.e. B and N) in the analytical expressions are discussed to determine the fitting formulae of B and N.Finally,the proposed analytical solutions are validated by the comparisons with the results of model test and numerical simulation.Results show that the fitting formulae provide a convenient and reliable way to evaluate the key parameters.Besides,the analytical solutions are reasonable and available in predicting the subsurface settlement induced by shield tunnelling in sandy cobble stratum.

1.Introduction

Sandy cobble stratum is a typical heterogeneous stratum composed of cobble,pebble,sand and a little clay.This stratum is widely distributed in Beijing,Chengdu,Lanzhou and other cities of China.Shield tunnelling in this stratum causes complicated ground movement,affecting not only the surface buildings but also the subsurface adjacent infrastructures (e.g.the existing tunnels,underground pipelines and pipe foundations).In general,the ground surface deformation caused by the excavation of new tunnels is a significant index to evaluate the safety of the surface buildings.However,due to the different magnitudes of ground movement at different depths,the ground surface deformation cannot accurately assess the safety of underground infrastructures.Therefore,it is necessary to investigate the characteristics of tunnelling-induced subsurface deformation and predict the magnitude of subsurface settlement.

For shield tunnels excavated in undrained clay,there is no plastic volumetric deformation (Li et al.,2019).The volume of subsurface settlement trough remains constant along depth,and soil volume loss at the ground surface is equivalent to tunnel volume loss (Marshall et al.,2012).Analytical solutions for surface settlement induced by ground loss in an elastic half-plane have been proposed by means of the virtual image technique(Sagaseta,1987;Verruijt and Booker,1996;Loganathan and Poulos,1998),complex variable method (Verruijt,1997,1998),Airy function(Bobet,2001;Chou and Bobet,2002;Park,2004,2005),and stochastic medium theory (Yang et al.,2004;Yang and Wang,2011).Based on the research of Peck (1969) settlement trough shape is always assumed to obey the Gaussian distribution function with two key parameters: the maximum settlementSmaxand the settlement trough widthi.According to the field monitoring data of ground settlement,it has been proved that both surface and subsurface settlement data can be well predicted by normal probability function (Fang et al.,2014).Besides,the subsurface settlement trough widthi(z) decreases with depthz,while the maximum settlementSmax(z) increases with depthz.Considering the data of field tests and centrifuge model tests in clays,Mair et al.(1993)revealed that the settlement trough widthi(z) was proportional toz0-z,wherez0is the depth of the tunnel axis.A linear relationship existed approximately between the trough widthi(z)and depthz,i.e.i(z)can be expressed asi(z)=i(0)-kz,wherei(0)=0.5z0andk=0.325 for clays.Lee et al.(1999)carried out a series of centrifuge model tests to investigate the ground movement around tunnels in soft soil and provided prediction formulae ofi(z) andSmax(z)for different tunnel volume loss rates.However,the results above are only available for tunnels in undrained clay.Shield tunnelling in drained soil (e.g.sand and granular soil) is also very common in engineering practice.

For shield tunnels in drained soil,the soil volume loss varies with depth due to the existence of tunnelling-induced volumetric deformation of the soil (Marshall et al.,2012).Laboratory model tests play an important role in studying ground movement due to ground loss.Chambon and Corte (1994) addressed the failure mechanism and limit support pressures in front of the face for shallow tunnels in cohesionless soil using centrifugal model tests.Yoo and Shin (2003) performed a series of reduced-scale model tests to investigate the effect of face reinforcement parameters on the deformation behavior of tunnel face.Shin et al.(2008) conducted a full-scale model test to reveal the effects of reinforcing patterns and pipe length on the ground deformation and face stability of tunnels in granular soil.Lee(2009)carried out three model tests of tunnels in granular soil with different cover depths to investigate the subsurface settlement profiles due to tunnel volume loss.It is found that the coefficientKof trough width fell within 0.15-0.35 for granular materials.Hu et al.(2020) reproduced the real construction of tunnels in cobble-rich soil using a reducedscale earth pressure balance (EPB) shield and analyzed the subsurface deformation.For tunnelling in sand,Wang et al.(2016)carried out a model test and concluded that the volume of subsurface settlement trough decreases significantly with the decrease of depth.Nevertheless,the results of centrifuge model tests(Marshall et al.,2012;Franza et al.,2019;Song and Marshall,2020)showed that the volume variation of subsurface settlement trough with depth is related to many factors(e.g.cover depth ratio,tunnel volume loss and soil density) and cannot be generalized.On the other hand,developing analytical solutions (Vu et al.,2015) and empirical methods (Wang,2021) is also significant in predicting and estimating the tunnelling-induced subsurface settlement.Considering the effects of tunnel volume loss,cover depth ratio and soil relative density,Franza and Marshall (2019) presented empirical and semi-analytical methods to evaluate the surface and subsurface ground movements in sand.Using stochastic medium theory,Wang et al.(2020) established a propagation model for subsurface settlement induced by tunnelling in soft soil.Based on the developed Gaussian function (Lu et al.,2020),Lin et al.(2021)proposed a concept of transmission ratio of ground loss to describe the variation of ground loss with depth.

For tunnelling in sandy cobble stratum,there have been some researches on the ground settlement through field monitoring data(He et al.,2012,2021;Zhang et al.,2013)and model tests(Fan and Zhang,2013;Wang et al.,2018;Hu et al.,2020).However,these researches only focused on the characteristics of surface or subsurface settlement for specific case histories,lack of a systematic study of ground movement characteristics.Besides,theoretical models to predict the subsurface settlement induced by tunnelling in sandy cobble stratum are rarely reported in the literature.Although Bai et al.(2013)and Lu et al.(2020)proposed prediction methods for ground settlement,the volumetric deformation response of the soil was not considered.Franza et al.(2019)pointed out that the propagation of ground movement in drained soil was closely related to the volumetric deformation response of the soil.Since sandy cobble stratum is a heterogeneous soil with high porosity,wide particle size distribution and significant soil arching effect (Huang et al.,2019),the variation of ground loss with depth and the volumetric deformation response of the soil are bound to be complicated.

The aim of this paper is to explore the volumetric deformation responses of sandy cobble soil and propose analytical solutions to predict subsurface settlement induced by shield tunnelling.Considering the ground loss variation with depth,the volumetric deformation modes are summarized from a series of numerical analyses.Meanwhile,the variations of subsurface settlement trough width,maximum subsurface settlement and soil volume loss rate with cover depth ratio,tunnel volume loss rate and volumetric block proportion are discussed.Then,analytical solutions for different volumetric deformation modes are derived using stochastic medium theory.Fitting formulae for the key parameters in the analytical solutions are provided based on the parametric studies.Finally,comparisons with the results of model test and numerical simulation are made to validate the analytical solutions.

2.Numerical analyses

A series of numerical analyses(about 240 groups)in plane strain condition was performed to understand the variation of ground movement (including subsurface settlement trough widthi(z),maximum subsurface settlementSmax(z) and soil volume loss rate ηs(z))with cover depth ratio(C/D),tunnel volume loss rate(ηt)and volumetric block proportion(VBP).CandDare the cover depth and diameter of a tunnel,respectively.C/D,ηtandVBPcharacterize the tunnel geometries,tunnel construction quality and geological conditions,respectively.They affect the level of confining pressure,magnitude of shear strain and content of coarse particles(Song and Marshall,2020).

2.1.Overall mechanical parameters for different VBP values

Due to the remarkable heterogeneity,sandy cobble soil was considered as a heterogeneous and two-component material composed of blocks and soil matrix (Zhang et al.,2018;Du et al.,2019).Here,“block” is a general term for gravels,cobbles and boulders with the particle sizes greater than 20 mm.The volumetric block proportion (VBP) has a remarkable impact on the overall mechanical properties of sandy cobble soil.Thus,five stochastic structural models representing the two-component material for differentVBPvalues (VBP=30%,40%,50%,60% and 70%)were built in plane strain condition,as shown in Fig.1.Here,the socalled“stochastic structural model”,which has been widely used in the research of soil-rock mixtures(You,2001;Li et al.,2004;Wang et al.,2005),emphasizes the significant influence of internal stochastic structure characteristics on the overall mechanical properties of sandy cobble soil.According to the results of field investigation,the common block sizes of sandy cobble soil are mainly in the range of 20-200 mm.When building stochastic structural models of sandy cobble soil,the minimum and maximum block sizes were set as 20 mm and 200 mm.Elliptical block shapes were assumed.In order to minimize the scale effect,the sizes of the numerical models were set as 1 m in width and 2 m in height.The position and orientation of blocks were generated randomly and placed without interpenetration through implementing a Python routine.The detailed generation procedure of stochastic structural models refers to Wang et al.(2023).

Fig.1.Stochastic structural models for different VBP values in plane strain condition: (a) VBP=30%,(b) VBP=50%,(c) VBP=70%,and (d) Representative blocks.

Using the respective mechanical parameters of soil matrix and blocks (Table 1) calibrated by Wang et al.(2023),five groups of numerical simulation in plane strain condition were performed to acquire the overall stress-strain curves for differentVBPvalues.Both soil matrix and blocks were assumed to follow the Mohr-Coulomb model (Xu et al.,2008;Napoli et al.,2018).Normal stress boundary conditions(i.e.confining pressure σ3)were applied to the lateral boundaries of the model.The bottom of the model was fixed in the normal direction.When performing numerical simulation,the stress equal to the confining pressure was applied to the top of the model at first.Then,the displacement boundary condition was applied for strain control loading.The overall stress-strain curves for differentVBPvalues are shown in Fig.2(see the scatters).

Table 1Respective mechanical parameters of soil matrix and blocks (Wang et al.,2023).

Fig.2.Comparisons of stress-strain curves.

Since the overall stress-strain behavior of sandy cobble soil shows a noticeable strain-hardening characteristic (Wang et al.,2023),the Hardening-Soil model proposed by Schanz et al.(1999)is adopted to describe this overall strain-hardening mechanical behavior (He et al.,2021).The Hardening-Soil model was formulated in an elastoplastic framework with stress-dependent stiffness modules based on the hyperbolic relationship between the axial strain and the deviatoric stress,overcoming the restrictions of the purely nonlinear elastic Duncan-Chang model.The parameters of the Hardening-Soil model include the failure parameters (i.e.friction angle φ,cohesioncand dilation angle ψ),the basic parameters for stiffness (i.e.reference secant stiffness modulus,reference tangent stiffness modulusand power for stress-level dependency of stiffnessm),and the advanced parameters (i.e.reference Young’s modulus for unloading and reloading,Poisson’s ratio for unloading-reloading μur,reference stress for stiffnesspref,K0-value for normal consolidationand failure ratioRf).The process of parameter inversion refers to Schanz et al.(1999) and Itasca (2017).The details are described as follows.

The calibrations of φ andchave no difference from the Mohr-Coulomb model.Using the stress-strain data under various confining pressures,the relationship betweenq=(σ1-σ3)/2 andp=(σ1+σ3)/2 is plotted to fit the Mohr-Coulomb envelope by a trend line.The slope of the trend line is sinφ and the intercept isccosφ.Due to the little clay content in sandy cobble soil,the soil is assumed to be cohesionless.The dilation angle ψ can be estimated through the maximum slopesmof the volumetric strain εvvs.axial strain ε1curve,i.e.sm=-Δεv/Δε1=2sinψ/(1-sinψ).However,for a sandy cobble stratum with two-component material,smis not easy to determine.Besides,at present,there is no clear and accurate method to calculate the dilation angle.Thus,ψ is obtained by trial and error so that the stress-strain curve develops stably and smoothly without fluctuations.ψ increases with the increase ofVBP.In addition,the densities ρ for differentVBPare calculated by the formula (Du et al.,2019): ρ=ρblockVBP+ρmatrix(1-VBP),where ρblockand ρmatrixare the densities for blocks and soil matrix,respectively.

Generally,the advanced parameters have default values and should not be necessary to specify.prefis commonly assumed as the standard atmospheric pressure (i.e.100 kPa).μurcan be estimated from the unloading-reloading slope of the εvvs.ε1curve (Itasca,2017).Schanz et al.(1999) pointed out that the realistic value of μuris about 0.2.Although μurmakes a difference on σ1-σ3and εv,the influence is slightly and not considered here.In the Hardening-Soil model,μuris assumed to be a constant value with a typical value of 0.2 (Itasca,2017).is a proper input parameter with a realistic correlation=1-sinφ as default.Rfis fixed to be 0.99.must be greater thanand use a value in the range of(3-5)for most soils.It is generally assumed to be.

Through parameter inversion above,using the Hardening-Soil model,the overall mechanical parameters of sandy cobble stratum for differentVBPare obtained in Table 2.Fig.2 shows the comparison of stress-strain curves between the homogeneous Hardening-Soil model and stochastic structural model under σ3=400 kPa.It can be seen that the stress-strain curves obtained by the Hardening-Soil model are in good agreement with the results of stochastic structural models,demonstrating the reasonability of the overall mechanical parameters listed in Table 2.

Table 2Overall mechanical parameters of sandy cobble stratum for different VBP values.

Note that,this paper focuses on the transverse ground movement induced by shield tunnels in sandy cobble stratum.Thus,all the numerical analyses are performed in two-dimensional (2D)plane strain condition.However,in three-dimensional (3D) space,the overlapping detection of blocks is more complex,especially for ellipsoidal blocks with highVBP(Wang et al.,2023).Although Wang et al.(2023)proposed a numerical modelling method for 3D condition,admittedly,the modelling efficiency is low and building a numerical model withVBP=70% is rarely possible.To solve this problem,a feasible alternative scheme may be that only the blocks with large sizes are built in the 3D numerical model,while the blocks with small sizes(e.g.smaller than 60 mm)are homogenized into the soil matrix(i.e.equivalent soil matrix).Simultaneously,the strength of soil matrix is improved due to the homogenization of small blocks.The mechanical behavior of equivalent soil matrix can also be characterized by the Hardening-Soil constitutive model proposed by Schanz et al.(1999).Similar to the 2D condition,the mechanical parameters of equivalent soil matrix can be calibrated through parameter inversion.Then,based on the mechanical parameters of equivalent soil matrix,the numerical modelling with highVBPin 3D condition can be realized by changing the content of large blocks.

2.2.Numerical models and boundary conditions

Fig.3 shows the numerical model withC/D=2 andD=6 m built by the finite difference software FLAC3D (Itasca,2017).Due to the symmetry,only half of the model was simulated to improve computational efficiency.The model consists of 9400 nodes and 4588 zones.The sizes of the model were set as 30 m in width and 33 m in height to reduce the boundary effects.The thickness of the model was set as 0.1 m to simulate the plane strain condition.The lateral displacement boundaries were fixed at the normal direction,while the bottom displacement boundary was fixed at both horizontal (x) and vertical (z) directions.The ground surface was free.

Fig.3.Numerical model for C/D=2.

Due to the factors such as physical clearance between the excavation face and external boundary of lining,3D elastoplastic deformation ahead of tunnel face,and over-excavation caused by workmanship problems,a gap occurs between the excavation face around the tunnel circumference and external boundary of lining(Lee et al.,1992).Considering the grouting technique,the tunnelling-induced gap and volume loss are reduced but still exist.Here,a typical non-uniform deformation pattern of the tunnel wall was adopted without considering the ovalization (Sagaseta,1998),as shown in Fig.4.The tunnel volume lossVtwas simulated by applying a displacement boundary condition around the tunnel circumference (Zhang et al.,2016).The reduced gapg(θ) can be calculated by Eq.(1).The maximum gapgmaxat the tunnel crown is related to the tunnel volume lossVtand calculated by Eq.(2).

Fig.4.Schematic diagram of ground settlement induced by tunnel volume loss.

whereRandrrepresent the radii of excavation face and external boundary of lining,respectively;and θ is the angle of a point at excavation face in the polar coordinate system.The tunnel volume loss rate ηtis denoted as

It should be pointed out that only the classical deformation pattern of the tunnel wall proposed by Sagaseta (1998) is considered in this paper.The lining moves downward due to the selfweight.In reality,field measurements in some cases showed that the tunnel lining will move upward due to the grouting pressure and buoyancy forces (Bezuijen et al.,2004).This causes another deformation pattern of the tunnel wall (i.e.another boundary condition).When the tunnel geometries,tunnel construction quality and geological conditions are constant,different deformation patterns of the tunnel wall may lead to different volumetric deformation responses of the soil above the tunnel crown.Here,this study only discusses the effects ofC/D,ηtandVBPon the ground movement.The effect of the deformation pattern of the tunnel wall on the ground movement needs further research.

2.3.Numerical results

Using the material parameters in Table 2 and the numerical model in Fig.3,a series of numerical analyses (about 240 groups)was performed to investigate the variation of tunnelling-induced ground movement with depth and discuss the volumetric deformation responses of sandy cobble soil above the tunnel crown.The numerical results are described as follows.

2.3.1.Subsurface settlement trough width i(z)

Fig.5 shows the effects of tunnel geometries (CandD) on the subsurface settlement trough width when ηt=2% andVBP=50%.Trough widthidecreases nonlinearly with depthz.At a certain depthz,trough widthiincreases with tunnel cover depthCor tunnel diameterD.Whereas,as shown in Fig.5c,dimensionless variablesi/DandC/Dare always positively correlated.Besides,althoughCandDvary,dimensionless variablei/Dkeeps constant as long asC/Dis constant.

Fig.5.Effects of tunnel geometries on the subsurface settlement trough width:(a)Variation of i with z for different C,(b)Variation of i with z for different D,and(c)Variation of i/D with z/C for different C/D.

The normalized trough widthi/z0is plotted against the relative depthz/z0in Fig.6 for differentVBPvalues.It can be seen that the normalized trough width at any depth decreases with ηtwhenC/Dis constant,which is consistent with the results of centrifuge tests by Bezuijen and Van Seters(2002).For shallow tunnels(Fig.6a and c),the higher theVBP,the smaller the trough width.Besides,the effect ofVBPis more apparent in the case of high ηt.For deep tunnels(Fig.6b and d),the influence ofVBPon trough width is not so apparent,especially for low ηt.It is worth noting in Fig.6d that the trough width far away from the tunnel crown in the case of highVBPis larger than that in the case of lowVBP.This can be explained by the soil arching effect (Franza et al.,2019).For the case of deep tunnels with highVBP,the height of the soil arch is lower than the case of deep tunnels with lowVBP.Ground deformation above the soil arch is smaller in magnitude but propagates upwards and outwards,causing a wider settlement trough.On the contrary,ground settlement beneath the soil arch is relatively larger in magnitude but narrower in width.

Fig.6.Normalized trough width i/z0 with z/z0 for different VBP: (a) C/D=1,ηt=0.5%;(b) C/D=4,ηt=0.5%;(c) C/D=1,ηt=4%;and (d) C/D=4,ηt=4%.

2.3.2.Width coefficient of subsurface settlement trough K(z)

Mair et al.(1993) investigated the field and centrifuge model test data of subsurface settlement above tunnels in clay and revealed that the settlement trough widthi(z)was proportional toz0-z.i(z) can be expressed as

whereKis the width coefficient of settlement trough.It was estimated,in clay,as

For tunnelling in sandy cobble soil,Fig.7 illustrates the variation of width coefficientK(z)of settlement trough with relative depthz/z0for different ηt.It is found thatK(z)decreases with the increase of ηtat any depth.For shallow tunnels (Fig.7a and c),the higher theVBP,the smaller the trough width coefficientK(z).For deep tunnels(Fig.7b and d),the influence ofVBPis not apparent.It is worth noting in Fig.7d thatK(z) increases firstly,and then decreases sharply with the increase ofz/z0when ηt2%.This also can be explained by the soil arching effect.For the case of deep tunnels with highVBP,when ηt2%,the soil arch is fully mobilized.Ground settlement beneath the soil arch is relatively larger in magnitude but narrower in width than that above the soil arch.Besides,the decrease of trough width is more radical than the increase of depth(Fig.6d).

Fig.7.Width coefficient K of settlement trough with z/z0 for different ηt: (a) VBP=30%, C/D=1;(b) VBP=30%, C/D=4;(c) VBP=70%, C/D=1;and (d) VBP=70%, C/D=4.

Marshall et al.(2012)extended the expression to sandy soils andK(z) was denoted as

whereK(0) is the width coefficient of settlement trough at the ground surface and ?i/?zis the variation rate of settlement trough width with depth.

Fitting the data of settlement trough width obtained from the numerical simulation in this section with Eq.(6),the variations of R-square with ηtfor differentC/DandVBPare depicted in Fig.8.It can be seen that all the values of R-square in the case of lowVBPare greater than 0.9,indicating that the subsurface settlement trough width in the case of lowVBPincreases approximately linearly with depth.With the increase ofVBP,the values of R-square for shallow tunnels are still greater than 0.9,while those for deep tunnels are greater than 0.9 only for low ηt.For the case of highVBP,apart from ηt=0.5%,the values of R-square are basically smaller than 0.9.Besides,the greater theC/Dand ηt,the smaller the values of Rsquare.Therefore,lower values ofVBP,C/Dand ηtfacilitate the linear variation of subsurface settlement trough width with depth.

Fig.8.Variations of R-square with ηt for different C/D and VBP: (a) VBP=30%,(b)VBP=50%,and (c) VBP=70%.

Figs.9 and 10 plot the variations ofK(0) and ?i/?zwith ηtfor differentC/DandVBP.It is found in Fig.9 that the range ofK(0)for sandy cobble stratum is 0.4-0.7,which is quite different from the typical range of 0.25-0.45 for sand or gravel stratum surveyed by Mair and Taylor(1997).Besides,the lower the ηt,or the higher theC/DandVBP(except forC/D=1),the greater theK(0).This indicates that the surface settlement in sandy cobble stratum for low ηt,especially for highC/DandVBP,is smaller in magnitude but wider in width than that in sand or gravel stratum.The soil behaviour at the ground surface may be approximately elastic (Marshall et al.,2012).In addition,in Fig.10,the range of ?i/?zfor sandy cobble stratum is-0.6 to-0.2.Besides,the lower the ηt(except forVBP=70%),or the higher theC/DandVBP,the smaller the ?i/?z,i.e.the variation of subsurface settlement trough width with depth is more sharply.

Fig.9.Variation of K(0)with ηt for different C/D and VBP:(a)VBP=30%,(b)VBP=50%,and (c) VBP=70%.

Fig.10.Variation of ?i/?z with ηt for different C/D and VBP: (a) VBP=30%,(b)VBP=50%,and (c) VBP=70%.

2.3.3.Subsurface maximum subsurface settlement Smax(z)

Fig.11 shows the effects of tunnel geometries (CandD) on the subsurface maximum settlementSmaxwhen ηt=2% andVBP=50%.It can be seen that maximum settlementSmaxincreases nonlinearly with depthz.At a certain depthz,the maximum settlementSmaxincreases with tunnel diameterDbut decreases with tunnel cover depthC.As shown in Fig.11c,the normalized maximum settlementSmax/Dis always negatively correlated withC/D.Besides,althoughCandDvary,Smax/Dkeeps constant basically as long asC/Dis constant.

The effects of ηtandVBPon the normalized maximum settlementSmax/Dat different relative depthsz/z0are displayed in Fig.12.It is found thatSmax/Dincreases with ηtbut decreases withVBPat any depth.Fig.13 compares the relationship betweenSmax/Dand ηtfor differentC/DandVBP.The values of relative depthsz/z0=0 and 0.5 are taken as a contrast.It can be seen that the increase of subsurface normalized maximum settlement is more sharply than that at the ground surface.In addition,the maximum settlement for shallow tunnels is more sensitive to the variation of ηtthan that for deep tunnels.Similarly,the maximum settlement for lowVBPincreases more sharply with ηtthan that for highVBP.

Fig.12.Effects of tunnel volume loss rate ηt and VBP on the normalized maximum settlement: (a) C/D=2, VBP=50%;and (b) C/D=2,ηt=2%.

2.3.4.Soil volume loss rateηs(z)

For tunnelling in undrained clay,the volume of the soil is constant and soil volume loss does not vary with depth (i.e.ηs(0)=ηs(z)=ηt).Whereas,for tunnelling in drained soil,the volumetric contraction or dilation of the soil occurs due to shear,meaning that soil volume loss varies with depth (i.e.ηs(0) ≠ηs(z) ≠ηt).

The expression ηs(z)/ηtis defined to determine whether the overall response of volumetric deformation of the soil at a certain depthzis contractive,constant or dilative.Note that the overall response of volumetric deformation refers to the cumulative volumetric deformation of the soil beneath the given depthz.If ηs(z)/ηt＞1,the overall response of volumetric deformation at the given depth is contractive.If ηs(z)/ηt=1,the overall response is constant.If ηs(z)/ηt＜1,the overall response is dilative.In addition,the first derivative of ηs(z) (i.e.η′s(z)) is defined to determine the localized response of volumetric deformation of the soil at a certain depthz.Here,the localized response refers to the volumetric deformation of an infinitesimal soil layer with the thickness of dzat the given depthz.If η′s(z)＞0,i.e.ηs(z)is positively correlated withz,the localized response of volumetric deformation at the given depth is dilative.If η′s(z)=0,i.e.ηs(z) does not vary withz,the localized response is constant.If η′s(z) ＜0,i.e.ηs(z) is negatively correlated withz,the localized response is contractive.The schematic diagram of volumetric deformation responses induced by tunnelling is displayed in Fig.14.

Fig.14.Schematic diagram of volumetric deformation responses induced by tunnelling.

Fig.15 shows the variation of ηs(z)/ηtwithz/Dfor differentC/D.It is found that for a givenVBPand ηt,the greater theC/D,the higher the ηsat any depth.The soil tends to have a more contractive(less dilative) response.This can be explained by the higher level of confining pressure(i.e.a deeper tunnel)or the lower magnitude of shear strain around the tunnel circumference(i.e.a smaller tunnel).In addition,for the case of lowVBPand ηt(Fig.15a),both the overall and localized volumetric deformation responses of the soil are contractive at any depth,regardless of shallow or deep tunnels.As theVBPor ηtincreases (i.e.higher content of coarse particles or higher magnitude of shear strain around the tunnel circumference),the soil above shallower tunnels tends to have a more dilative response.For the case of highVBPand ηt(Fig.15d),both the overall and localized volumetric deformation responses of the soil above shallow tunnels are dilative at any depth,while only the overall response of the soil above deep tunnels is dilative at any depth.Note that,for deep tunnels in Fig.15d,the soil close to the tunnel crown experiences considerably localized dilation(i.e.the value of η′s(z) is great).Therefore,although the soil away from the tunnel crown experiences localized contraction,it still indicates an overall dilative response due to the cumulative volumetric deformation of the soil beneath the calculated level.These results are consistent with the results of centrifuge model tests by Song and Marshall(2020).

Fig.16 describes the effect ofVBPon the variation of ηs(z)/ηtwithz/z0.It is found that for a givenC/Dand ηt,the lower theVBP,the higher the ηsat any depth.The soil tends to have a more contractive(less dilative)response.This can be explained by the lower content of coarse particles.In addition,for the case of deep tunnels with low ηt(Fig.16c),both the overall and localized volumetric deformation responses of the soil are contractive at any depth,regardless of low or highVBP.As theC/Ddecreases or as the ηtincreases (i.e.lower confining pressure or higher magnitude of shear strain around the tunnel circumference),the soil with higherVBPtends to have a more dilative response.For the case of shallow tunnels with high ηt(Fig.16b),both the overall and localized volumetric deformation responses of the soil with highVBPare dilative at any depth,while only the overall response of the soil with lowVBPis dilative at any depth.Note that,for the case of lowVBPin Fig.16b,the soil close to the tunnel crown experiences considerably localized dilation (i.e.the value of η′s(z)is great).Therefore,although the soil away from the tunnel crown experiences localized contraction,it still indicates an overall dilative response due to the cumulative volumetric deformation of the soil beneath the calculated level.

Fig.16.Effect of VBP on the variation of ηs(z)/ηt with z/z0: (a) C/D=1,ηt=0.5%;(b) C/D=1,ηt=4%;(c) C/D=4,ηt=0.5%;and (d) C/D=4,ηt=4%.

To sum up,a deeper and smaller tunnel with lower tunnel volume loss or lower volumetric block proportion tends to result in a more contractive(less dilative)response of the sandy cobble soil.On the contrary,a more dilative (less contractive) response of the sandy cobble soil needs a shallower and larger tunnel with higher tunnel volume loss or higher volumetric block proportion.

2.4.Volumetric deformation modes of sandy cobble soil

According to the localized response of volumetric deformation analyzed in Section 2.3.3,three volumetric deformation modes of sandy cobble soil are deduced,as shown in Fig.17:

Fig.17.Three volumetric deformation modes of the sandy cobble soil.

(1) Mode I: both the overall and localized volumetric deformation responses are contractive at any depth.

(2) Mode II is divided into two sub-modes according to the overall response.Mode II1: the overall response of the soil away from the tunnel crown is contractive,while that close to the tunnel crown is dilative.Mode II2:the overall response at any depth is dilative due to the cumulative effect of the volumetric deformation.For the localized responses of modes II1and II2,only a small region of the soil close to the tunnel crown has a considerable localized dilative response,while the rest of the soil has a localized contractive response.

(3) Mode III: both the overall and localized volumetric deformation responses are dilative at any depth.

Through a series of numerical analyses,the volumetric deformation modes for different conditions(i.e.differentC/D,ηtandVBP)are listed in Table 3 as a reference.With the variations ofC/D,ηtandVBP,three volumetric deformation modes can evolve with each other.

Table 3Volumetric deformation modes for different conditions.

3.Analytical solutions

Ground movement induced by tunnelling in heterogeneous sandy cobble stratum is so complicated that it is hard to describe it with classical elastic or elastoplastic mechanics.Thanks to the stochastic medium theory (Litwiniszyn,1957),the ground deformation induced by ground loss can be studied from the perspective of probability and statistics (Yang et al.,2004;Yang and Wang,2011).In stochastic medium theory,the soil is regarded as a stochastic medium,and the tunnelling-induced ground movement is regarded as a stochastic process.For the horizontally laminated medium in 2D condition,the governing equation is given as follows:

whereS(x,z)is the settlement of a point(x,z)beneath the ground surface;B(z) andN(z) are two generalized parameters characterizing the ground movement,taking into account the lack of homogeneity of the medium;B(z) characterizes the vertical propagation capacity of subsurface settlement trough;andN(z)defines the volumetric deformation performance (contraction or dilation rate)of the soil at a certain depthzabove tunnels.

Suppose that the subsurface settlement trough shape at any depthzis described by Gaussian function (Mair et al.,1993;Fang et al.,2014),S(x,z) is written as

The soil volume lossVs(z) at any depthzis calculated as

The soil volume lossVs(0) at the ground surface is denoted as

Combining Eqs.(8)-(10) gives

Substituting Eq.(11) into Eq.(7)and solving it gives

According to the numerical analyses in Section 2.3.1,the subsurface settlement trough widthi(z) always decreases with depth.Therefore,the relation ofB(z)＜0 is always true.Besides,the larger theB(z),the narrower and deeper the settlement trough width.Meanwhile,it is found from the numerical analyses in Section 2.3.4 that the localized volumetric deformation response of the soil is contractive whenN(z)＜0,but dilative whenN(z)＞0.Besides,the larger the absolute value ofN(z),the higher the volumetric deformation performance (contraction or dilation rate)of the soil.

Integrating Eq.(12),the settlement trough widthi(z) and soil volume lossVs(z) at the given depthzcan be calculated as

whereC1andC2are the constants and can be solved by the boundary conditions.

For simplicity,the soil is assumed to be a homogeneous and isotropic material to facilitate the exact analytical solutions,i.e.bothB(z)andN(z)are assumed to be constant and do not vary with the depthz.Eq.(13) can be simplified as

Next,the analytical solutions to predict the subsurface settlement for three volumetric deformation modes of sandy cobble soil are derived according to Eq.(14).

3.1.Modes I and III

For mode I and III,the localized volumetric deformation responses of the soil above the tunnel crown are contractive and dilative at any depth,respectively.

In general,it is easy to monitor the settlement data at the ground surface.Fitting the measured data with Eq.(8) and combining Eq.(10),the boundary conditions at the ground surface are obtained as follows: whenz=0,i(z)=i(0) andVs(z)=Vs(0).Substituting the boundary conditions into Eq.(14) gives

Substituting Eq.(15) into Eq.(11),the analytical solution of subsurface settlement for modes I and III is denoted as

This formula is consistent with that obtained by Wang et al.(2020) using the Fourier transform.Besides,for mode I:＜0 and＜0;for model III:＜0 and＞0.Here,the subscripts of the parametersBandNare only used to distinguish the volumetric deformation mode,while the superscripts ‘c’ and ‘d’represent the localized contractive and dilative responses of the soil,respectively.

For a project in engineering practice,if the observed data of subsurface settlement at a certain depth are available,the parametersBandNcan be obtained by back analysis according to Eq.(16).Then,the subsurface settlement at any depth can be predicted.

3.2.Mode II

For mode II,the localized volumetric deformation response of the soil close to the tunnel crown is dilative,while that close to ground surface is contractive.

The solution for the region with contractive response is identical with that for mode I.Besides,＜0 and＜0.

Suppose thatzcdenotes the critical depth between the regions with localized contractive and dilative response,the boundary conditions for the region with localized dilative response are as follows:whenz=zc,i(z)=i(zc)andVs(z)=Vs(zc).Substituting the boundary conditions into Eq.(14) gives

whereandare the parameters for the region with localized dilative response.Vs(zc) is calculated as

wherei(zc) andSmax(zc) can be calculated by the solution for the region with localized contractive response due to the displacement continuity.Combining Eqs.(8),(9),(17) and (18),the analytical solution of subsurface settlement for the region with localized dilative response is denoted as

where the parametersandcan be obtained by back analysis of the observed data of subsurface settlement in the region with localized dilative response.Besides,＜0 and＞0.Once the key parameters are obtained,the subsurface settlement at any depth can be predicted.

In Eq.(19),the critical depthzcis unknown.It may be relative to the cover depth ratio (C/D),tunnel volume loss rate (ηt) and volumetric block proportion (VBP).Through the fitting analyses of the numerical results in Section 2.3,it is found that whenVBPis constant,dimensionless variablezc/z0follows the relationship below:

wherea,bandcare constants and only related to theVBP.The fitting results are listed in Table 4.All the values of R-square are greater than 0.9,indicating the reliability and availability of Eq.(20).

Table 4Fitting results of a, b and c for different VBP.

Table 5Fitting formulae and coefficients of B/D for different VBP.

Table 6Fitting formulae and coefficients of Nz0 for different VBP.

4.Parametric study

If there is no available observed data of subsurface settlement,the key parameters(i.e.BandN)in the analytical solutions cannot be determined by back analysis,resulting in that the subsurface settlement at any depth cannot be predicted using the analytical solutions in Section 3.Therefore,the purpose of this section is to examine the influences ofC/D,ηtandVBPon the key parameters and study the fitting formulae forBandN.The values ofBandNare obtained by back analysis of the numerical results in Section 2.3 according to the analytical solutions in Section 3.

Dimensionless variablesB/DandNz0are taken to investigate the influences ofC/Dand ηtfor differentVBP,as shown in Figs.18 and 19.In Figs.18a-c and 19a-c,the upper part of each graph indicates a localized dilative volumetric deformation response of the soil,while the lower part indicates a localized contractive response of the soil.For the givenC/Dand ηt,if the data only appear in the upper part of the graph,the volumetric deformation mode of the soil belongs to mode III;if the data only appear in the lower part of the graph,the volumetric deformation mode of the soil belongs to mode I;if the data appear both in the upper and low parts of the graph,the volumetric deformation mode of the soil belongs to mode II.Figs.18d and e and 19d and e select some representative cases to display the variation relationship ofB,NwithC/D,ηtmore visually.

4.1.Propagation capacity (B) of settlement trough

As mentioned early,the parameterBcharacterizes the propagation capacity of subsurface settlement trough above the tunnel crown.Besides,the larger the value ofB,the narrower and deeper the settlement trough width.

Therefore,it is found from Fig.18 that for the case of soil contraction,Bc/Ddecreases approximately exponentially withC/Dbut increases approximately linearly with ηt.For the case of soil dilation,Bd/Dincreases firstly (only forVBP=30% and 50%),and then decreases approximately linearly withC/D.In addition,Bd/Dincreases approximately linearly with ηt.

Fig.18.Influences of C/D and ηt on the key parameters B for different VBP: (a) VBP=30%,(b) VBP=50%,(c) VBP=70%,(d) Bc/D vs. C/D and ηt,and (e) Bd/D vs. C/D and ηt.

4.2.Volumetric deformation performance (N) of sandy cobble soil

As mentioned early,the parameterNcharacterizes the volumetric deformation performance (contraction or dilation rate) of the sandy cobble soil above the tunnel crown.Besides,the larger the absolute value ofN,the higher the volumetric deformation performance (contraction or dilation rate)of the soil.

Therefore,it is found from Fig.19 that for the case of soil contraction(Nz0＜0),Ncz0decreases approximately linearly withC/Dbut increases approximately linearly with ηt.In addition,the higher theVBP,the larger the value ofNcz0,but the smaller the absolute value ofNcz0,i.e.the volumetric contraction performance of sandy cobble soil with lowVBPis greater than that with highVBP.For the case of soil dilation (Nz0＞0),Ndz0decreases firstly (for mode III),but then increases approximately linearly(for the dilative region of mode II)withC/D.The reason may be that the range of the dilative region of mode II decreases withC/D,resulting in that the dilation rate increases.In addition,Ndz0increases approximately linearly with ηt.The higher theVBP,the larger the value ofNdz0,i.e.the volumetric dilation performance of sandy cobble soil with highVBPis greater than that with lowVBP.

Fig.19.Influences of C/D and ηt on the key parameters N for different VBP: (a) VBP=30%,(b) VBP=50%,(c) VBP=70%,(d) Ncz0 vs. C/D and ηt,and (e) Ndz0 vs. C/D and ηt.

4.3.Fitting formulae for B and N

According to the variations ofBandNwithC/D,ηtandVBP,the fitting formulae and the coefficients for the three volumetric deformation modes are displayed in Tables 5 and 6.It can be seen that most of the values of R-square are greater than 0.9,indicating that the fitting formulae are reliable and available for the determination ofBandN.For a shield tunnel in sandy cobble stratum,it should be pointed out that back analysis is still the preferred way to determine the key parametersBandNas long as the monitoring data of subsurface settlement are available.Otherwise,if the monitoring data of subsurface settlement are unavailable(i.e.BandNcannot be obtained by back analysis),the fitting formulae provide a convenient way to determine the values ofBandN.

5.Comparison and validation

Limited by the complexity of measurement and the high cost of instrument installation,the field monitoring data of subsurface settlement induced by shield tunnelling in sandy cobble stratum are seldom found in engineering practice (Wang et al.,2016).In light of this,comparisons with the results of model test and numerical simulation are made to validate the proposed analytical solutions in this paper.In addition,in Eq.(16),ifN=0,the analytical solution for sandy cobble soil is reduced to that for undrained clays.Three field cases are taken to compare and validate the analytical solution in undrained clays.

5.1.Comparisons with model test results

5.1.1.Model tests conducted by Jiang (2014)

Using the developed EPB shield test system,Jiang(2014)carried out model tests to study the ground movement induced by tunnelling in sandy cobble stratum.The diameter and cover depth of the model tunnel are 800 mm and 1.2 m,respectively.Surface settlement (z=0) and subsurface settlement at the depth ofz=0.4 m were monitored.

The values ofi(z),Smax(z) and ηs(z) at the depths ofz=0 and 0.4 m are calculated by back analyses of the observed data as follows: whenz=0,i(0)=0.769 m,Smax(0)=4.016 mm and ηs(0)=1.54%;whenz=0.4 m,i(z)=0.505 m,Smax(z)=4.109 mm and ηs(z)=1.03%.Therefore,the localized volumetric deformation response of sandy cobble soil at the depths from 0 to 0.4 m is contractive.The analytical solution for mode I is adopted here.Due to the available observed data of subsurface settlement,back analysis is used to obtain the values ofand.Then,the prediction data of surface and subsurface settlement can be evaluated using Eq.(16).Fig.20 shows the comparisons between the predicted and observed data for surface and subsurface settlement,indicating that the analytical solution for mode I is reasonable and available.

Fig.20.Comparisons with the results of Jiang (2014).

5.1.2.Model tests conducted by Lee (2009)

Lee (2009) carried out 2D model tests of tunnels in granular mass to investigate the subsurface deformation.The diameter of the model tunnel is 100 mm.The depth of the tunnel axis is 270 mm(i.e.C/D=2.2).The tunnel volume loss rates are 5.83% and 10.94%.Subsurface settlements at the depths ofz=20 mm and 170 mm were monitored.

The values ofi(z),Smax(z) and ηs(z) at the depths ofz=20 mm and 170 mm are calculated by back analyses of the observed data as follows.For the case of ηt=5.83%,whenz=20 mm,i(z)=82.52 mm,Smax(z)=0.86 mm and ηs(z)=2.26%;whenz=170 mm,i(z)=58.61 mm,Smax(z)=1.63 mm and ηs(z)=3.05%.Looking up Table 3 shows that the volumetric deformation mode of this case may belong to mode II or III.If mode II is considered,the critical depthzcof localized contractive and dilative regions is about 170 mm according to Eq.(20).Then,the localized volumetric deformation response of the soil at the depths from 20 mm to 170 mm should be contractive.However,it is obvious that ηs(z)increases with the depths from 20 mm to 170 mm,indicating a localized dilative response.Therefore,mode III should be considered for the case of ηt=5.83%.For the case of ηt=10.94%,whenz=20 mm,i(z)=79.4 mm,Smax(z)=1.67 mm and ηs(z)=4.24%;whenz=170 mm,i(z)=54.4 mm,Smax(z)=3.47 mm and ηs(z)=6.03%.ηs(z) increases with the depths from 20 mm to 170 mm,indicating a localized dilative response.Besides,looking up Table 3 shows that the volumetric deformation mode of this case belongs to mode III.

Due to the available observed data of subsurface settlement,the values of,,i(0) andSmax(0) can be obtained through back analysis.Then,the prediction data are evaluated using Eq.(16).The comparisons between the predicted and observed data are displayed in Fig.21,indicating the reasonability of the analytical solution for mode III.In addition,the predicted values of settlement trough width at the depth ofz=70 mm were 75.39 mm and 72.04 mm for the case of ηt=5.83% and 10.94%,respectively,which are also in good agreement with the observed data (i.e.75.5 mm and 70.85 mm).

Fig.21.Comparisons with the model test results of Lee (2009): (a) ηt=5.83%,and (b) ηt=10.94%.

5.2.Comparisons with numerical results

5.2.1.Numerical simulation performed by Lee (2009)

In order to investigate the subsurface settlement profiles induced by shallow and deep tunnels in granular mass and make a comparison with the laboratory tests,Lee(2009)also performed a series of finite element analyses.The tunnel diameter is 100 mm.The depths of tunnel axis are 270 mm and 470 mm (i.e.C/D=2.2 and 4.2),corresponding to the shallow and deep tunnels,respectively.Comparison and validation for the shallow tunnel(C/D=2.2)have been carried out in section 5.1.2.Here,only the case of the deep tunnel (C/D=4.2) with ηt=5.83% is compared.Subsurface settlements at the depths ofz=70 mm and 370 mm were monitored.

The values ofi(z),Smax(z) and ηs(z) at the depths ofz=70 mm and 370 mm are calculated by back analyses of the observed data as follows:whenz=70 mm,i(z)=183.97 mm,Smax(z)=0.42 mm and ηs(z)=2.48%;whenz=370 mm,i(z)=96.33 mm,Smax(z)=0.78 mm and ηs(z)=2.34%.Looking up Table 3 shows that the volumetric deformation mode of this case belongs to mode II.According to Eq.(20),the critical depthzcof contractive and dilative regions is larger than 370 mm,regardless of theVBP.Besides,ηs(z)decreases with the depths from 70 mm to 370 mm.Therefore,the soil at the depths from 70 mm to 370 mm has a localized contractive response.Through back analysis,the values of,,i(0) andSmax(0) are obtained,and then the prediction data are evaluated using Eq.(16).The predicted data are compared with the observed data in Fig.22,showing good prediction performance.In addition,the predicted values of settlement trough width were 159.78 mm and 131.2 mm at the depth ofz=170 mm and 270 mm,respectively,also matching well with the observed data (i.e.161.08 mm and 138.39 mm).

Fig.22.Comparisons with the numerical results of Lee (2009).

5.2.2.Numerical simulation performed in this paper

Three cases,as shown in Table 7,are selected from the numerical simulation in Section 2 to validate the proposed analytical solutions and fitting formulae ofBandNfor modes I,II and III.The tunnel diameter is 6 m.Since the values ofC/D,ηt,VBPand the monitoring data of subsurface settlement are all available,the values ofBandNcan be obtained not only by back analysis but also by the fitting formulae in Tables 5 and 6.The corresponding depths of monitoring data taken to perform back analysis are also listed in Table 7.

Table 7Three cases to validate the analytical solutions and the fitting formulae.

Fig.23 compares the observed data with the predicted data by back analysis and fitting formula at different depths.As can be seen,almost all the predicted data at different depths agree well with the observed data,indicating that the proposed analytical solutions to predict subsurface settlement induced by shield tunnelling in sandy cobble stratum are reasonable and available.Meanwhile,the fitting formulae of the key parametersBandNin the analytical solutions are also validated to be reliable,providing a convenient way to determine the key parameters.

Fig.23.Comparisons with the numerical results in this paper: (a) Case 1,(b) Case 2,and (c) Case 3.

Fig.24.Comparisons between the observed and predicted surface settlement troughs: (a) Case 4,(b) Case 5,and (c) Case 6.

Fig.25.Comparisons between the observed and predicted subsurface maximum settlement above the tunnel centerline: (a) Case 4,(b) Case 5,and (c) Case 6.

5.3.Comparisons with the data in undrained clays

A total of three case histories in undrained clays are selected in this section to compare and validate the proposed analytical solution in engineering practice.Table 8 shows the tunnel details and soil types.The observed data include the surface settlement troughs and the subsurface maximum settlement above the tunnel centerline.All the observed data are taken from the published literature by Loganathan and Poulos (1998).

Table 8Case histories in undrained clays.

For tunnelling in undrained clay,the soil volume loss ηs(z)does not vary with depthz,i.e.N=0.The soil above the tunnel crown has a constant volumetric deformation response.The analytical solution of the subsurface settlement is simplified as

Based on the field monitoring data of cases 4-6 reported in Loganathan and Poulos(1998),the parameterBcan be obtained by back analysis.

Figs.24 and 25 compare the predicted and observed data of surface settlement troughs and subsurface maximum settlement,indicating the applicability of the analytical solution in this study.In fact,the applicability of Eq.(21) in undrained clays has been verified by Wang et al.(2020).In this section,comparisons between the predicted data obtained by this study and those obtained by Loganathan and Poulos(1998)are also presented in Figs.24 and 25.It can be seen that the analytical solution in this study provides a more accurate prediction than that by Loganathan and Poulos(1998).

6.Conclusions

By performing a series of numerical analyses on the subsurface settlement trough profiles and soil volume loss induced by tunnelling,this paper reveals three volumetric deformation modes of sandy cobble soil above the tunnel crown.Stochastic medium theory is used to derive the analytical solutions to predict the subsurface settlement for each volumetric deformation mode.The fitting formulae of key parameters in the analytical solutions are presented based on the numerical results.The conclusions are drawn,based on the presented results,as follows:

(1) The cover depth ratio(C/D),tunnel volume loss rate(ηt)and volumetric block proportion (VBP) have significant influences on the characteristics of subsurface settlement trough and the volumetric deformation modes of the soil above the tunnel crown.Dimensionless variablei/Dincreases withC/Dbut decreases with ηtat any depth.For a shallow tunnel,trough width decreases withVBPat any depth,while the variation of trough width withVBPfor a deep tunnel varies at different depth due to the soil arching effect.The normalized maximum settlementSmax/Ddecreases withC/DandVBP,but increases with ηtat any depth.With the variations ofC/D,ηtandVBP,three volumetric deformation modes can evolve with each other.

(2) A deeper and smaller tunnel (i.e.higher level of confining pressure and lower magnitude of shear strain)with lower ηt(i.e.lower magnitude of shear strain)or lowerVBP(i.e.lower content of coarse particle) tends to result in a more contractive (less dilative) volumetric deformation response of the soil.On the contrary,a more dilative(less contractive)volumetric deformation response occurs for a shallower and larger tunnel with higher ηtor higherVBP.

(3) The fitting formulae of the key parameters (i.e.BandN) in the analytical solutions are presented,providing a reliable and convenient way to determine the key parameters when the field monitoring data are unavailable.

(4) Comparisons with the results of model test and numerical simulation indicate that the proposed analytical solutions and the fitting formulae are reasonable and available in predicting the subsurface settlement induced by shield tunnelling in sandy cobble stratum.In addition,the analytical solution also provides a more accurate prediction in undrained clays than the solution by Loganathan and Poulos(1998).

It should be pointed out that this paper considered sandy cobble soil as a two-component material composed of blocks and soil matrix.Numerical analyses of tunnelling-induced ground movement characteristics were performed in the framework of continuum mechanics.In the following work,numerical analyses based on the discrete element method may be a good choice to explore the mechanical properties of sandy cobble soil and reveal the mechanism of tunnelling-induced ground movement.Meanwhile,the proposed analytical prediction only focused on the vertical ground settlement.In the future,it is worthwhile to investigate the characteristics of horizontal ground deformation induced by shield tunnelling in sandy cobble stratum and propose the corresponding prediction method.Finally,this paper only considered the case that the lining moved downward due to the self-weight.In reality,in some cases,the lining moves upward due to the grouting pressure and buoyancy forces.The difference of deformation pattern of the tunnel wall will affect the volumetric deformation response of the soil above the tunnel crown,deserving further research.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

This study was supported by the National Natural Science Foundation of China (Grant Nos.51538001 and 51978019).The authors appreciate their financial support.