Undrained stability of active and passive trapdoors

C

soil cover

D

tunnel diameter

E

extent of the failure surface

H

depth from the ground surface to the trapdoor opening

N

‘designed’ stability number

N_c

critical stability number

r²

correlation coefficient

S_u

undrained shear strength of soil

W

trapdoor width

γ

soil unit weight

σ_s

surcharge pressure

σ_t

support pressure

Φ

internal friction angle of soil

Sinkholes present environmental risk through subsidence or sudden ground collapse, leading to loss of life and infrastructure. The recent growth in the number of sinkhole occurrences due to human activities, such as urbanisation, mining and agricultural development, has highlighted the need for better understanding and prediction of the problem (Drumm et al., 2009).

Sowers (1996) outlined the sub-profile of karst soil and described the process of forming a sinkhole. It was suggested that in limestone areas, the gradual erosion of rock at a depth caused by the passing of underground water leads to subsidence of overburden deposited soil, resulting in a saucer-shaped depression. Field investigation studies (Newton, 1976; Sowers, 1996) also suggested that the underground voids, created either naturally or by humans, initiated in cracks between the underground rocks. As indicated by Tharp (2003), the initial size of the cavity does not reflect the actual size of the trapdoor at collapse because the size of the initial cavity will grow further due to internal erosion and will create a reverse-funnel shape.

A considerable number of studies have been published on the stability of trapdoors. A study was initiated by Terzaghi (1936), who experimentally investigated the effect of distributed stress in sand. The study categorised the failure as either an active or a passive mode and described active failure as occuring due to overburden pressure and passive failure occurring as an uplifting force such as an anchor.

Through laboratory experiments and field data collection, Broms and Bennermark (1967) stated that the support pressure required to maintain the stability of an opening on a vertical wall should equate to overburden pressures (surcharge and self-weight) and the undrained shear strength of the soil multiplied by a ‘factor’. The stability number (N) was therefore defined in the following equation

where σ_s is the surface surcharge pressure; σ_t is the support pressure; H is the depth of the opening; γ represents the soil unit weight; and S_u is the undrained shear strength of the soil. Broms and Bennermark (1967) also concluded that the opening would be unstable when the overburden pressure (σ_s + γH) is six times greater than the undrained shear strength of the soil.

One of the most influential studies of underground stability comes from Davis et al. (1980), who used an analytical approach to study tunnel heading stability. Davis et al. (1980) used the limit analysis theorem to determine upper-bound and lower-bound solutions to the problem. In their approach, unlike the original Broms and Bennermark (1967) stability number approach, the problem was approached differently by using a critical pressure ratio (σ_s − σ_t)/S_u, which is a function of the strength ratio (γD/S_u) and depth ratio (H/W), as indicated in Equation 2. Numerous studies have since been performed using this approach (Augarde et al., 2003; Drumm et al., 2009).

Koutsabeloulis and Griffiths (1989) used the displacement finite-element method to investigate soil displacement in active and passive modes of the trapdoor problem. Martin (2009) introduced a new slip-line solution for a shallow trapdoor. Craig (1990) utilised a centrifuge model to investigate the critical stability of a circular cavity. Using the centrifugal approach, many experimental investigations of trapdoor stability have been carried out by researchers, such as Abdulla and Goodings (1996) and Jacobsz (2016). Recently, Keawsawasvong and Ukritchon (2017) studied active trapdoor problems with a linear increase in undrained shear strength with depth using finite-element limit analysis (FELA).

Although extensive research has been carried out on the stability of trapdoors in the past, most of the studies have predominately focused on the use of the critical pressure ratio ((σ_s − σ_t)/S_u). Very few studies have used the original stability number approach (Broms and Bennermark, 1967) to study soil stability and further explore the relationship between the stability number and the factor of safety (FOS). In this paper, the shear-strength-reduction method (SSRM) is used with the finite-difference method (FDM) to obtain the FOS for a wide range of stability numbers (N) and depth ratios (H/W). This study also investigates the extent of sinkhole collapse on the ground surface. The results are validated by using FELA and other previous published studies. These FOS results are used to produce comprehensive design charts for the problem of trapdoor stability.

The development of cover-collapse sinkholes is a complex procedure due to the continuous expansion of the cavity size over the time. To simplify the problem, it is assumed that the cavity is in the critical stage where the failure is imminent. Figure 1 shows the problem definition of an idealised horizontal trapdoor underlying a homogeneous layer of cohesive soil. The undrained soil is modelled as uniform Mohr–Coulomb material with a zero soil internal friction angle (Φ = 0). S_u is the undrained shear strength, and γ is the soil unit weight. The trapdoor opening width is W, and the depth from the surface to the trapdoor opening is notated as H.

Note that the combination of surcharge pressure (σ_s), overburden pressure (γH) and support pressure (σ_t) can produce failure in either a collapse or a blowout. For undrained clay without volume loss during plastic shearing, stability results are independent of loading directions, and Broms and Bennermark’s original undrained stability number is a suitable design parameter (Shiau and Al-Asadi, 2018). For drained soils, the original stability number is not applicable.

A broad range of stability numbers (N = −15–15) and depth ratios (H/W = 1–10) have been chosen to cover all possible investigations of collapse and blowout. Note that the actual values of σ_s, σ_t, γ, S_u, H, and W used in the analyses are insignificant and are not to be reported here due to the nature of the dimensionless definition. The SSRM is adopted to solve the FOS, which is a function of the stability number N and the depth ratio H/W, as shown in the following equation

Due to the lack of previous literature on implementing the FOS approach as well as correlating the stability number N to the FOS, it is important to use two techniques for the comparison of results. For this reason, both the FDM and FELA were used to analyse the problem in this paper.

The FDM is one of the oldest techniques used in numerical studies and is a powerful method for analysing complex geotechnical stability problems involving non-linear solutions (Itasca, 2003). A typical grid for simulating the trapdoor stability problem in the FDM is shown in Figure 2. To improve the computational efficiency, a symmetrical condition is considered. This symmetrical condition is particularly important for deep cases that normally require more central processing unit time. An effective domain should be such that it is large enough to present the entire velocity field. Both the left and the right boundaries of the mesh were fixed in the x-direction, allowing the soil to move in a vertical direction. The lower boundary of the mesh was restrained in both the x- and y-directions except where the trapdoor opening is positioned. The trapdoor opening was not restrained so the soil body can freely move downwards into the cavity. A Fish script was developed to assist in auto mesh generation and problem solving of the trapdoor problem using SSRM (Shiau and Sams, 2019; Shiau et al., 2018). The Fish development was a particularly important tool in this study, as it allows parametric study to be conducted efficiently.

The numerical process of the FELA and SSRM is based on the limit theorems of classical plasticity, which were described by Lyamin and Sloan (2002a, 2002b) and Sloan (2013). The details of the formulation will not be repeated here. It is worth noting that the new technique utilises the finite-element discretisation to deal with complicated geometry and loading conditions and the plastic bounding theorems to bracket the true limit load using upper- and lower-bound solutions. Figure 3 shows the mesh used in the paper. The domain sizes of the models were carefully chosen by observing the non-zero velocity fields, so as to minimise boundary effects. Note that the symmetrical faces are fixed only in the normal direction (i.e. x-direction) to allow vertical movement, so as the outer face boundary. Similar to the FDM mesh, the lower face is restrained in the x- and y-directions and the top face is free to displace in all directions. The solutions are then triggered by involving the shear-strength-reduction technique stated in the paper by Krabbenhoft and Lyamin (2015).

The SSRM was utilised as early as 1975 by Zienkiewicz et al. (1975) and many other researchers to investigate various geotechnical engineering problems (Griffiths and Lane, 1999; Krabbenhoft and Lyamin, 2015; Matsui and San, 1992; Michalowski, 2002; Ugai and Leshchinsky, 1995; Yang and Drumm, 2002). Although the SSRM provides a straightforward solution to many geotechnical problems, such as slopes and retaining walls, this method has seldom been used in the analysis of underground stability problems (Shiau et al., 2017, 2018). In this study, the SSRM and FOS approach is adopted to analyse the stability of the trapdoor in collapse and blowout conditions. A total of 690 trapdoor cases are studied using the FDM, upper-bound FELA (FELA UB) and lower-bound FELA (FELA LB). Numerical results of the extensive investigation are presented in the form of design charts and equations.

Comprehensive findings of this study are presented in Tables 1–3 for a broad range of stability numbers (N = −15–15) and depth ratios (H/W = 1–10). Using the data in the tables, Figure 4 plots the FOS results of collapse and blowout from the analyses of lower bound, upper bound and finite differences for the depth ratio of H/W = 3. The results show that the curves are in hyperbolic form where FOS and N are the vertical and horizontal asymptotes, respectively. The general equation of the curve is presented as follows

Equation 4 suggests that for a given depth ratio (H/W = 3), any combination of FOS and N on the curve yields a unique value. This unique value is the critical stability number (N_c), which corresponds to an FOS of 1. By drawing an FOS = 1 horizontal line in Figure 4, the two intersection points give an N_C value of 4·925 for the collapse and −4·930 for the blowout.

When the supporting pressure ratio (SPR, σ_t/S_u) is greater than the overburden pressure ratio (OPR ((σ_s + γH)/S_u)), the negative value of N represents a blowout movement. Contrary to this, a positive value of N indicates that the soil moves in the collapse condition. This occurs when the overburden pressure ratio ((σ_s + γH)/S_u)) is greater than the supporting pressure ratio (σ_t/S_u). As N further increases, an incipient collapse is reached where FOS = 1 and the corresponding N is the critical N_c. When the supporting pressure ratio (σ_t/S_u) is equal to the overburden pressure ratio ((σ_s + γH)/S_u)), N is equal to zero and FOS is at a maximum (infinite) where a ‘stressless’ scenario exists on the asymptote line.

Broms and Bennermark’s original Equation 1 can be rearranged into a form that is more amenable to analysis, as shown in the following equation

Using Equation 5, a critical supporting pressure σ_t (when FOS = 1) can be determined as long as N_c is known. Note that there is only one unique N_c (±) for a particular depth ratio and N_c is a function of the depth ratio H/W regardless of the undrained shear strength of the soil. It is therefore important to study the effect of H/W on the critical stability number N_c. Figure 5 shows such a relationship between N_c and H/W. Note that the critical stability number (N_c) increases non-linearly as H/W increases, and the gradient of the curve decreases for large values of N_c. The area bounded by the collapse and the blowout curves represents the safe zone where FOS > 1. As the stability number (N) approaches zero (OPR = SPR), the FOS becomes infinite. Also, see the asymptote in Figure 4.

Note that the critical stability number (N_c) represents the ‘design’ N value when FOS = 1. This is the failure envelope for FOS = 1. Also, note that FOS ˃ 1 when a design N value is located within the failure envelope and FOS < 1 when a design N value is located outside the failure envelope.

The finite-difference results of N_c were chosen for the regression analysis. These are presented in Equations 6 and 7 for collapse and blowout, respectively, with a correlation coefficient r² = 0·998.

Using Equations 1 and 4, the FOS can be determined with the following equation

Substituting Equations 6 and 7 into Equation 8, Equation 9 can be used to determine the FOS for known design parameters (σ_s, σ_t, γ, H, W and S_u).

Equation 9 is also presented graphically in Figure 6. The design contour map of FOS was constructed based on the FDM numerical solutions.

By rearranging Equation 8, one can determine the required support pressure σ_t for a given FOS using the following equation

Figure 7 and Table 4 compare the N _c values obtained in this paper with those in the published literature. The comparison shows that the FDM results of N _c are considerably lower than the analytical upper-bound solutions of Davis (1968) for large H/W. Although Davis’s investigation had been improved by Gunn (1980) by using the three rigid block parameters, the analytical upper-bound solutions of Gunn (1980) are still 4·0–12·9% larger than the current FDM.

The comparison in Figure 7 also shows some 0·65–15·00% difference between the FDM results and the results by Sloan et al. (1990) for the range of depth ratios (H/W = 1–10) studied. Although Sloan et al. (1990) carried out extensive research on trapdoor stability, the use of linear programming was a major drawback on the solution accuracy. In particular, the lower bound shows a large variation in comparison with FDM. The slip-line solutions of Martin (2009) were limited to shallow depth ratios up to H/W = 2. Martin’s results agree well with the current FDM solutions of H/W = 1 and 2. Keawsawasvong and Ukritchon (2017) presented averaged upper- and lower-bound solutions for a wide range of depth ratios (H/W) for homogeneous and non-homogeneous clays. In general, their results of homogeneous clay agree well with the FDM results. Overall, the predicted trend of FDM results shows a good agreement with the FELA LB and UB solutions.

Results of the failure extent investigation are shown in Figure 8. The distance of failure extent was determined by inspection of the velocity vector plots produced in the program. Figure 8 suggests that the failure extent ratio (E/W) is linearly proportional to the depth ratio (H/W). The linear relationship is presented in the following equation

A practical conclusion can be drawn from Table 5, showing that an approximately 55° line can be drawn from the centre of the trapdoor to the outer boundary of the failure surface to estimate the failure extent. Note that this investigation is valid for all values of stability numbers (N).

An old vertical mining shaft has no internal pressure and no surcharge pressure. For the given parameters (S_u = 154 kPa, γ = 18 kN/m³, H = 36 m and W = 6 m), determine the FOS.

• Since there is no internal pressure, only the collapse failure should be considered.

• The stability number is N = γH/S_u = 4·21.

• Using H/W = 6 and N = 4·21, Equation 9 gives an FOS of 1·58 for the collapse.

• Using H/W = 6 and N = 4·21, Figure 6 gives an approximate FOS of 1·6.

• An actual computer analysis of this case gives an FOS of 1·63.

What is the critical support pressure σ_t when FOS = 1?

• Using Equation 6, N_c = 3·45 for collapse.

• From Equation 10, σ_t = σ_s + γH − (N_c × S_u/FOS) = 100 + (18 × 10) − (3·45 × 30/1) = 176·5 kPa.

An existing sinkhole has a diameter of 10 m. Estimate the depth of the sinkhole using the following parameters for cohesive soil: S_u = 54 kPa and γ = 19 kN/m³.

• Note that Equation 11 is independent of the design parameter N. The only needed information is E.

• From Equation 11, 10/W = 1·39 × (H/W) + 0·13.

• By ignoring the small value of (0·13W), the depth (H) is found to be 7·19 m.

An FOS of 4 is required for the design of an underground military bunker where the surcharge pressure is given as σ_s = 50 kPa. The following parameters are known: S_u = 25 kPa, γ = 18 kN/m³, H = 40 m and W = 30 m.

• Using Equation 6, the critical stability number is N_c = 2·89. Note that Figure 5 can also be used to find the N_c value.

• Substitute the N_c value into Equation 10, σ_t = σ_s + γH − (N_c × S_u/FOS) = 50 + (18 × 40) − (2·9 × 25/4) = 752 kPa.

• The required pressure to support the cavity for an FOS of 4 is 752 kPa.

This study successfully investigated the stability of trapdoor problems using Broms and Bennermark’s original stability number. Numerical results of this study were obtained by utilising the SSRM. Three numerical techniques were used – namely, the FDM, FELA UB and FELA LB.

The FOS was found to be a function of the depth ratio (H/W) and stability number (N). The numerical results suggest that the FOS increases when the stability number (N) decreases. The FOS becomes very large when the stability N is very small. Further investigation on failure mechanisms indicates a linear relationship between the failure extent ratio (E/W) and the depth ratio (H/W). The failure angle (θ) measured from the centre of the opening (W) to the outer boundary of the failure surface is approximately equal to 55° for all depth ratios (H/W).

This investigation has improved the understanding of trapdoor stability and associated surface failure extent. Further study is needed for a more realistic three-dimensional analysis of sinkhole failure.