Soil arching due to leaking of tunnel buried in water-rich sand

Tunnelling and Underground Space Technology 95 (2020) 103158

Contents lists available at ScienceDirect

Tunnelling and Underground Space Technology journal homepage: www.elsevier.com/locate/tust

Soil arching due to leaking of tunnel buried in water-rich sand Ying-Ying Long, Yong Tan

⁎

T

Department of Geotechnical Engineering, College of Civil Engineering, Tongji University, 1239 Siping Road, Shanghai 200092, China

A R T I C LE I N FO

A B S T R A C T

Keywords: Tunnel leaking FDM-DEM coupling analysis Water-rich sand Seepage erosion Stress arching

For tunnel buried in water-rich sand stratum, water intrusion into tunnel lining would bring sand particles into the tunnel and then caused ground subsidence or even cave-in of roadway (sinkhole). Until now, limited studies have been contributed to this topic. To investigate the geo-hazards and explore the associated failure mechanism due to tunnel leaking, both experimental test and comprehensive numerical simulations using the validated finite-difference-method and discrete-element-method (FDM-DEM) coupling method were carried out in this study. It was disclosed that displacements of the sand stratum surrounding the leaking tunnel could be divided into three regions (erosion, displacement, and static zones) by two ellipsoids; the stress changes caused by seepage erosion around the leaking tunnel were similar to that observed by the active door test. Finally, a stress arching model accounting for seepage erosion around a leaking tunnel was proposed. This model can help professionals evaluate stability of overburdening strata by giving soil stress distributions in sand formation surrounding a leaking tunnel with seepage erosion.

1. Introduction In recent years, sudden ground subsidence or cave-in of roadway (sinkhole) without obvious sign in advance (Fig. 1) have frequently happened in China (Wu et al., 2017a; Lyu et al., 2019) and worldwide (Davieset al., 2001; NILIM, 2006; Chai et al., 2018). These accidents led to adverse or detrimental effects on urban environment and even casualties sometimes. The post-failure investigations (Korff et al., 2011; Ni and Cheng, 2012; Sato and Kuwano, 2015; Tan and Lu, 2017; Chai et al., 2018) disclosed that many of these accidents occurred in waterrich sand and were closely associated with water intrusion into underground structures through openings (e.g., tunnel leaking; throughwall leaking). As a result, fine particles were eroded by seepage water and then flowed into the openings of underground structures, which severely disturbed the strata and eventually incurred abrupt ground subsidence or even sinkhole (cave-in of roadway). For tunnels that have been in service for many years, openings on tunnel linings were generally caused by open dislocation of joints due to differential settlement of tunnel segments or by failure of waterproofing at the bolt holes. Loss of soil particles due to tunnel leaking could not only incur ground subsidence but also change stress distribution in ground. As observed in the trapdoor tests (Terzaghi, 1943; Adachi et al., 2003; Chevalier et al., 2012; Rui et al., 2019), loss of particles from the opening produced soil arching effect on the soil layers above the opening. Soil arching is common in underground engineering, which governs the stability of

⁎

ground (Lai et al., 2014; Wu et al., 2017b, 2019a). Because it is closely related to the stability of tunnel excavation face and surrounding strata, soil arching in tunnel excavation has already been extensively studied (e.g., Lee et al., 2006; Chen et al., 2013; Kong et al., 2018). In contrast, soil arching due to tunnel leaking has received much less attention and few of relevant studies have been known in literature. In light of this, it is worthwhile to explore the mechanism of soil arching development due to tunnel leaking; the relevant results will be useful references for timely assessing associated adverse influence and then preparing or adopting appropriate countermeasures or remedial measures. At present, some studies have been devoted to seepage erosion of soil particles, but they mainly focused on internal erosion (Kenney, 1985; Sterpi, 2003; Ke and Takahashi, 2014; Ouyang and Takahashi, 2016). Internal erosion defines that fine particles are carried away through the voids among coarse particles by seepage flow, which forms a continuous flow path; thus, fine particles continuously get lost and eventually leave soil skeleton (Ke and Takahashi, 2012). Distinct from the internal erosion, movement of soil particles due to tunnel leaking (seepage erosion) includes not only suffusion of fine particles in coarse particles, but also erosion of soil at the opening. Therefore, the severity of seepage erosion is related to not only particle grading and hydraulic gradient, but also the size of the damaged spot (opening) on tunnel. Furthermore, the force of fluid acting on soil particles is directed to the damaged spot on tunnel, which causes flow path of soil particles to be different from that of the internal erosion. Currently, some researchers

Corresponding author. E-mail addresses: [email protected] (Y.-Y. Long), [email protected] (Y. Tan).

https://doi.org/10.1016/j.tust.2019.103158 Received 11 June 2019; Received in revised form 18 October 2019; Accepted 19 October 2019 0886-7798/ © 2019 Elsevier Ltd. All rights reserved.

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) Lanzhou City, Gansu Province, China

(b) Dongguan City, Guangdong Province, China

(c) Zhengzhou City, Henan Province, China

Fig. 1. Typical accidents of cave-in of roadways above tunnels in China.

diameter). The water level was maintained at the ground level during testing. Since ground water is semi-infinite for buried tunnels, it was reasonable to assume that the water level was constant during the entire process of seepage erosion. To achieve a dense state of the sand in the experimental test, the following steps were adopted. First, compaction tests were carried out on the sand samples and an optimum water content of 10% for soil densification was determined. Then, the sands were mixed with distilled water at this water content and then placed into the model box layer by layer (10 cm thick per layer); each layer was compacted densely before placing another layer above. For facilitating observation on loss path of soil particles (Fig. 2b), colored sands were placed at three different elevations of the stratum. Once the sands were filled to the required height, distilled water was added into the two water tanks on the two sides of the model box to saturate the sands in the model box until the water flooded the soil surface. This process took approximately 48 h to ensure full saturation of the sands. Before the tunnel leaking test, the water above the soil surface was drained to keep the water line steady at the surface. The maximum void ratio, the minimum void ratio, and void ratio of the sands used in this test were 0.83, 0.52 and 0.61, respectively; then, a relatively density of 70% can be obtained for the sands in the model box. Once the test was completed, the images were processed using an image processing program developed by Lu et al. (2018a, 2020).

have studied seepage erosion due to leaking of buried pipelines (Mukunoki et al., 2009; Guo et al., 2013; Tang et al., 2017). However, these studies just focused on seepage velocity and particle loss, without considering the effect of soil particles loss on the stress of ground. As for a leaking tunnel buried in water-rich sand, uneven ground subsidence due to particle erosion would lead to redistribution of stress in the ground and then affect the stability of overburdening strata. Thereby, it is essential to investigate stress arching resulting from leaking of tunnel. In order to explore stress changes in water-rich sand formation due to tunnel leaking, both experimental model test and FDM-DEM coupling numerical simulations were conducted in this study. Based on the displacement development law obtained by the model test, the rationality of the 2D FDM-DEM coupling method adopted in this study was validated first. Thereafter, comprehensive parametric studies using the validated FDM-DEM coupling method were carried out to analyze stress changes around a leaking tunnel and development of soil-arching following the process of seepage erosion. The arching mechanism obtained in this study was called as stress arching thereafter. Then, a stress arching calculation method derived from the 2D simulations was incorporated into the 3D FDM-DEM coupling simulations. Finally, mathematical equations describing stress arching in water-rich sand around leaking tunnel were given from both the 2D and 3D perspectives. These equations will be helpful for engineers to evaluate ground stability associated with tunnel leaking in water-rich sandy strata.

2.2. Experimental test results 2. Experimental test at 1 g The recorded seepage erosion process of the sand due to tunnel leaking is shown in Fig. 4, which can be divided into three distinct stages. The first stage (stage 1) is the initial erosion stage (Fig. 4a-b), which lasted for about 12 s. During this period, ground loosening zone mainly occurred within a certain range above the model tunnel, and did not reach the ground surface. Thereafter, the loosening zone gradually expanded to the surface and apparent subsidence at the ground level was observed (stage 2), refer to Fig. 4 (c-d). In stage 2, two inclined cracks appeared near the ground surface, probably due to the friction between the model box and the sand particles. Meanwhile, a fixed seepage area was formed above the tunnel and the soil particles within this area continuously got lost, which resulted in aggravation of surface subsidence. After a certain period of time, the particles near the surface gradually stabilized and intrusion of soil particles into the tunnel stopped (stage 3), refer to Fig. 4(e–h).

2.1. Test procedures The model box used to simulate tunnel leaking is schematically illustrated in Fig. 2. In order to observe the seepage erosion process during tunnel leaking, one side wall of the box consisted of transparent plexiglass. Two water tanks were placed on the two sides of the model box for controlling water level. Tunnel in practice has a large diameter, usually around 6 m for metro tunnels. To observe the development law of seepage erosion due to tunnel leaking, the tunnel diameter was reduced by 25 times and replaced by a plexiglass tube with a diameter of 0.2 m in this experimental test. Since only qualitative observations were made, the experimental test was performed at 1 g. An opening of 10 mm × 50 mm was cut on the crest of the plexiglass pipe to simulate cracks at the joints of the tunnel segments. The pipe was buried inside the model box, 200 mm from the bottom of the model box. Two cameras were used to record the seepage erosion process, which were placed in front of and above the top of the model box, respectively (Fig. 2a). Since this test just focused on the water-rich sand stratum, Fujian standard sand was used as the soil medium. The particle size distribution of the sand is continuous, refer to Fig. 3. The sand layer placed in the model box was 0.6 m thick, which was 0.3 m above the crest of the tunnel (the buried depth of the model tunnel was 1.5 times the tunnel

2.2.1. Stage 1 The image processing program developed by Lu et al. (2018a, 2020) was used to process the photos taken at different times for obtaining the displacement field of the soil particles at stage 1 and the vector diagram in Fig. 5 illustrates the corresponding particle displacements. As shown in this figure, particles displacement mainly occurred within a coneshaped zone directly above the tunnel. The particles in this zone tended to flow into the tunnel through the opening. As time went by, the loose 2

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) Three-dimensional schematic (unit: mm)

(b) Model box Fig. 2. Model devices used in this study.

area began to expand and particles near the ground level moved faster than particles near the tunnel opening. This was an obvious sign that the ground surface was about to settle.

VN =

4π AN BN2 3

(1)

where AN = the maximum semi-axis length of the ellipsoid; BN = the minimum semi-axis length of the ellipsoid; the horizontal section of the ellipsoid was assumed to be circular. The eccentricity of the ellipsoid can be expressed as:

2.2.2. Stage 2 In stage 2, cracks appeared on the ground surface, marking the beginning of surface subsidence. Subsequently, an elliptical erosion zone was formed over the tunnel. The soil near the ground surface slipped into the elliptical area, and then flowed into the tunnel through the tunnel opening. This process was similar to the flow process of the granules in the silo (Grudzien et al., 2017; Fullard et al., 2018). Janelid and Kvapil (1966) proposed a concept of flow ellipsoid for loss of granular materials in hopper or silo. As illustrated in Fig. 6(a), when an opening at the chamber bottom was opened, the particles in the chamber started to flow outward. After a certain period of time, all particles flowed out of this approximate ellipsoidal region which was called as the ellipsoid of motion. For the particles located between the ellipsoid of motion and a corresponding limit ellipsoid, they loosened and displaced. While, particles outside the limit ellipsoid remained stationary (Brady and Brown, 2004). The volume of the ellipsoid of motion can be calculated by

ε=

1 (AN2 − BN2 )1/2 AN

(2)

Janelid and Kvapil (1966) pointed out that in practice ε varies between 0.90 and 0.98 in most cases. By assuming VN equal to the volume EN of the outgoing particle, the semi-axis length of the ellipsoid can be obtained and then was substituted into the eccentricity within the range of values:

BN =

3

EN (1 − ε 2) 3π /4

(3)

Then, a limit ellipsoid corresponding to this ellipsoid can be obtained. Janelid and Kvapil (1966) had given a relationship between this limit ellipsoid and the ellipsoid of motion using a parameter β: 3

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

really simulate a seepage erosion problem. Discrete element method DEM is based on discrete mechanics, which was first proposed by Cundall (1971) for the application of rock mechanics and later applied to the study of soil mechanics (Cundall and Strack, 1979). Combined with computational geometry, DEM can reasonably simulate the true shape of particles which enforces interlocking forces as well as frictional forces among particles (Lu et al., 2017). To date, it has been successfully used for analysis of slope stability (Lu et al., 2018b) and soil-structure interaction (Tan et al., 2018) in soil-rock mixtures. Additionally, seepage erosion is a process of complex interaction between fluids and particles; hence, the combination of FDM and DEM enables simultaneous simulation of fluids and solids. The coupling FDM-DEM is an ideal method to simulate the seepage erosion of sand particles and has been adopted in some studies (Cai et al., 2007; Saiang, 2010; Itasca Consultant Group, 2012; Ma et al., 2018; Zhang et al., 2019). Based on the considerations above, comprehensive numerical simulations using a FDM-DEM coupling method were carried out in the following sections. In order to ensure the reliability of the numerical simulations, the adopted FDM-DEM coupling method has been validated by the experimental results prior to performing extensive parametric studies.

Fig. 3. Particle grading curve of the sand used in the experimental test.

β=

VG VG − VN

3.1. Formulation of the coupled numerical model (4)

The DEM analysis was performed with the commercial DEM software PFC5.0 program (Itasca Consulting Group, 2015). The soil was composed of a certain number of rigid particles, the movement between the particles was independent of each other, and the force between the particles was transmitted through the contact points between the particles. The Newton's second law was satisfied between the force or moment of the particle and the particle motion, i.e., the motion of the particle was generated by the contact force and the physical force acting on the particle. The contact force between the particles followed the force-displacement law, which was used to update the contact force generated by the relative motion at the contact (Itasca Consulting Group, 2015). Many studies have demonstrated that linear contact models can simulate the mechanical behavior among sand particles (Minh and Cheng, 2013; Dai et al., 2014; Lopez et al., 2016). For the linear contact model, the moment is always equal to 0, and the forcedisplacement law updates the contact force as follows:

where, VG is the volume of the limit ellipsoid. In the process of seepage erosion, the opening on the model tunnel can be regarded as a discharge point, through which soil particles got lost. The particles in the vicinity of the opening were displaced, while the particles outside the limited range remained stationary. According to the method above, the flow ellipsoid around the tunnel can be roughly estimated for the experimental model test, refer to Fig. 6(b). The two ellipsoids divided the particles around the model tunnel into three regions, namely erosion, displacement and static zones. “Erosion” means that particles inside the ellipsoid of motion got lost through the opening on the tunnel; “Displacement” refers to looseness of the particles between the ellipsoid of motion and the limit ellipsoid; “Static” means that particles outside the limit ellipsoid remained stationary. 2.2.3. Stage 3 Once the stratum was severely disturbed, it evolved into stage 3. The particles around the model tunnel reached a new stable state and did not suffer erosion any longer. Meanwhile, groundwater was still overflowing. As the soil particles have sealed off the tunnel opening, no intrusion of soil particles into the tunnel took place and the exposed groundwater level went upwards to the surface rapidly. At the end of the seepage erosion, the front of the erosion pit (Fig. 4g) showed an inverted triangle profile with equal slopes on either side; the tunnel was not exposed to air or water and there still existed certain thickness of soils overlying the tunnel. It implies that there may be a soil arching formed above the tunnel. From the top view of the test result in Fig. 4(h), it can be seen that the slump featured a semicircle shape in plane, and the diameter of the ground subsidence caused by tunnel leaking was about twice the tunnel diameter.

Fc = F l + F d

(5)

Fl

Fd

where Fc is the contact force, is the linear force, and is the dashpot force. When considering the fluid force, the equation of motion of the particle can be expressed as:

→ → fmech + ffluid → u ∂→ = + g m ∂t

(6)

→ where → u is the particle velocity, m is the particle mass, fmech is the sum → of additional forces acting on the particle, ffluid is the total force applied → by the fluid on the particle, and g is the acceleration of gravity. The total force of the fluid acting on the particles consists of two parts: drag force and hydrodynamic forces. Related studies (Morsi et al., 1972; Zhu et al., 2007) have shown that the hydrodynamic forces are negligible compared to the resistance between the fluid and the particles. Therefore, the drag force of the fluid as the only force acting on the particles can be calculated as below:

3. FDM-DEM coupling simulations From the experimental test results, only the development law of particle displacement can be obtained; while the relevant mesoscopic mechanism cannot be explored. Compared with the experimental test, numerical simulation is an effective and convenient method to explore the problem of seepage erosion. The process of particle seepage erosion is a highly nonlinear behavior accompanied with loss and migration of soil particles. Inherently, it is difficult for a continuous medium method (e.g., finite element method – FEM; finite difference method – FDM) to

Fi = Vi ∙γw i

(7)

where Vi is the volume of the particles, γw is the gravity of the water, and i is the hydraulic gradient. i = qi / k , where, qi=the specific discharge and k=the permeability coefficient (Itasca Consultant Group, 2012). Generally, it is difficult for the CFD (computational fluid dynamics) module that comes with the PFC program to simulate a circular tunnel 4

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) t = 0 s

(b) t = 10s

(c) t = 18 s

(d) t = 60 s

(e) t = 90 s

(f) t = 110 s

(g) Front view at end of erosion

(h) Top view at end of erosion

Fig. 4. Process of seepage erosion observed in the experimental test.

where Cs is the shape factor of granular material; Ss is the surface area per unit volume of soil solids; T is the tortuosity factor; μ is the dynamic viscosity; and e is the void ratio. In porous media containing approximately uniform pore sizes, Cs can be taken as 2.5, and T can be taken as 2 (Yao et al., 2018). For spherical particles, Ss can take 6/ dp (Tang et al., 2017), where, dp is the particle diameter. Porosity can be converted to a void ratio by a formula, which can be measured in PFC; then, formula (8) can be converted to:

boundary, because it can only support a coarse mesh composed of tetrahedrons or hexahedrons. In contrast, FLAC, which belongs to the same company as PFC, can realize modeling and calculation of circular boundary (2D) or cylindrical boundary (3D). Moreover, FLAC is a software based on finite difference method (FDM), which can simulate flow field change in the process of seepage erosion. Based on these considerations, FLAC was chosen for seepage analysis. By FLAC-PFC coupling analysis, the seepage erosion in sandy strata due to tunnel leaking can be simulated reasonably. FLAC models the flow of fluid by Darcy’s law, an effective method for predicting fluid behavior. When FLAC only performs fluid calculation, the key parameter is the permeability coefficient of the soil, k . It is generally obtained by the method of porosity conversion for granular material, the most common of which is the Kozeny-Carman equation (Carman, 1956; Kozeny, 1927):

k=

γw e3 1 Cs Ss2 T 2 μ 1 + e

k=

γw 1 n3 2 2 Cs Ss T μ (1 − n)2

(9)

where n is a void ratio. The FDM-DEM coupling process by PFC and FLAC can be described briefly as the following major steps. First, the corresponding stratigraphic model was established in each of PFC and FLAC. Then, the void ratio in the PFC was converted to a permeability coefficient and passed to FLAC via the socket Input / Output (I/O). After FLAC ran stably, the obtained flow field was converted into seepage force of the particles

(8) 5

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) t = 1 s

(b) t = 5 s

(c) t = 9 s

(d) t = 12 s

Fig. 5. Particle displacement vectors in stage 1. Note: The red dash line in each figure is used to divide the area, whose location is dependent on whether the particles move or not. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

accordingly to ensure the reasonability of the numerical simulation for tunnel leaking. Thereby, a tunnel opening of 15 mm corresponding to 10 mm in the experimental test was used in the DEM model. In this paper, the linear contact model was chosen as the contact model between particles. The main model parameters were particle density, friction coefficient, damping, normal stiffness, shear stiffness and friction coefficient (Table 1). The particle contact stiffness and friction coefficient were obtained by biaxial-test calibration. Two damping mechanisms were set. One was local damping, which accelerated dissipation of kinetic energy by applying a damping force with magnitude proportional to the unbalanced force to each particle; the other was viscous damping, which generated energy consumption by setting a damping opposite to the direction of motion at the contact. Local damping was used for balance of numerical models and viscous damping was used to simulate energy dissipation of particles in seepage erosion process. Once DEM model was established, the void ratio of the model was derived and converted to the permeability coefficient by Eq. (9). A mesh model was built in FLAC, and the permeability coefficient was introduced via socket O/I to start seepage analysis. During analysis, the water line was set to be constant on the model surface, and the pressure at the tunnel opening was set to 0 for achieving leakage at the tunnel opening. The hydraulic gradient of the cell grid in FDM was extracted, which was converted into seepage force by Eq. (7) and then imported into DEM model. As shown in Fig. 7, the arrows represented the seepage forces of the particles, and the sizes of the arrows represented the magnitudes of seepage forces.

and transmitted to PFC via the socket I/O. Finally, PFC updated the force of particles and began a simulation of the seepage erosion process. 3.2. Validation of FDM-DEM coupling model by experimental test 3.2.1. Numerical model In the previous experimental test, a rectangular opening was used to simulate the tunnel opening. The tunnel itself extended in the longitudinal direction; thereby, the seepage erosion problem due to tunnel leaking can be treated as a plane-strain problem. It can be simulated by a two-dimensional (2D) FDM-DEM coupling method. A numerical model of the same size as experimental test was established, and a circle was used to simulate the buried tunnel, see Fig. 7. The difference in particle colors in Fig. 7 was to separate the particles and facilitate observation of particle loss. If soil particles were established according to their actual sizes, a huge amount of particles would be generated. This would make a computer unable to perform calculation. In light of this, a representative particle-size range (0.2–1.3 mm) for the sands used in the experimental test was adopted for numerical simulation. In this range, only 25% sand were finer than 0.2 mm or greater than 1.3 mm. Considering huge amounts of particles to be generated in the DEM program, it was necessary to enlarge particle size of a certain multiple. By this way, the calculation rate of DEM can be improved significantly. Based on this consideration, the soil particles were enlarged by 1.5 time in the DEM model, i.e., the particle sizes distributed between 0.3 and 2.0 mm corresponding to 0.2–1.3 mm in the experimental test were adopted in the numerical model. When the particle sizes were enlarged, the corresponding size of the tunnel opening had to be enlarged 6

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) Flow ellipsoid concept (redrawn based on Janelid and Kvapil 1966)

(b) Flow ellipsoids in the experimental test

(c) Flow ellipses in the 2D numerical model Fig. 6. Flow ellipsoids/ellipses.

7

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

Fig. 7. 2D FDM-DEM coupling model (unit: mm).

particles got lost layer by layer under dragging of groundwater and gravity; as a result, ground subsidence began to take place. Eventually, the particles formed a certain slope and remained stable, and the angle of inclination was smaller than the inclination of the experimental test when it was stable. This discrepancy was probably due to that the rolling resistance at particle contact had been ignored in the DEM simulation, which is normal in a DEM simulation (Jiang et al., 2006; Lu et al., 2018c). With the beginning of tunnel leaking, soil particles directly above the opening began to flow into the tunnel; with the progress of leaking over time, soil particles flowed into a particular elliptical region and then flowed into the tunnel; as a consequence, significant ground subsidence were tracked at the ground level. According to the concept of flow ellipsoid, flow ellipses similar to those of the experimental test can be obtained in Fig. 6(c). Although these comparisons were qualitative rather than quantitative, it demonstrated that the adopted FDM-DEM coupling simulation can well restore the displacement development law of the particles observed in the experimental test.

Table 1 Parameters adopted in numerical simulation. Soil parameters Particle density (Kg/m3) Friction coefficient Cohesion (kPa) Particle normal stiffness (N/m) Particle shear stiffness (N/m) Damping Normal critical damping ratio Porosity of granular material Wall normal stiffness (N/m) Wall shear stiffness (N/m) Wall friction Fluid parameters Density (kg/m3) Fluid bulk modulus (Pa)

2650 0.50 0.00 1.00 × 108 1.00 × 108 0.70 0.20 0.30 1.00 × 108 1.00 × 108 0.50 1000 2.00 × 109

3.2.2. Validation of FDM-DEM coupling model Fig. 8 illustrates typical seepage erosion processes of sand particles in this numerical simulation at different time steps. As shown in this figure, the direction of particle loss in the numerical simulations well captured the trends of particle loss in the experimental test (Fig. 4). The

(a) Initial erosion

(b) Particles erosion layer by layer

(c) Particles erosion layer by layer

(d) End of erosion

Fig. 8. Process of seepage erosion in the numerical simulation. 8

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) Phase 1

(b) Phase 2

(c) Phase 3

(d) Phase 4 Fig. 9. Development of contact force chains (left column) and particle displacements (right column). Note: The black and red balls in (d) represent the particle positions in phase 1 and phase 4, respectively. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

9

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) C/D = 1.5

(b) C/D = 2.0

(c) C/D = 2.5

(d) C/D = 3.0

(e) C/D = 3.5

(f) C/D = 4.0

Fig. 10. Distributions of vertical stresses above the tunnel. Note: he = q/(γD) where, q is vertical stress, γ is soil unit weight, and D is tunnel diameter.

4. Soil arching

The model boundary was a rectangle of 8-m long and 10-m wide. A circle with a diameter of 1 m was used to simulate the tunnel with a buried depth ratio, C/D, of 1.5 to 4.0, where C denotes buried depth of tunnel and D denotes tunnel diameter. The water level was maintained at the ground level during the numerical simulation. According to the FDM-DEM coupling method described previously, a seepage erosion test was started. After the seepage force was introduced, the soil began to be eroded and the particles were washed into the tunnel. The four figures on the left side of Fig. 9 are the diagrams of the force chains of contact forces between particles, and the thicknesses of the force chain curves were proportional to the magnitudes of contact forces. The four figures on the right side of Fig. 9 are

4.1. Force chain during seepage erosion process In order to study the arching effect above the leaking tunnel after the seepage erosion ended, a series of tunnel-stratigraphic models with different buried depths of the tunnel were established. Considering the limited computational power of a computer, it would be difficult to establish a full-scale numerical model. Based on the similarity principle of the geotechnical centrifuge test, the numerical model was reduced by a factor of five while applying a gravitational acceleration of 5g, so that the initial stress level of the strata was close to the actual stress level. 10

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

Fig. 11. Natural balance arch (redrawn based on Miao, 1990).

disregarding the effects of soil relative density, as the buried depth of the tunnel increased, the size of the stress arching increased accordingly. The numerical simulation of tunnel leakage was analogous to the trapdoor test. The local opening in the tunnel was similar to the trapdoor in the trapdoor test. The sand particles were gradually lost from the opening or the trapdoor and then tended to stabilize under the friction between the particles. Through the trapdoor tests, Terzaghi (1943) first discovered the soil arching effect. With reference to the idea of the trapdoor test, the vertical stress and stress arching above the leaking tunnel opening were explored herein. The vertical stress above the tunnel after seepage erosion was measured by the measuring circles in the PFC and then converted to the equivalent depth by Eq. (13):

the displacement maps of the particles. The color of the particles varied from red to blue, representing the displacement of the particles from large to small. When the seepage erosion began, the force chains at the tunnel opening were broken by the particles and fluid (Fig. 9a), indicating loss of particles. At some points, new archings of force chains were formed above the opening (Fig. 9b). However, if the arching of force chains was not stable enough to withstand the impact of the upper particles and fluid, the new arching of force chains would be broken again (Fig. 9c), and then particles continued to get lost. With the formation of arching and repeated occurrence of breakage, the amount of particle loss increased. Once the loosening area of the particles expanded to the ground surface, the strata collapsed suddenly. If a new stable arching of force chains (Fig. 9d) sufficient to resist the upper impact force was formed above the tunnel opening, the particles stopped escaping and the seepage erosion ended completely. In order to clearly observe the soil arching phenomenon above the tunnel, only highly contact forces larger than 0.04Fmax were shown in Fig. 9(d), where, Fmax = the maximum contact force in the entire model. It was divided into two zones, depending on whether the directions of the contact forces were disturbed or not. The directions of the force chains in zone I rotated and even became horizontal when they were above the tunnel; adjacent contact forces were connected to each other to form an arching of force chains. The force chains in Zone II were almost undisturbed, and their directions were still mostly vertical. Since the stability of the ground was mainly related to the position of the outermost arching (the red line in Fig. 9d), the particles on this arching (red balls) are marked in the particle displacement in Fig. 9(d). Meanwhile, the positions of these particles before the seepage erosion (black balls) were also marked on the same diagram, in which the white arrows indicated the directions of particle movement. It can be seen that the particles were distributed on the periphery of the arching before the seepage erosion; with the beginning of seepage erosion, the particles at different locations started to flow towards the tunnel opening; eventually, the particles were joined together to form a stable arching of force chains.

he =

q γD

(13)

where q= vertical stress, γ =unit weight of soil, D= tunnel diameter, and he= an equivalent burial depth corresponding to q . As shown in Fig. 10, the different line types represented vertical stresses tracked at different depths (ha ) in terms of D . It can be seen that the vertical stress above the tunnel center was greatly reduced while the stresses on the left and right sides of the tunnel were slightly increased and the vertical stress near the surface was almost unchanged. Within a certain zone above the tunnel, the vertical stresses were convexly distributed in a symmetric manner. In the trapdoor test, the vertical stresses at different depths in the loose area above the trapdoor were also observed to be convexly distributed (Adachi et al., 2003). This consistency indicated occurrence of a soil arching above the tunnel after seepage erosion. In this study, stress arching was used to describe the soil arching effect in the tunnel leaking. It was an arching mechanism based on the distribution of contact force chains between particles. The shape of the stress arching was determined according to the principle of natural equilibrium arch proposed by Protodyakonov (1907). He pointed out that after a deep tunnel was excavated, a natural balanced soil arching would be formed in the surrounding rock above the tunnel and the edge of the balance arch was the arch with a thin thickness like a grain. This opinion was based on two assumptions: (1) there was no bending moment on any section and (2) the arch could remain stable without slipping. According to the theory of the natural balance arch, the natural balance arch formed under the confining pressure of longitudinal pressure q and lateral pressure λq can be drawn in Fig. 11 (Miao, 1990), where λ is the lateral earth pressure coefficient. On the arbitrary section

4.2. 2D stress analysis The previous analyses showed that in the process of seepage erosion, the force chains formed some archings above the tunnel, an obvious evidence of soil arching. As pointed out by Franza et al. (2019), the arching above the tunnel was not only affected by C/D but also by soil relative density. The latter has not been considered in this study. If 11

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

the half-axis length of the ellipse in the Y direction. The critical buried depth (boundary) between the uniform and the non-uniform distributions of the vertical stresses in Fig. 10 was taken as the vault of the stress arching. According to the arching diagram of force chains (Fig. 9d), the principal contact force between the soil particles and the tunnel lining closest to the opening was determined as the intersection of the stress arching and the tunnel boundary. The lateral earth pressure coefficient, λ, could be obtained from the ratio of the horizontal stress to the vertical stress above the tunnel, which was tracked by the measuring circles in PFC. When λ was obtained, the expression of the stress arching was calculated by inputting the coordinates of the tunnel vault and the intersection point of the stress arching and the tunnel boundary into Eq. (16). The lateral earth pressure coefficient and expression of the stress arching at different buried depth ratios are summarized in Table 2. As shown in Table 2, the magnitudes of λ measured at different C/D were greater than 1.0, which were different from the case where the atrest lateral earth pressure coefficient is less than 1.0. This discrepancy was due to that during seepage erosion, the particles got lost along the opening, which caused the vertical stress to decrease and expand downward and then resulted in reduction in the horizontal restraining force. The degrees of reductions in the horizontal stress and the vertical stress were different; thereby, the principal stresses inside the soil rotated and the horizontal stress was greater than the vertical stress. This finding matched with the previous studies in literature (Terzaghi, 1936; Hewlett et al., 1988; Chevalier et al., 2012; Vivanco et al., 2012). The outermost archings of force chains at different C/D were drawn and then compared with the calculated stress archings in Fig. 12. It can be noticed that the calculated stress archings and the outermost archings of force chains almost coincided. In case the effects of relative soil density were not considered, as the buried depth of the tunnel increased, the size of the stress arching increased accordingly. It implies

Table 2 The lateral earth pressure coefficients and expressions of the stress arching. C/D

Lateral earth pressure coefficient

expression of the stress arching

1.5

1.286

(y − 0.275)2 0.5252

+

2.0

1.291

(y − 0.330)2 0.5702

+

2.5

1.272

(y − 0.490)2 0.7102

+

3.0

1.673

(y − 0.808)2 0.9922

+

3.5

1.529

(y − 0.960)2 1.1402

+

4.0

1.646

(y − 1.112)2 1.2882

+

x2 0.5952 x2 0.6472 x2 0.8002 x2

1.2842 x2 1.4522 x2 1.6532

=1 =1 =1

=1 =1 =1

of a natural balance arch, there only existed axial force without bending moment or shearing force. According to the symmetry, one quarter of the ellipse was intercepted for analysis with ∑ Mm = 0 :

λqb2 (b2 − y ) −

1 2 1 qx − λq (b2 − y )2 = 0 2 2

(14)

Eq. (14) can be simplified as below:

y2 x2 + =1 b22 λb22

(15)

In this test, the position of the coordinate axis moved from the center of the ellipse to the center of the tunnel (Fig. 9d), and the expression of the arching was converted into the following formula:

(y − m)2 x2 + =1 b22 λb22

(16)

where m is the distance from the center of the ellipse to X-axis; and b is

(a) C/D = 1.5

(d) C/D = 3.0

(b) C/D = 2.0 (a)C/D= 1.5

(b)C/D= 2.0

(c)C/D= 2.5

(d)C/D= 3.0

(e)C/D= 3.5

(f)C/D= 4.0

(e) C/D = 3.5

(c) C/D = 2.5

(f) C/D = 4.0

Fig. 12. Stress arching and the largest arching of force chain with C/D = 1.5–4.0. Note: The red-color dash lines represent the calculated stress arching, and the black-color bold lines represent the outermost archings of force chains. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) 12

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

stress arching can be described by an ellipsoid, whose expression is:

y2 (z − m)2 x2 + + =1 2 2 b2 λ x b2 λ y b22

(17)

where m is the distance from the center of the ellipsoid to X-axis, b is the half-axis length of the ellipsoid in the Z direction, λ x is the lateral earth pressure coefficient in the X direction, and λ y is the lateral earth pressure coefficient in the Y direction. The vault of the stress arching was obtained from the critical buried depth (boundary) of the uniform and the non-uniform distributions of the vertical stresses in the 3D numerical model. Then, the stress along the longitudinal direction of the tunnel was measured, and the intersection points of the stress arching and the tunnel were assumed to be the locations featuring the maximum contact force between the particles and the tunnel. Finally, by measuring λ x and λ y above the tunnel, the expression of the 3D stress arching was obtained. In this research work, the 3D stress archings after the seepage erosion were analyzed at three different buried depths, namely C/D = 1.0, 1.5 and 2.0. The water level was also maintained at the ground level during numerical simulation. The vertical stress distribution above the tunnel after the seepage erosion was measured and plotted with contour lines of the stress distribution within the longitudinal and the transverse planes to the tunnel axis, as shown in the four figures on the left side of Fig. 14. The four figures on the right side of Fig. 14 are the calculated 3D stress arching. The definitions of X-axis, Y-axis and Z-axis are shown in Fig. 14. It can be seen that after the seepage erosion, the vertical stress above the tunnel opening was significantly reduced, and the vertical stresses on the left and right sides of the tunnel were slightly increased. Compared with the arching shape formed by the vertical stress distribution, the 3D stress arching featured a shape in an oblate spheroid. The trend of the 3D stress arching was consistent with the trend of the vertical stress distribution. When the tunnel depth was increased, the volume of the 3D stress arching increased accordingly.

(a) 3D numerical model in DEM

(b) Fluid force on the particles Fig. 13. 3D FDM-DEM coupling numerical model for simulation.

4.4. Discussion that when the buried depth was greater, the strata were more stable and the damage caused by tunnel leaking was smaller.

4.4.1. The difference between 2D and 3D numerical simulations In this study, both 2D and 3D numerical simulations have been performed to investigate the soil-arching arising from tunnel leaking in the water-rich sand stratum. It should be pointed out that the forms of damage to the tunnel were different in the 2D and 3D numerical models. 2D numerical simulation is suitable for seepage erosion phenomena in the case of narrow cracks along the longitudinal tunnel alignment. To ensure the consistency of the 2D and the 3D numerical results, a 3D numerical model with a narrow crack on the tunnel crest extending in the longitudinal direction was established in Fig. 15(a) for analysis. After the seepage erosion, the particle displacement on the X-Z plane was extracted for Y/D = 1.0, 1.5, and 2.0 (Fig. 15b–d). Obviously, the loss of particles in terms of particle displacements on each tunnel section was the same and consistent with the particle loss path observed in the 2D numerical simulation (Fig. 8). Different from 2D numerical simulation, 3D numerical simulation is suitable for seepage erosion phenomena in the case of local failure of the tunnel, such as bolt holes or damages at the grouting holes. Therefore, the shapes of stress archings in 2D and 3D models would be different. In the 2D model, the stress arching was part of the ellipse, intersecting near the arching waist of the tunnel, refer to Fig. 12; in the 3D model, the stress arching was part of the ellipsoid, intersecting near the arching waist and the crest of the tunnel, refer to Fig. 14.

4.3. 3D stress analysis 2D simulation can intuitively observe the process of seepage erosion, which is convenient for summarizing and discovering the law of seepage erosion. 3D simulation is more in line with the actual force and constraint state of the investigated object. Therefore, the stress analysis method applied in the 2D numerical simulation was implemented into the 3D numerical simulation to verify whether the above analysis method is reasonable or not. A 3D numerical model was built to simulate the seepage erosion of tunnel leaking, as shown in Fig. 13(a). The numerical model was 10 m long, 6 m wide and 8 m high, with a cylinder surface for simulating tunnel lining. The cylinder surface had a diameter of 2 m and was provided with a circular hole having a diameter of 250 mm for simulating tunnel leaking. Like the 2D numerical simulation method, the process of seepage erosion in this 3D numerical simulation was simulated by the fluid–solid coupling method of FDM-DEM. In addition to the DEM model, the seepage field calculated by CFD is also very important (Liu et al., 2018, 2019; Wu et al., 2019b; Xu et al., 2019). Fig. 13(b) shows the fluid force that the particles were subjected to after fluid-solid coupling, where the sizes of the arrows reflected the force magnitudes. It can be seen that the fluid force at the opening was the largest, followed by the fluid force above the tunnel, and the fluid force under the tunnel was the smallest. By applying the previously introduced method for solving 2D stressarching problem to the solution of 3D stress arching, Eq. (16) was calculated based on X-Z plane and Y-Z plane, respectively; then, the

4.4.2. Stress changes caused by seepage erosion As demonstrated previously, seepage erosion due to tunnel leaking had a great impact on the stress distribution in the ground. Fig. 16 shows the contour maps of the stress ratio R with C/D = 2.0–4.0. R is defined as below: 13

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) C/D = 1.0

(b) C/D = 1.5

(c) C/D = 2.0 Fig. 14. Vertical stresses (left column) and 3D stress archings (right column) at C/D = 1.0–2.0.

R=

σn − σ0 σ0

As shown in Fig. 16, at different C/D, the zone with the maximum variations in the soil stress featuring the maximum absolute R magnitude (the zone surrounded by red lines in Fig. 16) almost occurred directly above the tunnel. If taking a close look at Fig. 16, it can be seen that R had negative magnitudes of the vertical stress in the zone directly above the tunnel, i.e., during tunnel leaking, the vertical stress had a tendency of decreasing. For the horizontal stresses in the zone above

(18)

where σ0 represents the horizontal or vertical stress before seepage erosion; σn represents the horizontal or vertical stress after seepage erosion. Positive and negative R values represent increase and decrease in soil stress, respectively. Due to the symmetry of the model, only the stress changes of the left half model were analyzed. 14

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) 3D numerical model with a narrow longitudinal crack

(b) Y/D = 2.0

(c) Y/D = 1.5

(d) Y/D = 1.0

Fig. 15. The case of the 3D numerical model with a narrow longitudinal crack.

the flow ellipses in Section 3.2. Compared with the stress arching, the displacement arching was relatively narrower. Since the formation of the stress arching requires a certain thickness of the overburdening soils, there is a certain correlation between the stress arching and the ground subsidence. In order to explore this correlation, the distance from the vault of the stress arching to the center of the displacement at the ground level was defined as hC / D . It can be observed in Fig. 17 that h4.0 > h3.0 > h2.0 , i.e., the larger the burial depth of the tunnel, the thicker the stratum where the stress is not disturbed. This result means that if a deep buried tunnel leaked, the water-rich sand stratum would have a large stress arching enough for supporting the upper soil stress and then the ground subsidence would not be obvious. This result was consistent with the finding of Franza et al. (2019). The volume loss caused by tunnel leakage has two meanings: the ground volume loss, VG , due to leakage at the tunnel to the surface and subsurface soil volume loss, VS . VG was obtained by calculating the area of settlement trough at the ground level after tunnel leaking and VS was obtained by the sum of the volume of lost particles caused by tunnel leaking. Due to the voids and the volume of soil expanded in shear during the tunnel leakage process, their values are not equal in the granular soil. Table 3 summarized VS and VG in 2D numerical simulation at C/D = 1.5–4.0, in which Vtun refers to the volume of the tunnel per unit length. The ratio of subsurface soil volume loss to ground volume

the tunnel, R featured negative values at C/D = 2.0–3.0; however, R had positive value as C/D increased to 4.0. This inconsistency disclosed that reduction in vertical stresses was greater than the reduction in the horizontal stresses. Thereby, the lateral earth pressure coefficients in Table 2 were greater than 1.0, a clear sign of soil-arching formed above the leaking tunnel. 4.4.3. Relationship between soil movement and stress arching Soil movement is one of the consequences of tunnel leakage and critical for assessing the risk of damage to nearby structures. Moreover, soil movements can be used as an input for SSI (soil-structure dynamic interaction) two stage analyses based on continuum/Winkler models or modification factor approaches (Franzius et al., 2004; Klar et al., 2005; Liu et al., 2018a, 2018b; Franza and DeJong, 2019). In seepage erosion, the distribution of stress arching is interrelated with soil movement. On the one hand, soil movement caused stress redistribution in the strata and a stress arching was gradually formed; on the other hand, the stability of the stress arching restrained movement of the soil. In order to explore the relationship between soil movement and stress arching, Fig. 17 plots the contour lines of soil movement and stress arching at different C/D. The red, the blue, and the black lines in Fig. 17 represent the stress arching, the displacement contour, and the displacement at the ground level, respectively. Obviously, the displacement contours (displacement archings) were elliptical in shape, similar to the shape of 15

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

(a) C/D = 2.0

(b) C/D = 3.0

(c) C/D = 4.0 Fig. 16. Distributions of stress ratios.

Fig. 17. Displacement arching and stress arching in the 2D numerical simulation. Note: hC / D represents the distance from the vault of the stress arch to the center of the surface collapse.

loss was defined as RVL , i.e, RVL = VS / VG . As shown in Table 3, the value of RVL decreased as C/D increased. This can be interpreted by that when the buried depth was relatively large, the particles formed stable stress archings; thereby, the amount of lost particles was reduced. Moreover, since the particles near the tunnel opening were subjected to large overburdening pressure from the upper stratum, the volume of shear expansion caused by the loss of particles was small. Consequently, VG was reduced accordingly. At C/D = 1.5, the tunnel was too close to the ground surface and hence leaking caused less particle loss followed by smaller ground volume loss. As a result, relatively smaller RVL was

Table 3 The volume losses caused by tunnel leakage at different C/D. C/D

VS

VG

RVL

4.0 3.5 3.0 2.5 2.0 1.5

0.50Vtun 0.69Vtun 0.73Vtun 1.53Vtun 1.62Vtun 0.34Vtun

1.06Vtun 1.36Vtun 1.32Vtun 2.29Vtun 2.37Vtun 0.65Vtun

0.47 0.51 0.55 0.67 0.69 0.52

16

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

numerical models, the formation of soil arching during seepage erosion was observed. Using the principle of natural balance arch, the soil arching above the tunnel (stress arching) was derived. Stress arching was the result of changes in the soil stress arising from tunnel leakage, which was closely related to the buried depth of the tunnel. Ground subsidence above the leaking tunnel was inherently related to the distribution of the stress arching. Generally, the wider the distribution of the stress arching, the smaller the ground subsidence. (3) Unlike 2D stress arching in an ellipse shape, 3D stress arching featured a shape in ellipsoid and can better capture the trend of 3D vertical stress distribution after seepage erosion. Since it was inconvenient to compare the distribution of the contact force chain in a 3D model with the 3D stress arching, no discussion on the force chains in the 3D model was conducted in this study. The increase of the 3D stress arching with the buried depth of the tunnel indicated that the stress arching was closely related to ground subsidence directly above the leaking tunnel. (4) The expression of the stress arching due to tunnel leaking in waterrich sand stratum was built upon the shape of the natural balance arch under the loose medium proposed by Miao (1990). Combined with the force chain distributions observed in the numerical model, the coordinate origin of the expression of the natural balance arch was relocated above the tunnel for description of the stress-arching due to tunnel leaking. By assuming that the natural balance arch intersected the tunnel, the arching above the tunnel was taken as the stress arching and the calculated stress arching was close to the outermost archings of force chains observed in the numerical model. (5) This is a new attempt to apply the soil arching effect to the study of the seepage erosion around a leaking tunnel buried in water-rich sand. With the distribution of stress arching, it is possible to judge the stability of the ground above a leaking tunnel and the severity of ground subsidence.

(a) Water-rich sand

(b) Dry sand Fig. 18. Relationship between surface subsidence and particle loss.

obtained. 4.4.4. Influence of groundwater on surface subsidence As shown in Fig. 17, the ground had the maximum ground subsidence at C/D = 2.0. Hence, the ground settlements at C/D = 2.0 for two different cases (particle losses through the opening on the tunnel in case of water-rich sand and dry sand) were selected for analyses. Fig. 18 recorded the process in which the surface gradually sank with the loss of particles in the two cases. The amount of particle loss was defined as the ratio of the mass of the particle loss to the total particle mass above the tunnel (RPL ). It is apparent that under the action of groundwater, the amount of particle loss was greatly increased, which resulted in much larger ground settlement than that of dry sand. In addition, the force of the groundwater acting on the particles played a key role in the loss path of the particles, which caused the particles to move toward the opening. Therefore, the magnitude of the ground settlement directly above the tunnel center in the case of ground water (Fig. 18a) was almost three times that in the case of dry sand (Fig. 18b). Hence, it can be concluded that the presence of groundwater would incur more severe damage.

Declaration of Competing Interest The authors (Ying-ying Long; Yong Tan) acknowledge that there is no conflict of interest for this paper (Title: Soil Arching due to Leaking of Tunnel Buried in Water-Rich Sand; MS. No. TUST-2019-744) Acknowledgement Financial supports from the National Natural Science Foundation of China (No. 41877286) and the National Basic Research Program (973 Program) (Grant 2015CB057800) are gratefully acknowledged. The comments and suggestions from the three anonymous reviewers and Editor-in-Chief Prof. Qianbing Zhang are sincerely appreciated.

5. Summary and conclusions

Appendix A. Supplementary material

Through the experimental test and both 2D and 3D FDM-DEM coupling numerical simulations, the seepage erosion due to tunnel leaking in water-rich sand stratum was investigated and special attention was paid to the soil stress variation and ground displacement. In general, the following major conclusions can be drawn from this study:

Supplementary data to this article can be found online at https:// doi.org/10.1016/j.tust.2019.103158. References

(1) The experimental test show that with the development of seepage erosion due to tunnel leaking, the ground settled gradually. At the end of the test, seepage erosion of the soil particles did not completely destroy the stratum above the tunnel, an evident sign of soil arching effect above the tunnel. The consistency between the 2D numerical simulation results and the experimental results demonstrated that FDM-DEM coupling is an effective and reliable method for simulating seepage erosion due to leaking of tunnel buried in water-rich sand stratum. (2) By analysis of soil stress and contact force chains in the 2D

Adachi, T., Kimura, M., Kishida, K., 2003. Experimental study on the distribution of earth pressure and surface settlement through three-dimensional trapdoor tests. Tunn. Undergr. Space Technol. 18 (2–3), 171–183. Brady, B.H.G., Brown, E.T., 2004. Rock Mechanics for Underground Mining. Allen and Unwin. Cai, M., Kaiser, P.K., Morioka, H., Minami, M., Maejima, T., Tasaka, Y., Kurose, H., 2007. FLAC/PFC coupled numerical simulation of AE in large-scale underground excavations. Int. J. Rock Mech. Min. Sci. 44 (4), 550–564. Carman, P.C., 1956. Flow of Gases through Porous Media. Butterworths Scientific Publications. Chai, J., Shen, J.S., Yuan, D.J., 2018. Mechanism of tunneling-induced cave-in of a busy road in Fukuoka city, Japan. Undergr. Space 3 (2), 140–149.

17

Tunnelling and Underground Space Technology 95 (2020) 103158

Y.-Y. Long and Y. Tan

Lu, Y., Li, X., Wang, Y., 2018a. Application of a flexible membrane to DEM modeling of axisymmetric triaxial compression tests on sands. Eur. J. Environ. Civ. Eng. 22 (S1), s19–s36. Lu, Y., Sun, H.Q., Li, H., 2018b. Application of modified DIC technology in press-in pile model test. J. Shanghai Univ. 24 (06), 993–1001 (in Chinese). Lu, Y., Tan, Y., Li, X., 2018c. Stability analyses on slopes of clay-rock mixtures using discrete element method. Eng. Geol. 244, 116–124. Lu, Y., Tan, Y., Li, X., Liu, C., 2017. Methodology for simulation of irregularly shaped gravel grains and its application to DEM modeling. J. Comput. Civil Eng. 31 (5), 04017023. Lyu, H.M., Shen, S.L., Zhou, A.N., Yang, J., 2019. Perspectives for flood risk assessment and management for mega-city metro system. Tunn. Undergr. Space Technol. 84, 31–44. Ma, S.L., Wei, Z.Y., Chen, X.L., 2018. A systematic approach for numerical research of realistic shaped particle-fluid interactions. Powder Technol. 339, 377–395. Miao, X.X., 1990. Stability of natural balance arch and roadway surrounding rock. J. Min. Saf. Eng. 1990 (2), 55–57. Minh, N.H., Cheng, Y.P., 2013. A DEM investigation of the effect of particle-size distribution on one-dimensional compression. Géotechnique 63(1), 44–53. Morsi, S.A., Alexander, A.J., 1972. An investigation of particle trajectories in two-phase flow systems. J. Fluid Mech. 55 (02), 193–208. Mukunoki, T., Kumano, N., Otani, J., Kuwano, R., 2009. Visualization of three dimensional failure in sand due to water inflow and soil drainage from defective underground pipe using X-ray CT. Soils Found. 49 (6), 959–968. National Institute for Land and Infrastructure Management (NILIM), 2006. Annual Report of Waste Water Management and Water Quality Control. Ministry of Land, Infrastructure and Transport, Japan. Ni, J.C., Cheng, W.C., 2012. Characterising the failure pattern of a station box of Taipei Rapid Transit System (TRTS) and its rehabilitation. Tunn. Undergr. Space Technol. 32, 260–272. Ouyang, M., Takahashi, A., 2016. Influence of initial fines content on fabric of soils subjected to internal erosion. Can. Geotech. J. 53 (2), 299–313. Protodyakonov, M.M., 1907. Rock pressure on mine support (theory of mine support). Tipografiya Gubernskogo Zemstva 23–45. Rui, R., Han, J., van Eekelen, S.J.M., Wang, Y., 2019. Experimental investigation of soilarching development in unreinforced and geosynthetic-reinforced pile-supported embankments. J. Geotech. Geoenviron. Eng. 145 (1). https://doi.org/10.1061/ (ASCE)GT.1943-5606.0002000. Saiang, D., 2010. Stability analysis of the blast-induced damage zone by continuum and coupled continuum-discontinuum methods. Eng. Geol. 116 (1–2), 1–11. Sato, M., Kuwano, R., 2015. Influence of location of subsurface structures on development of underground cavities induced by internal erosion. Soils Found. 55 (4), 829–840. Sterpi, D., 2003. Effects of the erosion and transport of fine particles due to seepage flow. Int. J. Geomech. 3 (1), 111–122. Tan, Y., Lu, Y., 2017. Forensic diagnosis of a leaking accident during excavation. J. Perform. Constr. Facil 31 (5), 04017061. Tan, Y., Lu, Y., Xu, C., Wang, D., 2018. Investigation on performance of a large circular pit-in-pit excavation in clay-gravel-cobble mixed strata. Tunn. Undergr. Space Technol. 79, 356–374. Tang, Y., Chan D.H., Zhu, D.Z., 2017. A coupled discrete element model for the simulation of soil and water flow through an orifice. Int. J. Numer. Anal. Methods Geomech. 41(14), 1477–1493. Terzaghi, K., 1936. Stress distribution in dry and in saturated sand above a yielding trapdoor. Proceedings of the 1st International Conference on Soil Mechanics, Harvard University, Cambridge, Massachusetts. Terzaghi, K., 1943. Theoretical Soil Mechanics. John Wiley and Sons, New York. Vivanco, F., Rica, S., Melo, F., 2012. Dynamical arching in a two dimensional granular flow. Granular Matter 14 (5), 563–576. Wu, J., Liao, S.M., Liu, M.B., 2019a. An analytical solution for the arching effect induced by ground loss of tunneling in sand. Tunn. Undergr. Space Technol. 83, 175–186. Wu, H.N., Shen, S.L., Yang, J., 2017a. Identification of tunnel settlement caused by land subsidence in soft deposit of Shanghai. J. Perform. Constr. Facil 31 (6), 04017092. Wu, J.J., Wang, L.J., Cheng, Q.G., 2017b. Soil arching effect of lattice-shaped diaphragm wall as bridge foundation. Front. Struct. Civil Eng. 11 (4), 446–454. Wu, Y.X., Lyu, H.M., Han, J., Shen, S.L., 2019b. Dewatering-induced building settlement around a deep excavation in soft deposit in Tianjin, China. J. Geotech. Geoenviron. Eng., ASCE 145 (5), 05019003. Xu, Y.S., Yan, X.X., Shen, S.L., Zhou, A.N., 2019. Experimental investigation on the blocking of groundwater seepage from a waterproof curtain during pumped dewatering in an excavation. Hydrogeol. J published online. Yao, T., Chan, D.H., Zhu, D.Z., Guo, S., 2018. An analytical solution for steady seepage into a defective pipe. Water Sci. Technol.-Water Supply 18 (3), 926–935. Zhang, D.M., Gao, C.P., Yin, Z.Y., 2019. CFD-DEM modeling of seepage erosion around shield tunnels. Tunn. Undergr. Space Technol. 83, 60–72. Zhu, H.P., Zhou, Z.Y., Yang, R.Y., Yu, A.B., 2007. Discrete particle simulation of particulate systems: theoretical developments. Chem. Eng. Sci. 62 (13), 3378–3396.

Chen, R.P., Li, J., Kong, L.G., Tang, L.J., 2013. Experimental study on face instability of shield tunnel in sand. Tunn. Undergr. Space Technol. 33, 12–21. Chevalier, B., Combe, G., Villard, P., 2012. Experimental and discrete element modeling studies of the trapdoor problem: influence of the macro-mechanical frictional parameters. Acta Geotech. 7 (1), 15–39. Cundall, P.A., 1971. A computer model for simulating progressive large-scale movements in block rock mechanics. In: Proceedings of the Symposium of International Society of Rock Mechanics, Nancy, II, Article 8. Cundall, P.A., Strack, O.D.L., 1979. Discrete numerical-model for granular assemblies. Géotechnique 29 (1), 47–65. Dai, B.B., Yang, J., Luo, X.D., 2014. A numerical analysis of the shear behavior of granular soil with fines. Particuology 21, 160–172. Davies, J.P., Clarke, B.A., Whiter, J.T., Cunningham, R.J., 2001. Factors influencing the structural deterioration and collapse of rigid sewer pipes. Urban Water 3 (1–2), 73–89. Franza, A., DeJong, M.J., 2019. Elastoplastic solutions to predict tunneling-induced load redistribution and deformation of surface structures. J. Geotech. Geoenviron. Eng. 145 (4), 04019907. Franza, A., Marshall, A.M., Zhou, B., 2019. Greenfield tunnelling in sands: the effects of soil density and relative depth. Géotechnique 69, 297–307. Franzius, J.N., Potts, D.M., Addenbrooke, T.I., Burland, J.B., 2004. The influence of building weight on tunnelling-induced ground and building deformation. Soils Found. 44 (1), 25–38. Fullard, L.A., Davies, C.E., Neather, A.C., Breard, E.C.P., Godfrey, A.J.R., Lube, G., 2018. Testing steady and transient velocity scalings in a silo. Adv. Powder Technol. 29 (2), 310–318. Grudzien, K., Chaniecki, Z., Babout, L., 2017. Study of granular flow in silo based on electrical capacitance tomography and optical imaging. Flow Meas. Instrum. 62, 186–195. Guo, S., Shao, Y., Zhang, T.Q., Zhu, D.Z., Zhang, Y.P., 2013. Physical modeling on sand erosion around defective sewer pipes under the influence of groundwater. J. Hydraul. Eng. 139 (12), 1247–1257. Hewlett, W.J., Randolph, M.F., 1988. Analysis of piled embankments. Ground Eng. 21 (3), 12–18. Itasca Consultant Group, 2012. Fast Lagrangian Analysis of Continua in 3 Dimensions Online Manual Version 5.0. Consulting Group Incorporated, Minneapolis, Minnesota. Itasca. Itasca Consultant Group, 2015. Particle Flow Code in 2 Dimensions Online Manual Version 5.0. Consulting Group Incorporated, Minneapolis, Minnesota. Itasca. Janelid, I., Kvapil, R., 1966. Sublevel caving. Int. J. Rock Mech. Min. Sci. 3 (2), 129–153. Jiang, M.J., Yu, H.S., Harris, D., 2006. Discrete element modelling of deep penetration in granular soils. Int. J. Numer. Anal. Meth. Geomech. 30 (4), 335–361. Ke, L., Takahashi, A., 2012. Strength reduction of cohesionless soil due to internal erosion induced by one-dimensional upward seepage flow. Soils Found. 52 (4), 698–711. Ke, L., Takahashi, A., 2014. Experimental investigations on suffusion characteristics and its mechanical consequences on saturated cohesionless soil. Soils Found. 54 (4), 713–730. Kenny, T.C., 1985. Internal stability of granular filters. Can. Geotech. J. 22 (2), 215–225. Klar, A., Vorster, T.E.B., Soga, K., Mair, R.J., 2005. Soil-pipe interaction due to tunneling: comparison between Winkler and elastic continuum solutions. Géotechnique 55 (6), 461–466. Kong, X.X., Liu, Q.S., Zhang, Q.B., Wu, Y.X., Zhao, J., 2018. A method to estimate the pressure arch formation above underground excavation in rock mass. Tunn. Undergr. Space Technol. 71, 382–390. Korff, M., Mair, R.J., Van Tol, A.F., Kaalberg, F., 2011. Building damage and repair due to leakage in a deep excavation. Proc. Inst. Civil Eng-Forensic Eng. 164(4), 165–177. Kozeny, J., 1927. Überkapillareleitung des wassersimboden: (aufstieg, versickerungundanwendung auf die bewässerung). Hölder-Pichler-Tempsky. Lai, H.J., Zheng, J.J., Zhang, J., Zhang, R.J., Cui, L., 2014. DEM analysis of “soil”-arching within geogrid-reinforced and unreinforced pile-supported embankments. Comput. Geotech. 61 (9), 13–23. Lee, C.J., Wu, B.R., Chen, H.T., Chiang, K.H., 2006. Tunnel stability and arching effects during tunneling in soft clayey soil. Tunn. Undergr. Space Technol. 21 (2), 119–132. Liu, X., Fang, Q., Zhang, D., 2018a. Mechanical responses of existing tunnel due to new tunnelling below without clearance. Tunn. Undergr. Space Technol. 80, 44–52. Liu, X.X., Shen, S.L., Xu, Y.S., Yin, Z.Y., 2018b. Analytical approach for time-dependent groundwater inflow into shield tunnel face in confined aquifer. Int. J. Numer. Anal. Meth. Geomech. 42, 655–673. Liu, X.X., Shen, S.L., Zhou, A., Xu, Y.S., 2019. Evaluation of foam conditioning effect on groundwater inflow at tunnel cutting face. Int. J. Numer. Anal. Meth. Geomech. 43, 463–481. Lopez, R.D., Silfwerbrand, J., Jelagin, D., Birgisson, B., 2016. Force transmission and soil fabric of binary granular mixtures. Géotechnique 66 (7), 578–583. Lu, Y., Jiang, Y., Wang, X.Y. 2020. Investigation on soil displacements caused by press-in process of close-ended model piles using imaging technique. Acta Geotechnica Slovenia, in press.

18