Sinkhole Stability in Elliptical Cavity under Collapse and Blowout Conditions

Road subsidence and sinkhole failures due to shallow cavities formed by defective water main have increased in recent decades and become one of the important research topics in geotechnical engineering. The present paper numerically studies the stability and its associated failure mechanism of ellipse-shaped cavity above defective water mains using the finite element limit analysis technique. For a wide range of geometrical parameters, the pressure ratio method is used to formulate the stability solutions in both blowout and collapse scenarios. Even though there is no published solution for elliptical cavities under blowout failure conditions, the obtained numerical results are compared with available circular solutions. Several conclusions are drawn based on the failure mechanism study of the various ellipse shape transformations in this study, whilst design charts and equations proposed for practical uses.


Introduction
Natural and human-induced sinkholes are the two main categories of sinkholes that can be instantaneous and catastrophic. The human-induced sinkhole events are speedily increasing in recent decades, and they are the news headline of most national and international media [1]. Sinkholes are a significant hazard in many areas worldwide and are studied scientifically by Waltham et al. [2], Gutierrez et al. [3] and Parise [4]. The recent increase in sinkhole events due to water main defects necessitates the problem to be better understood, predicted and prevented. The water main damage-driven sinkhole is a significant concern in many metropolises across the globe. Several factors such as maintenance lack, water pressure increased, differential settlement, root damage and corrosion are responsible for causing the water main damage [1]. Table 1 presents some of the latest sinkhole events that occur across the globe [5][6][7][8].
The trapdoor problem is a classical geotechnical stability example for the pipe leakage-driven cavities. Collapse failure (downward movement) and blowout failure (upward movement) are the two main modes of trapdoor failures. The self-weight of soil mass and ground surcharge pressures are two critical components associated to collapse failures, whereas blowout failures are only due to external forces, such as the water main pressure, exerted against the soil weight and ground surcharge. There have been numerous studies concerning "collapse" and "blowout" stability, though much of them are limited to a flat planar trapdoor [9].
To the authors' best knowledge, very few studies were reported on the stability effects of elliptical cavity shape transformation, especially for the "blowout" stability problems. Dutta and Bhattacharya [24] studied the stability of dual elliptical tunnels in soft clay using lower bound finite element limit analysis with second-order conic programming. Yang et al. [25] investigated the effects of surcharge loading on an elliptical tunnel using the upper-bound finite element method with rigid translatory moving elements in cohesive-frictional soils. Zhang et al. [26] analyzed the stability of elliptical tunnels in cohesionless soils using the upper bound solution approach. Moreover, the stability analysis of unlined elliptical tunnel using finite element upper-bound method with plastic deformation element method was studied [27,28]. The presence of a soil cavity due to subsurface soil erosion can cause a sinkhole. Interestingly, most sinkhole formations are circular at the ground surface [2][3][4][20][21][22][23]29,30]. Figure 1 presents the potential stages of water main related soil failures. The passive (blowout) and active (collapse) failures are based on the limit loads acting on the sinkhole's upward and downward orientations, respectively. Despite the uncertainty of the shape transformation, it has been found that most of the previous research were centered on the fundamental trapdoor shape transformation, i.e., circular, square and rectangle [31]. Very few studies were performed on elliptical shape transformation. Furthermore, there has also been a lack of research in relation to sinkhole stability above the pressurized drinking water mains, although sinkhole events are often seen on main roads worldwide. This paper applies advanced finite element limit analysis (FELA) and adaptive mesh techniques to study sinkhole stability due to defective water mains. It is noted that low pressure in the pipe leads to dragging the loose soil particle into the sewer pipe through the crack in the pipe, creating a cavity in the soil. In contrast, high water pressure is accountable for blowout failure. It is therefore hypothesized that the initiation of a cavity is caused by a fracture in the pipe, which leads to ultimate failure in one of two scenarios: blowout or collapse. Elliptical cavity shapes transformation is studied for a wide range of depth ratios, ellipse width to height ratios and soil shear strength ratios. Design charts and equations are presented for practical uses in the preliminary stability assessment under both collapse and blowout scenarios.

Problem Definition and FELA Model
An idealized cavity shape is assumed following Equation (1), which describes the shape of an ellipse.
where a = B/2, b = D/2 and x and y are coordinates of an ellipse. Using Equation (1), a typical adaptive FELA mesh is shown in Figure 2 for a width to height ratio B/D = 2, whilst in Figure 3 for width to height ratio B/D = 0.5. In both cases, the trapdoors have a cover (C), depth (D) and width (B). The inner perimeter of the opening is exposed to a pressure acting perpendicular throughout the inner face (σt), whereas the ground surface is exposed to a surface pressure (σs). The soil is considered as a rigid-perfectly plastic Tresca material with the soil unit weight γ, and the undrained shear strength of the soil Su. The domain size was chosen carefully to diminish boundary effects of all sides in view of the overall development of velocity field. Note that both left-hand and right-hand sides are kept stationary in the x-direction and the bottom of the domain is fixed in both x and y directions. The nodes on the ground surface are free to move in all directions.
In contrast with the Finite Element Analysis (FEA) or the Finite Different Analysis (FDA), the FELA technique employed rigorous upper and lower bound theorems, where the true solutions can always be bracketed by the two solutions. In the current investigation, both upper and lower bound theorems with finite element technique are used [32]. With the latest advances, the automatically adaptive mesh refinement is utilized in both UB and LB simulations to diminish the solution's discretization error. The upper bound (UB) limit formulation produces kinematically admissible velocity fields, whereas the lower bound (LB) formulation is based on statically permissible stress fields. Note that five iterations of adaptive meshing were employed in this study, with the number of elements increasing from 5000 to 10,000 elements.
Definition of stability number (N) and the application of it was initiated by Broms and Bennermark [33] using a vertical trapdoor problem. The effect of surface surcharge (σs), soil self-weight (γH) and the supporting pressure (σt) are combined into a single dimensionless stability number (N), and it is shown in Equation (2).
Instead of using the combined features in Equation (2), David et al. [34] proposed a new way to represent the results using a pressure ratio {PR = (σs -σt)/Su}, which is a function of soil strength ratio (SR = γD/Su), width to height ratio (WR = B/D) and depth ratio (DR = C/D) in our study. This is shown in equation (3).
Following Davis's approach, the current paper formulates numerical solutions using the critical pressure ratio for wide ranges of dimensionless parameters. A range of depth ratios (C/D = 1-5), width to height ratio (B/D= 0-2) and shear strength ratio (γD/Su = 0-3) are studied for both blowout and collapse scenarios. The objective function is to determine the lower and upper bound limits of the critical supporting pressure (σt) with selected input parameters (C, D, B, σs, γ and Su). The obtained (σt) from the numerical analyses are substituted into equation (3) to calculate the critical pressure ratio (PR).

Results and Discussion
The blowout and collapse pressure ratios {PR = (σs − σt)/Su} versus the width to height ratio (B/D = 0.5 to 2) for various values of the depth ratio (C/D = 1-5) are presented in Figure 4a (for blowout) and Figure 4b (for collapse), respectively, for a weightless soil γD/Su = 0.
Due to the PR definition, a negative value of {PR = (σs − σt)/Su} means that the compressive normal supporting pressure (σt) is greater than the compressive surcharge pressure (σs). It is to be noted that negative PR values are seen for all blowout results. This is because (σt) must be greater than (σs) to cause a blowout failure. The larger the absolute value of PR, the greater the (σt). Numerical blowout results have shown that the pressure ratio (PR) increases, in negative PR, as (B/D) increases for all values of C/D. The increase in negative PR literally means a less critical blowout pressure (σt) is required to cause a blowout failure as (B/D) increases. Figure 4a. has also shown that the soil stability decreases (in negative PR) as the depth ratio (C/D) increases. Similar trends can be found in the blowout results (see Figures 5a, 6a and 7a) for dimensionless strength ratios (γD/Su = 1-3). As (γD/Su) increases, the system becomes heavier and therefore greater (σt) is required to cause a blowout.
In contrast, for the weightless collapse scenario, the PR values are always positive, i.e., (σs) must be greater than (σt) to cause a collapse. Numerical results in Figure 4b have shown that the pressure ratio (PR) decreases with increasing (B/D) for all the depth ratios (C/D). Furthermore, the soil stability increases as the depth ratio (C/D) increases. A symmetrical result can be found from both the collapse and blowout solutions.
Similar observations are noticed with the shear strength ratio (γD/Su = 1-3), as shown in Figures 5b, 6b and 7b. The "heavier" the system the larger the supporting pressure (σt) is needed to prevent collapse failure. As (γD/Su) increases, the value of PR becomes negative indicating the need of a supporting pressure for the active failure mechanism.   Using the same data, the effects of (C/D) are presented in Figures 8-12. In general, for the blowout scenario (Figures 8a-12a), the absolute values of the pressure ratio {PR = (σsσt)/Su} increases as (C/D) increases. Note that the rate of increase (in negative PR) is greater as (γD/Su) increases. The finding is applicable to the various values of (B/D) of the blowout cases. Similar trends can be found in the collapse scenario (Figures 8b-12b), where larger σt (in negative PR) is expected to support the soil as (γD/Su) increases. Furthermore, the larger the (C/D), the larger the σt (in negative PR) is. Design charts for assessing the blowout and collapse stability for various strength ratios (γD/Su = 0-3) are presented in Figures  13-16, in which the x-axis represents the width to height (B/D) and the y-axis represents the depth ratio (C/D). One would need to input the "designed" parameters (B/D) and (C/D) to obtain a critical pressure ratio (PR) directly from the charts.        Contour plots of the absolute velocity fields of the elliptical cavity shape transformation are presented in Figures 17-19 for various width to height ratios (B/D = 0.5, 1 and 2). These plots are both blowout and collapse modes with a shear strength ratio (γD/Su = 2) and two depth ratios (C/D = 1 and 3). In general, the failure extent in the collapse mode is greater than that in the blowout mode. This observation is supported by noting the location of the initial slip surface points from the inner cavity (see Figure 19). A chimneytype of failure is presented for the shallow depth ratio C/D = 1 in both the blowout and collapse scenarios and larger lateral expansions are seen in the deeper case C/D = 3, which supports the fact that greater σt is needed to cause soil blowout failure or to prevent soil collapse. Noting that the absolute values of the coloured velocity fields are not real for such a perfectly plasticity soil model, and therefore they are not presented here.

Comparison and Examples
Since there was no previous study on the soil stability of elliptical-shaped cavity under blowout condition, the comparison in Figure 20 is made by using previous circular cavity results (i.e., B = D) in either blowout or collapse scenarios. To compare the current pressure ratio (PR) with the critical stability number (Nc) in Shiau and Al-Asadi [35], a weightless soil condition (γD/Su = 0) is needed. In this way, our pressure ratio (PR) is equal to the (Nc) reported in Shiau and Al-Asadi [35]. The comparison shows a good agreement between Shiau and Al-Asadi [35] and the present study for the blowout study. Similarly, for the collapse study, it also shows a good arrangement amongst Shiau and Al-Asadi [35], Wilson et al. [36] and the present study. This comparison has greatly enhanced the confidence with the numerical results presented in the paper. An example is given next to explain how to use the design charts.
A simple example: Evaluate the ground stability with the following given parameters: γ = 20 kPa, Su = 40 kPa, σs = 100 kPa, C = 6 m, D = 2 m and B = 1 m. Refer to Figure 3 for the notations.
Blowout Check 1. Both the cover depth ratio and the width depth ratio are C/D = 3 and B/D =0.5, respectively.

Conclusions
Ground stability due to water main leakage-related sinkhole in blowout and collapse scenarios was investigated in this study using finite element limit analysis. Elliptical cavity shapes transformation was studied for a wide range of dimensionless ratios, namely the depth ratios, width to height ratios and soil strength ratios. The study was aimed to compute the rigorous upper and lower bounds of pressure ratios (PR). It was found that the thinner cavity (B/D = 0.5) provides better resistance than the wider one (B/D = 2), and the extent of surface failure increases with the increase in depth ratio irrespective of the opening shape. Numerical results obtained were compared with published literature, and comprehensive design charts and tables produced for practical uses with illustrated examples. It should be noted that the simulations in this paper are limited to 2D plane strain cavities which are close to the settlement phenomenon of the ground due to the excavation of unreinforced tunnels. The immediate future work will be the study on sinkholes with elliptical shape by using which the finite element limit analysis under axisymmetric conditions, which can provide a more realistic phenomenon of 3D sinkholes. In addition, future work recommendations may also include an elliptical-shaped transformation in cohesionless soil.