Experimental–numerical insights into large-diameter shield tunneling through sand–clay and sand–pebble composite strata

Abstract

Large-diameter shield tunneling in composite strata presents significant challenges due to the heterogeneous mechanical properties, posing substantial risks to excavation safety. Existing research, however, provides limited insights into how specific layer configurations influence soil strength, failure mechanisms, and the soil arching effect. This study aims to bridge these knowledge gaps by systematically investigating the mechanical behaviors of sand–clay and sand–pebble composites, as well as their excavation-induced responses. Methodologically, we integrate laboratory triaxial tests on reconstituted composite specimens with 3-dimensional (3D) finite element simulations of the tunneling process. Key contributions of this work include: 1) the identification of a critical clay layer thickness (20–40 mm) that governs the transition from interfacial slippage to bulging failure in sand–clay composites; 2) the quantification of a progressive strength reduction (up to 14.52%) in sand–pebble composites as the sand layer shifts downward; and 3) a novel comparative analysis using a soil stress coefficient method, which reveals a more extensive destruction zone and a more fragile soil arching effect in sand–pebble strata, thereby indicating an elevated risk of collapse. These findings provide crucial insights into the performance of large-diameter shield tunnels in complex geological conditions and offer guidance for safer design and operational control.

Graphical abstract

Keywords

large-diameter shield tunneling / composite strata / laboratory triaxial test / numerical simulation / mechanical property / excavation response

1 Introduction

With rapid urbanization and the growing demand for underground space, large-diameter (generally larger than 10 m [1]) shield tunneling has become the preferred method for major tunneling projects, such as cross-river and subsea tunnels, due to its high construction efficiency and enhanced safety [2]. However, the excavation of these large tunnels often traverses composite strata [3] characterized by distinctly different geotechnical properties, such as sand–clay or sand−pebble layers. The heterogeneity across these stratified soils leads to complex shield-strata interactions, posing critical engineering challenges, including uneven ground settlement and excavation face instability [4]. Therefore, understanding the excavation-induced soil response in such conditions has become a pivotal research focus.

To address these challenges, numerous studies have been carried out to examine the geotechnical behavior of composite soils. Studies on the mechanical properties of composite soils have primarily utilized triaxial tests to investigate the influence of factors like saturation [5], particle content [6–8], and density [9] on soil behavior. Specifically, Soroush and Soltani-Jigheh [6] conducted unconfined monotonic, cyclic, and post-cyclic triaxial compression tests on clay-sand and clay-gravel composite specimens to examine the effects of variables such as sand or gravel content on the mechanical behavior of the mixtures. Huang et al. [9] examined the relationship between relative density and shear strength parameters of sandy pebbles under various disturbed conditions through indoor triaxial testing. While these studies have successfully investigated the effects of these intrinsic properties, a clear research gap exists: the critical influence of the spatial arrangement of the soil layer, specifically the stratigraphic configuration, such as layer thickness and vertical position, on the overall strength and failure mechanisms of composite soils remains underexplored.

In parallel, research on tunnel face stability and ground response has employed theoretical analysis, physical modeling, and numerical simulations. Theoretical approaches, such as limit equilibrium [10–12] and limit analysis [13,14] methods, have established fundamental models for face stability. Physical model tests have offered valuable insights into failure mechanisms under various conditions [15–19]. More recently, numerical simulations have become a primary tool for analyzing deformation and failure in composite strata. Recent studies have investigated surface deformation behavior by examining the effects of construction parameters and geological conditions [20,21]. Moreover, face stability in soft–hard composite strata has been extensively studied, with particular emphasis on the effects of parameters such as stiffness and strength ratios on failure mode [22–26]. Specifically, in coastal soft–hard layered formations, Wu et al. [20] utilized numerical simulations to explore ground deformation mechanisms, identifying optimal configurations of stratum ratio and earth pressure to mitigate subsidence and uplift. Wei et al. [26] examined the influence of relative stiffness and unconfined compressive strength ratios on excavation face failure mechanisms in upper-hard-lower-soft composite strata through numerical and theoretical analyses.

However, despite these advances, several critical issues remain unresolved. 1) The failure mechanisms of the tunnel face in heterogeneous strata, particularly the influence of the soil interface position relative to the tunnel, have not been adequately investigated; 2) a comparative analysis of the soil arching effect and collapse risk in different common composite strata (i.e., sand–clay vs. sand–pebble) has not yet been thoroughly conducted.

This study aims to address the aforementioned research gaps. The primary objectives and innovations are: 1) to experimentally quantify the influence of layer configuration, including thickness and position, on the strength, deformation, and failure characteristics of sand–clay and sand–pebble composite specimens, an area that remains underexplored in existing studies; 2) to systematically elucidate the effects of tunnel burial depth and stratigraphic interface position on ground surface deformation and the failure patterns at the tunnel face; 3) to propose and implement a novel approach based on soil stress coefficients to visualize and quantitatively compare the soil arching effect between sand–clay and sand–pebble composite strata, providing new insights into collapse risk.

This paper is structured as follows. Section 2 describes a series of systematic triaxial tests conducted on composite specimens. Section 3 provides a comprehensive description of the three-dimensional numerical simulations conducted to investigate tunnel face stability and excavation-induced ground responses. Finally, the main conclusions are summarized to offer practical guidance for the safer design and operation of tunnels in heterogeneous geological conditions.

2 Experimental investigation on mechanical behavior of composite strata

This section investigated the physical properties of composite strata specimens obtained from a large-diameter tunneling project in Beijing. The subsurface conditions are characterized by Quaternary alluvial deposits typical of the North China Plain, primarily consisting of silt, silty clay, fine sand, medium sand, and pebble layers. The geological longitudinal section of the strata and typical geological interface traversed by the shield tunnel is presented in Fig. 1. It is a schematic geological profile based on the integration of engineering borehole logs and geotechnical reports. The profile aims to reflect the typical stratigraphic conditions and composite layer interfaces encountered along the tunneling path, and serves as the basis for the numerical simulation model. The tunnel face primarily traversed sand–clay and sand−pebble composite strata. Accordingly, this study focused on a systematic analysis of two representative composite strata: sand–clay and sand−pebble strata.

2.1 Sand–clay composite strata

2.1.1 Materials and specimen preparation

The test soils were obtained from the aforementioned large-diameter tunneling project in Beijing. The mechanical behavior of clay is primarily governed by its plasticity characteristics. Therefore, the Atterberg limits are provided as the key classification parameters, in accordance with standard geotechnical practice. The physical properties of the clay and sand are summarized in Table 1, and the particle size distribution curve of the sand is presented in Fig. 2.

The configuration of the composite specimens, specifically the thickness and position of the clay layer, was designed to reflect the stratigraphic variability observed in the geological data from the Beijing large-diameter tunneling project (Fig. 1(b)). As depicted schematically in Fig. 1, the tunnel frequently traversed interbedded sand–clay strata where the relative thickness of the clay layers varied significantly. To systematically investigate this variability in the laboratory, the specimen design was based on the ratio of the clay layer thickness to the total specimen height (t/H). Given the standard specimen height of 80 mm, two representative scenarios were established. A 20 mm clay layer was selected, representing a thickness ratio of 25%. This configuration models a common field condition of a distinct, relatively thin clay seam within a predominantly sandy stratum. A 40 mm clay layer was chosen, corresponding to a thickness ratio of 50%. This scenario represents a more substantial clay layer where sand and clay are co-dominant components of the composite ground. These two ratios provide a systematic basis for parametrically evaluating the influence of the clay layer’s relative thickness on the overall strength and deformation behavior. The specimen size is φ39.1 mm × 80 mm, with the specimen grouping schematic illustrated in Fig. 3. By varying the position and thickness of the clay layer, the effects on the strength and deformation behavior of the composite specimens were systematically investigated.

A layered compaction method was employed using a saturator: for homogeneous specimens, sandy soil was compacted into 5 layers with equal mass distribution per stratum, while clay specimens required 7 layers of identical mass configuration. Composite specimens were fabricated using different layering strategies: specimens with a 20 mm clay interlayer consisted of 2 clay layers alternating with 3–4 sandy soil layers, while those containing a 40 mm clay interlayer exhibited a configuration of 3 clay layers and 3 sandy soil layers. Following compaction of each layer to the target height, the interfacial surface was scarified to ensure bonding integrity before proceeding with subsequent layer placement. This iterative process produced fully compacted specimens. All specimens were vacuum-saturated at low pressure for 15 h and stored in distilled water to prevent moisture evaporation prior to testing.

2.1.2 Triaxial testing systems and experimental design

The experimental apparatus employed was a standard stress path triaxial testing system, comprising a Bishop and Wesley triaxial pressure chamber, a pressure controller (GDS instruments, UK), and a data acquisition system. The technical specifications of the system were as follows: specimen dimensions of φ39.1 mm × 80 mm; maximum axial load capacity of 7 kN; maximum confining pressure of 2 MPa; and a maximum axial displacement of 25 mm. Each fabricated specimen was enclosed in a rubber membrane and placed on the base of the pressure chamber. Back pressure was applied in incremental steps starting from 50 kPa, while confining and axial pressures were simultaneously applied, each maintained at 20 kPa above the corresponding back pressure and increased in 50 kPa increments. Once the pore pressure stabilized, the B-value was checked to ensure it exceeded 0.98 for achieving full saturation. Notably, the total confining pressures for sand and sand–clay specimens were 550, 650, and 750 kPa, whereas the confining pressures for clay specimens were 150, 250, and 350 kPa, respectively. The reason for the different confining pressure settings was that a back pressure of 450 kPa was required for the sand and sand–clay composite specimens, whereas 50 kPa was sufficient for the clay specimens due to their different pore structures and permeability. After consolidation was complete, all specimens were sheared under undrained conditions at a constant strain rate of 0.07%/min until the axial strain reached 15%–20%.

2.1.3 Results analysis

Figure 4 demonstrates distinct shear-induced deformation patterns across specimen types. The pure sand specimen developed pronounced barrel-shaped bulging at mid-height due to differential lateral constraints. In contrast, pure clay exhibited faint bulging despite similar deformation characteristics, a consequence of its inherent high compressibility. For sand–clay composites, bulging deformation localized predominantly within the clay stratum, highlighting the critical role of material heterogeneity in governing strain distribution.

Figure 5 presents the deviatoric stress-strain curves for each specimen under varying confining pressures. The pure sand specimen exhibited strain softening behavior, reaching peak deviatoric stress at approximately 8% axial strain. In contrast, the pure clay specimen showed strain hardening characteristics, attributed to its high compressibility. Sand–clay composite specimens with a 20 mm clay layer exhibited a three-stage response: 1) rapid initial deviatoric stress increase (0%–2% strain); 2) a moderated strengthening phase (2%–13% strain); and 3) a stabilized co-deformation phase. Conversely, the specimen with a 40 mm clay layer displayed only a two-stage response: 1) rapid initial deviatoric stress increase (0%–1% strain); and 2) a moderated strengthening phase beyond 2% strain. Notably, in the 40 mm clay specimen, incomplete interfacial consolidation persisted even at 18% strain, leading to early interfacial slippage and a lower peak strength compared to the thinner-clay counterpart, thereby preventing the development of a third deformation phase. Based on the analysis of strength curves and failure modes of specimens with 20 and 40 mm clay layers, it was inferred that a critical clay thickness exists, governing the deformation and failure mechanisms of sand–clay composites. This critical thickness was within the range of 20–40 mm.

Figure 6 illustrates the variation in peak deviatoric stress across the tested specimens. As the deviatoric stress-strain curves of the composite specimens did not exhibit a distinct peak, the deviatoric stress at 15% axial strain was adopted as the peak value for comparison. The use of 15% axial strain is a widely accepted convention in large-strain geotechnical testing for strain-hardening materials [27–29]. The results indicated that the peak deviatoric stress of sand–clay composite specimens decreased with both the downward displacement of the clay layer and the increase in its thickness. This trend suggested that, within sand–clay composite strata, a lower clay layer position and greater clay thickness were associated with reduced overall specimen strength.

The tests in this study used the p–q' (deviatoric stress-average effective stress) curve of the soil to determine the shear strength indices c' and φ' using Mohr−Coulomb theory. The obtained values of the effective angle of internal friction and effective cohesion for each specimen are shown in Table 2.

2.2 sand−pebble composite strata

2.2.1 Materials and specimen preparation

The maximum size of the particles in the triaxial test needs to be less than 20% of the specimen diameter, so a large triaxial tester is required for testing coarse-grained soils. According to the actual size distribution characteristics of sand and pebble strata, specimens were divided into four groups, with the size of φ300 mm × 600 mm, and the schematic diagram is shown in Fig. 7.

Soil samples were collected from the excavation face of the Beijing Large-Diameter Shield Tunnel Project. In situ samples were obtained by utilizing a pressurized chamber opening during the maintenance of the shield machine. The collected samples were sealed on site, air-dried, and subsequently sieved using six standard sieves with aperture sizes of 60, 40, 20, 10, 5, and 2 mm. The resulting particle size distribution of the original pebble stratum is presented in Table 3. Sieve analysis indicated that particles larger than 60 mm accounted for more than 5% of the total mass. To prepare large-scale triaxial test specimens, the oversized particles were replaced using the equal mass substitution method, in which particles within the 5–60 mm range were proportionally substituted by mass. The particle size distribution curves of both the original and remolded test soils are shown in Fig. 8.

The target unit weight (γ) was 22 kN/m³ for the pebble specimens and 19 kN/m³ for the sand specimens. Pebble specimens were prepared based on the remolded gradation curve, while sand specimens were composed of 75% sand, 10% gravel (20–40 mm), 10% gravel (10–20 mm), and 5% gravel (5–10 mm). To ensure uniform density across specimens, each sample was compacted in four layers, with a fill height of 150 mm per layer. Each layer was tamped using a heavy hammer from a loose state to the target height, and the height deviation of the final layer was strictly controlled within 2 mm.

2.2.2 Triaxial testing systems and experimental design

The test equipment was the TAJ-2000 large-scale dynamic and static triaxial apparatus. The equipment was capable of automatically collecting vertical load and displacement, and the specimen loading was controlled by displacement. The technical indexes and parameters of the system were as follows: the specimen size was φ300 mm × 600 mm; the maximum axial test force was static 2000 kN; the maximum confining pressure in the confinement system was 10 MPa; the control system adopted the stress and strain control, and the vertical and circumferential loading synchronized or controlled by phase.

Prior to testing, specimens were vacuum-saturated and subsequently saturated with water. Following saturation, isotropic consolidation under drained conditions was performed. Upon completion of consolidation, compression tests were conducted under strain-controlled conditions at a constant shear rate of 1 mm/min. Failure was defined as the point at which axial displacement reached 15% of the specimen height (i.e., 90 mm), at which time data and stress−strain curves were recorded. For multi-stage loading tests, each loading stage involved an axial displacement corresponding to 5% of the specimen height. The specimen was considered to have failed when the cumulative displacement reached 90 mm (15% of the height of the specimen).

2.2.3 Results analysis

In triaxial testing, the ends of the specimen were constrained, resulting in less lateral deformation compared to the middle section. Consequently, the specimen exhibited bulging at the center, where the constraint was weakest. The maximum diameter observed in the bulging region was defined as the bulging diameter. Table 4 summarizes the magnitude and characteristics of deformation observed in each specimen. The homogeneous pebble specimen exhibited a minimal bulge deformation of 10%. In contrast, sand−pebble composite specimens demonstrated significantly greater deformation, ranging from 18.7% to 20%. Notably, upper sand specimens subjected to variable confining pressures (200–600 kPa) exhibited a confining pressure-dependent reduction in bulge deformation, decreasing from 20% to 16.7%. A positive correlation was observed between deformation magnitude and the depth of the sand layer: specimens with basal sand exhibited the highest deformation (19.3%), followed by those with mid-layer sand (18.7%), and the lowest strain occurred in specimens with upper-layer sand configurations (16.7%). Bulging diameter measurements (ranging from 350 to 358 mm) corroborated this trend, with composite specimens showing 8–9% greater deformation compared to homogeneous pebble specimens.

Upon completion of the test, internal damage characteristics were assessed by removing the rubber membrane from the specimens. As shown in Fig. 9, the upper pebble-lower sand specimen exhibited a reduction in sand layer thickness from the initial 150 to 130 mm, corresponding to a 13.3% deformation. This value was lower than the total axial deformation of the specimen (15%), indicating that the pebble layer experienced greater compression than the sand layer. Figure 10 illustrates the deviatoric stress−strain responses of the tested specimens. All specimens exhibited strain-hardening behavior without a distinct peak; accordingly, the deviatoric stress at 15% axial strain was adopted as the representative peak value for subsequent analysis. The peak deviatoric stress values under varying confining pressures are summarized in Table 5.

At a confining pressure of 200 kPa, the homogeneous pebble specimen exhibited the highest peak deviatoric stress. For the composite specimens, the peak deviatoric stress decreased progressively with the downward migration of the sand layer, ranging from 87.59% (527.12 kPa) to 70.07% (421.68 kPa) of the pebble specimen’s strength. This trend indicated that the overall strength of the sand−pebble composite stratum diminished as the sand layer was positioned deeper within the specimen.

Figure 11 presents the strength envelopes for all specimens, indicating zero cohesion across the board. The homogeneous pebble specimen exhibited the highest internal friction angle at 38.8°. In contrast, the composite sand−pebble specimens demonstrated a progressive decline in internal friction angle, decreasing from 38.0° to 34.5° as the sand layer position shifted downward. This trend suggested that lower sand layer placement adversely affected the shear strength of the composite stratum.

Table 6 summarizes the mechanical parameters and shear strength of each specimen, using the whole pebble specimen as a reference. The composite specimens exhibited strength reductions corresponding to the position of the sand layer: 2.83% for the upper layer, 6.28% for the middle layer, and 14.52% for the lower layer. These results highlighted a clear inverse relationship between sand layer depth and overall shear strength in sand−pebble composite specimens.

3 Numerical simulation of shield tunnel excavation

The numerical simulations presented in this section are designed to explore the large-scale excavation-induced responses of composite strata. To ensure the reliability of the simulation results, the constitutive models for the soils were calibrated based on fundamental mechanical parameters directly obtained from the systematic triaxial laboratory tests described in Section 2.

3.1 Simulation scheme

Considering the geological and engineering characteristics of the large-diameter tunneling project in Beijing, the shield tunnel excavation process was simulated based on the structural composition of the typical strata encountered. The simulation investigated the mechanical behavior and deformation patterns of soils in composite strata under varying stratum configurations and tunnel burial depths, with a focus on analyzing the failure mechanisms at the excavation face. The detailed numerical modeling schemes for soil excavation are summarized in Table 7.

3.2 Numerical models

3.2.1 Constitutive model

The shear strength parameters of sand, clay, as well as pebble layers were obtained from the consolidated undrained triaxial tests on sand–clay composite specimens and the drained triaxial tests on sand−pebble composite specimens in Section 2, respectively. The Mohr-Coulomb constitutive model was used for the calculations. The equations for the yield surface of this intrinsic model are as follows:

(1)Rq-p\tan \phi -c=0,

where p is the mean stress, q is the shear stress, \phi is the slope of the Mohr−Coulomb yield surface on the q−p plane, and c is the cohesive force. R controls the shape of the yield surface in the \text{π}-plane and it is calculated by:

(2)R=\frac{1}{\sqrt{3}\cos \phi }\sin \left(\mathit{\Theta }+\frac{\text{π}}{3}\right)+\frac{1}{3}\cos (\mathit{\Theta }+\frac{\text{π}}{3})\tan \phi ,

where \mathit{\Theta } is calculated by:

(3)\cos \left(3\mathit{\Theta }\right)=\frac{I_3^3}{q^3},

where I₃ is the third bias stress invariant.

The yield surface of the Mohr−Coulomb model is a straight line in the p−q plane. However, in the π-plane or meridional plane, the conventional associated flow rule becomes numerically unstable when applied to yield surfaces with vertex singularities. The simulation algorithm adopted a non-associated flow rule formulation, thereby ensuring the unique determination of plastic strain increments.

3.2.2 Computational models and parameters

The three-dimensional finite element model of the sand–clay composite stratum, when the soil layer interface is located in the center of the tunnel, is illustrated in Fig. 12. The computational domain extended 81 m along the tunnel advancement direction (Z-axis) and 110 m perpendicular to the excavation face (X-axis). Vertical dimensions (Y-axis) varied with different tunnel burial depth ratios (C/D). The model heights were taken as 54 m (C/D = 0.5), 60 m (C/D = 1), and 72 m (C/D = 2).

All elements in the model were 3 dimensional (3D) 8-node hexahedral elements, with solid elements for the soil and uncoordinated elements for the shield, lining, and grouting. A non-uniform meshing scheme was implemented in current simulations. The mesh was significantly refined in the region adjacent to the tunnel face.

The soil strata were modeled as homogeneous, isotropic, elastic–perfectly plastic materials based on the Mohr−Coulomb failure criterion. The parameters used for each soil layer in the numerical model (Table 8) were directly derived from the corresponding homogeneous triaxial tests (pure clay, pure sand and pure pebble tests). The shield machine, segmental lining, and grouting material were simulated using a linear elastic constitutive model. The simulation was assumed to be quasi-static, approximating continuous excavation by a series of discrete steps and neglecting dynamic and inertial effects.

A geostatic step was performed prior to excavation to establish the initial stress field from self-weight and ensure the model reached a stable equilibrium. This was verified by monitoring the maximum unbalanced force ratio within the model until it became negligible. This same convergence criterion was applied after each incremental excavation step to ensure that stress redistribution and deformation were calculated from a quasi-static equilibrium state. The key parameters of the present numerical simulations are shown in Table 8.

The numerical simulation of shield tunneling advancement adopted a staged construction approach, where the excavation step length (ΔL) was defined as the width of a single precast segment ring (1.5 m). To mitigate boundary interference effects, the tunnel face was positioned at a distance of 3D (D = tunnel diameter) from the front boundary of the finite element domain. In this context, the “time step” corresponded to increments within each analysis step. The ABAQUS/Standard solver automatically adjusted these increments to ensure that a converged solution was found at the end of each step, effectively managing the simulation’s stability and accuracy. The excavation process was simulated through sequential element operations: 1) deactivation of excavated soil elements corresponding to the shield position; 2) simultaneous activation of lining elements and pressurized grout elements behind the shield tail; 3) application of equivalent face support pressure.

3.2.3 Determination of limit support pressure

The stability of the tunnel excavation face is critically dependent on the applied support pressure. The objective of this analysis was to determine the lower and upper bounds of the limit support pressure range, known as the active and passive limit support pressures.

A stress-controlled approach was implemented in the simulations for achieving face instability. First, an initial support pressure (P₀) was calculated, which was taken as the average of the static and active earth pressures at the excavation face:

(4)P_0=\frac{\left(K_0{\sigma }_{\nu }^{\prime }+K_a{\sigma }_{\nu }^{\prime }\right)}{2},

where K₀ is the static earth pressure coefficient, K_a is the active earth pressure coefficient. As shown in Table 9, this initial pressure is primarily dependent on the tunnel burial depth and is not significantly affected by the composite strata interface position for the same depth ratio. Therefore, under conditions with the same burial depth ratio, the initial support force P₀ was assumed to be the same for the tunnel excavation face in sand–clay composite strata. The support force behavior in the sand−pebble composite stratum followed a similar pattern to that in the sand–clay composite stratum. The support pressure ratio f (defined as the ratio of current support pressure to P₀) was incrementally decreased or increased until tunnel face instability occurred. The instability of the tunnel face in current simulations was defined by two criteria: 1) a sudden increase of maximal horizontal displacement or 2) the maximal horizontal displacement exceeded 100 mm, which are widely adopted in tunnel excavation simulations [30–33].

3.2.4 Model validation

To verify the accuracy and reliability of the established three-dimensional finite element model, a validation study was conducted by comparing the simulation results with a published case study. The work by Cui et al. [34] was selected as a benchmark, as it provides both field monitoring data and numerical simulation results for a large-diameter shield tunnel excavated in a comparable geological setting in Beijing, specifically sandy pebble soil. The transverse distance from the tunnel axis (x) was normalized by the tunnel diameter (D) as x/D, and the vertical surface settlement (S) was normalized by the maximum settlement (S_max) at the axis as S/S_max.

Figure 13 presents the comparison between the normalized transverse settlement troughs obtained from our simulations for sand−pebble strata (at burial depth ratios C/D = 1 and C/D = 2) and the field data and simulation results reported by Cui et al. [34]. The results show very good agreement in the overall trend and shape of the settlement troughs. The width and curvature of the ground deformation profile obtained by our model are highly consistent with both the in situ measurements and the benchmark simulation. This strong correlation confirms that the numerical model developed in this study can accurately capture the fundamental ground response mechanisms to shield tunneling, thus validating its use for the subsequent parametric analyses.

3.3 Shield tunnel excavation soil response analysis

3.3.1 Limit support pressure

Figure 14 illustrates the relationship between f and maximum horizontal displacement for clay-sand composites under varying cover-to-diameter ratios (C/D). Five distinct behavioral phases were observed: 1) Catastrophic failure (f = 0.1–0.4): Minor reductions in f triggered rapid inward displacement (> 100 mm); 2) progressive failure (f = 0.4–0.9): Gradual inward displacement increased with decreasing f; 3) stable phase (f = 0.9–1.5): Linear displacement response, indicating face equilibrium; 4) passive failure onset (f = 1.5–3.0): Slow outward displacement growth, escalating with f; 5) passive collapse (f > 3.0): Sudden displacement surged in the tunneling direction, exceeding 100 mm. Critical support ratios were derived: 1) active failure: f = 0.2 (C/D = 0.5), 0.3 (C/D = 1), 0.4 (C/D = 2); 2) passive failure: f = 3.0, consistent across all C/D ratios.

3.3.2 Factors influencing the surface deformation pattern

The deformation characteristics of the ground surface were analyzed under different stratum combinations and burial depth ratios when the shield tunnel was pushed forward to 52.2 m. The settlement and uplift curves of the ground surface were investigated separately for active damage (f = 0.2, 0.3, 0.4 for C/D = 0.5, 1, 2, respectively) and passive damage (f = 3.0 for all C/D) to the excavation face.

1) Tunnel burial depth ratio

Figure 15 illustrates the surface deformation curves at various burial depth ratios for sand–clay and sand−pebble composite strata, where the soil layer interface crosses the tunnel axis. For sand–clay composite strata, surface deformation under active failure conditions directly correlated with the burial depth ratio in both upper clay-lower sand and upper sand−lower clay configurations. The analysis of surface heave reveals that shallow tunnels are particularly susceptible to passive failure, resulting in pronounced ground uplift near the tunnel axis. This underscores the importance of accurate real-time control of face support pressure to prevent over-pressurization, especially when tunneling beneath sensitive infrastructure. As the burial depth ratio increased, the maximum settlement and the width of the settlement trough both expanded (Fig. 16). Under passive failure conditions, surface deformation exhibited spatial dependence: within 1.5D of the tunnel axis, heave decreased with increasing burial depth ratio, while beyond 1.5D, the trend reversed, with heave increasing at larger burial depth ratios. In sand−pebble composite strata, the settlement trough width also increased with burial depth ratio, similar to the sand–clay strata, indicating a broader ground disturbance at greater depths. Heave behavior, however, differed at critical lateral positions: within 1.7D from the axis, heave decreased as the burial depth ratio increased, whereas beyond this threshold, heave amplified with larger burial depth ratios. Notably, in lower burial depth ratio conditions, close monitoring of near-axis uplift is critical due to the pronounced heave effects. The finding that trough width increases with the burial depth ratio (Fig. 16(a)) implies that for deeper tunnels, the surface settlement monitoring must cover a wider lateral area to capture the full extent of ground movement.

2) Stratigraphic interface

Figure 17 illustrates the surface deformation curves of the tunnel cross-section (z = 52.2 m) and longitudinal section (x = 0) under different stratigraphic interface positions in sand–clay composite strata, with C/D = 2. Under active failure conditions, the surface settlement curves for the cross-section closely followed the Peck curves, with the maximum settlement occurring directly above the tunnel axis. In contrast, the surface uplift curve (under passive failure conditions) deviated from this pattern, with the peak uplift occurring approximately 15–20 m (1.3D–1.7D) from the tunnel axis. The uplift directly above the tunnel was comparatively smaller due to the hardening effect of the grouting unit—at the time of excavation, the grout was not yet fully cured and possessed low strength, allowing the overlying soil to move downward slightly, which limited surface uplift at the centerline. In the longitudinal section, the maximum settlement appeared directly above the excavation face, stabilizing beyond a distance of 22 m (1.8D) from the face. Conversely, the maximum heave was observed approximately 15 m (1.3D) ahead of the excavation face, and remained essentially constant once the distance exceeded 20 m (1.7D).

Figure 18 presents the surface deformation curves for the tunnel cross-section (z = 52.2 m) and longitudinal section (x = 0) under varying stratigraphic interface positions in sand−pebble composite strata, with C/D = 2. Tunnel excavation in sand−pebble composite strata produced a broader zone of surface deformation compared to sand–clay strata. In the cross-section, settlement effects became negligible beyond 40 m (approximately 3.4D) from the tunnel axis. The maximum surface heave occurred at a lateral distance of around 25 m (2.1D), with heave magnitudes stabilizing beyond 45 m (3.8D). Longitudinal profiles indicated that both settlement and heave magnitudes reached a stable state once the distance from the excavation face exceeded 25 m (2.1D) during shield tunneling.

Maximum surface settlement occurred when soil interfaces intersected the tunnel centerline in upper clay-lower sand and upper sand−lower pebble strata. As demonstrated in Fig. 19, settlement magnitudes reduced with increasing interface proximity to the ground surface (i.e., thinning overlying clay/sand layers). Additionally, peak heave was observed under centerline-aligned interfaces, while upward-shifted interfaces demonstrated negligible positional sensitivity in heave magnitude and disturbance range. Notably, upper sand−lower clay and upper pebble-lower sand strata exhibited inverse deformation trends governed by interface elevation. The strong correlation between the interface position and peak ground settlement (Fig. 19) highlights its potential as a predictive indicator. Specifically, in upper sand−lower clay stratum, settlement tends to increase significantly as the tunnel approaches the clay stratum.

3) Stratigraphic composition

Figure 20 represents the surface deformation curve of the cross section when the stratigraphic interface passes through the tunnel center under different stratigraphic combinations (C/D = 2). In sand–clay composite strata, upper sand−lower clay configurations exhibited greater heave magnitudes compared to upper clay-lower sand counterparts, while the latter generated larger settlements with broader settlement troughs. For sand−pebble composites, stratigraphic composition exerted limited influence, with slight differences in settlement magnitudes. Within ±1.3D, upper pebble-lower sand strata showed higher heave values with minimal variation (maximum differential: 0.75 mm), while beyond this range, the trend reversed. These findings indicate that the stratigraphic compositions in sand–clay composite strata are more influential on surface deformation than in the sand−pebble case.

3.3.3 Factors influencing the failure patterns

For the sand–clay and sand−pebble composite strata, the damage modes of the tunnel excavation surface under different burial depth ratios, stratigraphic combinations, and soil stratigraphic interfaces were investigated by analyzing the displacement changes of the soil layer under the limit state of the excavation surface.

1) Tunnel burial depth ratio

Taking the upper clay-lower sand composite stratum as a benchmark case, the failure mechanisms of the tunnel excavation face were analyzed under varying burial depth ratios, with the soil interface aligned at the tunnel centerline. Upon destabilization of the excavation face, the longitudinal displacement fields of the surrounding soil at three different burial depth ratios are illustrated in Fig. 21. In the figure, arrows denote the displacement vectors: their orientation reflects displacement direction, while their length indicates relative magnitude. To further evaluate the failure mode of the excavation face, the potential slip surfaces were delineated based on the distribution and intensity of the displacement vectors observed after tunnel advancement.

Active destabilization of tunnel excavation faces induced wedge-shaped failure zones in the soil mass ahead of the face, characterized by differentiated slip surface geometries: the lower segment presents gentler inclination angles, whereas the upper region displays broader spatial propagation. Comparative analysis indicated that the burial depth ratio (C/D) significantly influenced these failure mechanisms. Specifically, configurations with C/D = 1 and 2 exhibited more extended and gradually inclined slip surfaces compared to the steeper, more localized failure observed at C/D = 0.5, accompanied by increased soil disturbance. Under passive failure conditions, the slip surface was constrained by the sand–clay stratigraphic interface, with the lower portion exhibiting gentle slopes and the upper portion forming steeper gradients. As the burial depth ratio increased, the lower segment of the slip surface became increasingly shallow, expanding the disturbance zone ahead of the excavation face. Concurrently, the upper portion steepened, reducing the excavation-induced influence on the overlying soil. Notably, for C/D = 2, the failure region did not propagate to the ground surface, indicating effective confinement by the overburden.

2) Stratigraphic interface

The damage modes of the excavation surface under different soil strata interface conditions at C/D = 2 were analyzed as an example of the upper clay-lower sand stratum. In this way, the influence of the stratigraphic interface on the destabilization mode of the excavation surface was investigated.

Figure 22 presents the displacement vectors of the soil under various stratigraphic interface conditions during active failure of the excavation face. Tunnel excavation caused disturbance in the surrounding soil, inducing inward movement toward the tunnel axis. This resulted in the formation of a wedge-shaped failure zone ahead of the excavation face, extending upward in a basin-like pattern toward the ground surface. When the interface intersected the excavation face, its location had minimal influence on the curvature of the slip surface, which remained continuous across the excavation front. In contrast, when the interface lay above the tunnel crown, a distinct change in the curvature of the slip surface was observed at the interface, resulting in discontinuity along the slip surface at that location.

The numerically observed behavior is consistent with the physical model tests on layered soils reported by Ma et al. [16]. The discontinuity of the slip surface arises directly from the pronounced mechanical contrast between the clay and sand layers. As the slip surface crosses the interface, it refracts owing to the abrupt change in shear strength and stiffness. Within the weaker clay layer, a larger plastic zone must be mobilized to provide sufficient resistance, whereas upon entering the stiffer sand layer, the failure path steepens due to the higher shear strength and stiffness of the material.

Figure 23 illustrates the displacement vectors of the soil under various stratigraphic interface conditions during passive failure at the excavation face. Unlike the basin-shaped surface deformation observed in active failure, passive destabilization was characterized by a localized slip surface that did not extend to the ground surface; the failure zone remained confined to the vicinity ahead of the excavation face. When the stratigraphic interface intersected the tunnel, the failure mechanism displayed a combination of wedge-shaped damage in the lower section and “silo-type” failure in the upper section. This provides a clear illustration of well-established theoretical models. The observed failure mode is consistent with the classic limit equilibrium mechanisms originally proposed by Horn [10] and subsequently refined by Zhang and Hu [12]. Moreover, the silo-shaped failure observed in the upper, more confined stratum reflects vertical shearing and soil arching effects, which constitute central aspects of the face stability analyses by Anagnostou and Kovari [11]. If the interface lay above the tunnel crown, failure became localized and spherical, with minimal disturbance to the surrounding soil. In this case, damage was concentrated near the excavation face. However, when the interface was positioned closer to the top of the tunnel, the soil above the shield structure experienced greater disturbance. The extent of the failure zone strongly correlated with the proximity of the stratigraphic interface. Once the vertical distance between the interface and the tunnel crown exceeded approximately 1.2D, the influence of interface position on failure patterns became negligible, suggesting a critical threshold beyond which stratigraphic effects on passive destabilization diminish.

3) Stratigraphic composition

The influence of different soil combinations on the instability modes of tunnel excavation faces was systematically analyzed under conditions of a burial depth ratio C/D = 1, with the stratigraphic interface intersecting the tunnel centerline.

Under active failure conditions (Fig. 24), the clay−sand composite stratum demonstrated a prominent wedge-shaped failure zone ahead of the excavation face, coupled with a basin-like disturbance zone above it. This configuration featured a narrower base and a broader upper extent, indicating the largest disturbance range among all tested stratigraphic combinations. In contrast, other soil configurations consistently exhibited wedge-shaped failure zones beneath the interface and silo-shaped failure zones above. Among these, sand–clay strata displayed greater susceptibility to instability compared to sand−pebble strata, highlighting the significant role of stratigraphic composition in modulating excavation-induced failure mechanisms.

The observed difference in stability can be attributed to the contrasting geomechanical properties of clay and pebble strata. Clay, with its lower strength and greater deformability, provides limited resistance to the formation of a failure wedge. Its plastic behavior induces a wider and deeper basin-shaped subsidence zone. By contrast, pebbles develop a dense, interlocking granular skeleton with a high friction angle. This interlocking enhances internal stability and promotes a pronounced soil arching effect, which confines the failure zone within a compact, silo-shaped region. The arching action in the pebble layer redistributes stresses more effectively than the plastic deformation in the clay layer, leading to a smaller and more localized failure region.

Under passive failure conditions (Fig. 25), the clay−sand composite stratum exhibited a gently inclined wedge-shaped slip surface below the interface and a basin-shaped disturbance zone above. The presence of clay beneath the interface significantly amplified deformation, contributing to broader and deeper failure zones. In contrast, sand−pebble combinations induced more localized instability, characterized by a wedge-shaped failure zone beneath and a silo-shaped disturbance zone above the interface, with limited sensitivity to the vertical arrangement of layers. Importantly, the spatial distribution of the clay layer strongly influenced the geometry of the slip surface and the extent of ground disturbance in clay-sand strata. Conversely, sand−pebble strata demonstrated negligible dependence on layer configuration.

The pronounced influence of the clay layer’s position, particularly when located beneath the interface, can be attributed to its role as a mechanically weak bearing stratum. Under passive loading, the tunnel face transmits stresses to the surrounding soil mass. When these stresses are carried by a deformable, low-strength clay layer, the stability of the supporting stratum is compromised, leading to extensive plastic flow and a widely distributed failure zone. In contrast, when sand or pebbles form the lower layer, they act as a stiff, high-strength bearing stratum that resists deformation, thereby restricting failure to a more localized region. The negligible sensitivity to layer configuration in sand–pebble strata further underscores the inherent stability of interlocking granular materials. In such systems, failure is primarily governed by internal friction and particle rearrangement, which are comparatively insensitive to stiffness contrasts between sand and pebbles—unlike the pronounced disparity observed between sand and clay layers.

These findings highlight that clay-dominated interfaces require prioritized reinforcement due to their propensity for large-scale instability, whereas sand−pebble strata offer relative stability with localized failure risks, informing optimized support strategies for heterogeneous ground conditions.

3.3.4 Soil arching effect

The soil arching effect is a critical stress redistribution mechanism triggered by shield tunneling, arising from the reorientation of principal stress directions and the evolution of interparticle force chains within heterogeneous strata [35]. This phenomenon facilitates the transfer of overburden loads to surrounding stable regions, significantly influencing the stability of the excavation face and the resulting surface deformation patterns. Chen et al. [36] proposed a method for estimating the extent of soil arching above tunnel excavations by analyzing the redistribution of vertical stress during tunneling. They suggested that the inflection points in vertical stress-depth curves at different locations above the tunnel approximately delineate the boundary of the soil arch. In analyzing the soil arching effect during tunnel excavation, the arch was typically divided into a loose zone and an arch zone based on vertical stress distribution [37–39]. The influence of the soil arching effect was then evaluated by examining the height of the soil arch [40,41].

In this study, the soil stress evolution was quantified through normalized vertical and horizontal stress coefficients: α_v = σ_v / σ_v0, α_h = σ_h / σ_h0, where σ_v0 and σ_h0 denote initial geostatic stresses. Compared with traditional approaches that rely solely on vertical stress redistribution, the presented method used the normalized stress coefficients (α_v, α_h) in the vertical and horizontal directions to distinguish the three zones of the arch (destruction zone, arch foot and arch crown) and allowed for a more refined understanding of the load transfer mechanism within the soil arch. Numerical simulations focused on transverse (A-A, z = 52.2 m) and longitudinal (B-B, x = 0 m) monitoring sections (Fig. 26), analyzing stress zoning under active failure conditions (C/D = 2). Figures 27 and 28 demonstrate the soil stress coefficients at different locations in the two profiles in the upper clay-lower sand and upper sand−lower pebble composite stratum.

Three distinct arching zones were identified across the transverse (A-A) and longitudinal (B-B) sections surrounding the tunnel excavation face. 1) Destruction Zone (α_v < 1, α_h < 1): Characterized by reductions in both vertical and horizontal stress due to stress release and soil deformation toward the tunnel face. 2) Arch Foot (α_v > 1, α_h < 1): Exhibited increased vertical stress due to load transfer from upper soil layers and reduced horizontal stress due to loosening effects. 3) Arch Crown (α_v < 1, α_h > 1): Showed decreased vertical stress from downward soil movement and increased horizontal stress due to lateral convergence and interparticle wedging. During active failure, vertical and horizontal stresses decreased in the destruction zone as soil deformed inward, forming a distinct disturbed region above the excavation face. In the arch crown, vertical stress decreased as the upper soil loosened and moved downward, while lateral compression from surrounding soil raised the horizontal stress. At the arch foot, vertical stress increased as load was transferred downward, while horizontal stress decreased due to peripheral loosening. From the stress coefficient distribution, the soil arch boundary was defined. In the sand–clay composite stratum, the destruction zone extended up to 4.00 m (0.34D) horizontally and 2.80 m (0.23D) vertically, with the full soil arch reaching a height of 18.78 m (1.57D). In the sand−pebble composite stratum, the destruction zone expanded further, to 4.71 m (0.39D) horizontally and 3.03 m (0.25D) vertically, with a slightly greater arch height of 19.00 m (1.59D). In the longitudinal section (B-B), the destruction zone ahead of the excavation face was significantly more extensive in sand−pebble strata, reaching 6.13 m (0.51D), compared to 4.51 m (0.38D) in the sand–clay strata. This disparity was attributed to the greater disturbance sensitivity and deformation tendency of the coarse-grained sand−pebble medium. While the soil arch height was similar between the two composite strata, the arch effect was more fragile in sand−pebble formations. The soil arching effect in these strata was observed to be more fragile, resulting in a larger destruction zone and thus a higher susceptibility to face instability and collapse.

4 Conclusions

This study systematically investigated the mechanical behavior and excavation responses of large-diameter shield tunnels in sand–clay and sand−pebble composite strata, yielding critical insights for tunnel design and construction.

The laboratory triaxial tests revealed a strong dependence of the mechanical behavior of soil on its stratigraphic configuration. For sand–clay composites, a critical clay layer thickness (20–40 mm) was identified that governs the transition from interfacial slippage to bulging failure. For sand−pebble composites, the overall shear strength progressively decreased (up to 14.52%) as the sand layer’s position shifted downward, highlighting the critical influence of the weaker layer’s location.

The numerical simulations clarified the influence of geological conditions on excavation-induced responses. The settlement trough width consistently expanded with increasing burial depth. For strata similar to upper clay-lower sand configuration, shallower interfaces reduce maximum settlement, while surface heave is weakly sensitive to interface position. However, the failure patterns were highly dependent on the stratigraphic interface position, which dictated the shape and extent of the failure zone, including transitions between wedge-shaped, silo-shaped, and localized spherical failures.

A novel analysis using soil stress coefficients revealed that sand−pebble strata, despite their higher strength, form a more fragile soil arching effect and a more extensive destruction zone compared to sand–clay strata. This finding indicates a higher susceptibility to collapse in sand−pebble formations. These results provide a refined risk assessment framework, emphasizing that engineers must consider specific geotechnical combinations and their spatial configurations to ensure safer tunnel design.

References

[1] Dai HJi Y. Statistical analysis of Chinese large-diameter shield tunnel and state-of-art and prospective of comprehensive technologies. Tunnel Construction, 2022, 42(5): 757 (in Chinese)

[2] Broere W . Urban underground space: Solving the problems of today’s cities. Tunnelling and Underground Space Technology, 2016, 55: 245–248

[3] Zhang J , Hu J , Zong C , Feng T , Xu T . A tunneling speed enhancement method for super-large-diameter shield machines considering strata heterogeneity. Tunnelling and Underground Space Technology, 2025, 159: 106496

[4] Do N A , Dias D . Tunnel lining design in multi-layered grounds. Tunnelling and Underground Space Technology, 2018, 81: 103–111

[5] Wang RJia BDeng AWang J. Dynamic triaxial testing study on dynamic characteristics of sandy pebble soil. Chinese Journal of Rock Mechanics and Engineering, 2006, 25: 4059–4064 (in Chinese)

[6] Soroush A , Soltani-Jigheh H . Pre- and post-cyclic behavior of mixed clayey soils. Canadian Geotechnical Journal, 2009, 46(2): 115–128

[7] Dine B S E , Dupla J C , Frank R , Canou J , Kazan Y . Mechanical characterization of matrix coarse-grained soils with a large-sized triaxial device. Canadian Geotechnical Journal, 2010, 47(4): 425–438

[8] Soltani-Jigheh H , Zarei C . A comprehensive experimental study on the monotonic behavior of clay-sand mixtures. Arabian Journal of Geosciences, 2020, 13(13): 540

[9] Huang M , Yu D , Li X , Zhao Y . Study on the mechanical characteristics of sand pebble surrounding rock considering the disturbance effect of tunnel excavation. Frontiers in Earth Science, 2024, 12: 1341922

[10] Horn M. Horizontal earth pressure on perpendicular tunnel face. In: Proceedings of the Hungarian national conference of the foundation engineer industry. Budapest: Academic Publishing House, 1961, 7–16

[11] Anagnostou G , Kovari K . The face stability of slurry-shield-driven tunnels. Tunnelling and Underground Space Technology, 1994, 9(2): 165–174

[12] Zhang ZHu W. Investigation on excavation face support pressure calculation methods of shield tunnelling in clayey soil. Chinese Journal of Rock Mechanics and Engineering, 2014, 33(3): 606–614 (in Chinese)

[13] Dai ZHu Z. Upper bound limit analysis of limit support pressure for shield excavation face in composite ground. Advanced Engineering Sciences, 2021, 53(2): 95–102 (in Chinese)

[14] Kong D , Teng S , Shi J , Zhao M , Deng M . Analytical model for tunnel face stability in composite strata considering the arching effect. International Journal of Geomechanics, 2024, 24(12): 04024277

[15] Song YWang WDu C. Model tests on stability and ultimate support pressure of shield tunnel in sand−gravel composite stratum. Chinese Journal of Geotechnical Engineering, 2020, 42(12): 2206–2214 (in Chinese)

[16] Ma S , Duan Z , Huang Z , Liu Y , Shao Y . Study on the stability of shield tunnel face in clay and clay-gravel stratum through large-scale physical model tests with transparent soil. Tunnelling and Underground Space Technology, 2022, 119: 104199

[17] Cui X , Li P , Ge Z , Li S , Chen Y . Experimental study on face stability of shield tunnel in water-rich inclined composite strata considering different inclination angles. Applied Ocean Research, 2025, 154: 104323

[18] Cui X , Li P , Wang C , Ge Z , Li S , Di Q . Effect of seepage flow on face stability for a tunnel in water-rich silty clay overlying sandy cobble strata. Tunnelling and Underground Space Technology, 2025, 161: 106539

[19] Liang C , Wang Y , Yu H , Zhao C , Zhang Z , Zou H , Zheng S , Guo Y . Model experiment and numerical simulation study on the instability of shield tunnel face in upper-clay and lower-sand composite strata. Environmental Earth Sciences, 2025, 84(7): 174

[20] Wu X , Xu J , Wang S , Sha P , Han Z , Chen X , Shu S , Qiao W , Zeng X . Ground deformation of shield tunneling through composite strata in coastal areas. Buildings, 2024, 14(5): 1236

[21] Mei Y , Zhou D , Shi W , Zhang Y , Zhang Y . Laws and numerical analysis of surface deformation caused by excavation of large diameter slurry shield in upper-soft and lower-hard composite stratum. Buildings, 2022, 12(9): 1470

[22] Wang QLi QZhang YLi J. Study on stability of shield tunnel excavation face with soft upper and hard lower strata. In: Proceedings of the 8th International Conference on Hydraulic and Civil Engineering: Deep Space Intelligent Development and Utilization Forum (ICHCE). Xi’an: IEEE, 2022, 333–337

[23] Xiao HSun Y. Stability analysis of excavation surface of shield tunnel in upper-soft and lower-hard composite strata. China Railway Science, 2022, 43(5): 70–77 (in Chinese)

[24] Zhang X , Song Q , Yao Z , Su D , Zhang Y , Li Q . Limit equilibrium models for tunnel face stability in composite soft-hard strata. Applied Sciences, 2023, 13(19): 10748

[25] Duan Z , Yang J , Ma S , Huang Z . Study on instability mechanism of tunnel face in clay-gravel strata under confined water conditions based on DEM-CFD coupling method. Computers and Geotechnics, 2023, 164: 105769

[26] Wei LYang CLiu LSong QSu D. Research on the stability of the working face of large-diameter shield tunnelling in upper hard and lower soft composite strata. Modern Tunnelling Technology, 2024, 61(5): 21–28,41 (in Chinese)

[27] Zhang H , Xu W , Yu Y . Triaxial tests of soil–rock mixtures with different rock block distributions. Soil and Foundation, 2016, 56(1): 44–56

[28] Liu J , Pan J , Wang B , Hu C , Liu Q . Study on the shear and deformation characteristics of geogrid-reinforced gravelly soils based on large-scale triaxial tests. Frontiers in Earth Science (Lausanne), 2024, 12: 1287718

[29] Cetin H , Gökoğlu A . Soil structure changes during drained and undrained triaxial shear of a clayey soil. Soil and Foundation, 2013, 53(5): 628–638

[30] Chambon P , Corté J F . Shallow tunnels in cohesionless soil: Stability of tunnel face. Journal of Geotechnical Engineering, 1994, 120(7): 1148–1165

[31] Chang Y , Cao P , Zhang J , Fan Z , Xie W , Liu Z , Deng H , Zhao Q . Face stability of tunnel in multi-stratum: Limit analysis and numerical simulation. Geotechnical and Geological Engineering, 2023, 41(5): 3203–3215

[32] Weng X , Dang B , Li X , Fei Y , Ma Y . Study on the instability mode of a tunnel face under variable seepage conditions in sandy soil shield tunnels: Centrifuge tests and numerical simulation. Tunnelling and Underground Space Technology, 2025, 159: 106515

[33] Chen X , Hu Y , Cheng P , Fang P , Yao K . Effect of soil arching evolution on tunnel face stability considering spatially variable sandy soils. Transportation Geotechnics, 2025, 51: 101496

[34] Cui J , Yao Z , Yu T , Wang J , Ying K , Liu B , Zhu S , Yan X . Soil displacement of slurry shield tunnelling in sandy pebble soil based on field monitoring and numerical simulation. Buildings, 2024, 14(10): 3043

[35] Chappell B A . Deformational response in discontinua. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts, 1979, 16(6): 377–390

[36] Chen C N , Huang W Y , Tseng C T . Stress redistribution and ground arch development during tunneling. Tunnelling and Underground Space Technology, 2011, 26(1): 228–235

[37] Chen R , Song X , Meng F , Wu H , Lin X . Analytical approach to predict tunneling-induced subsurface settlement in sand considering soil arching effect. Computers and Geotechnics, 2022, 141: 104492

[38] Liu M , Meng F , Chen R , Cheng H , Li Z . Numerical study on the lateral soil arching effect and associated tunnel responses behind braced excavation in clayey ground. Transportation Geotechnics, 2023, 40: 100970

[39] Khandouzi G , Khosravi M H . An analytical investigation of soil arching induced by tunneling in sandy ground. Tunnelling and Underground Space Technology, 2023, 140: 105242

[40] Lin X , Chen R , Wu H , Cheng H . Three-dimensional stress-transfer mechanism and soil arching evolution induced by shield tunneling in sandy ground. Tunnelling and Underground Space Technology, 2019, 93: 103104

[41] Song X , Meng F , Chen R , Wu H . Effect of seepage on soil arching effect in deep shield tunnel. Underground Space, 2023, 12: 218–233

RIGHTS & PERMISSIONS

Higher Education Press