Bearing Capacity of Cylindrical Caissons in Cohesive-Frictional Soils Using Axisymmetric Finite Element Limit Analysis

This paper investigates the use of stability factors for estimating the ultimate bearing pressures on cylindrical caissons in cohesive-frictional soils. Rigorous upper and lower bound limit analyses with finite elements in axisymmetric (AX) condition are used for a series of numerical studies. The bearing capacity factors (Nc0, Nq0 and Nγ0) for a surface circular footing are firstly revisited. This is followed by a study on the effect of caisson’s embedded depth ratio (L/D). A comprehensive set of depth factors (Fcd, Fqd and Fγd) is then reported as a function of caisson’s embedded depth ratio (L/D) and soil internal friction angle (ϕ). The obtained results are compared with published solutions in the literature. Several examples are given to validate the principal of superposition as well as to illustrate on how to use the produced factors to estimate the ultimate bearing pressures on cylindrical caissons in cohesive-frictional soils. The study should be of great interests to practitioners.


Introduction
Cylindrical caissons, also known as skirted or bucket foundations, have been proven to be an effective foundation system for bridges, piers, floating platforms and other offshore foundations in deep water Mello et al. 2021;Sales et al. 2021). Caisson is generally made from a large steel thin-walled cylindrical structure that is open at the bottom and closed at the top. A comprehensive review of the offshore foundation was presented in Randolph and Gourvenec (2011), who reported that the depth ratio between the diameter and the depth of caissons is in the range of 1-6.
Early studies on the behaviour of caissons in cohesive soils were conducted by using field experiments and centrifuge model tests Dyvik et al. 1993;Clukey and Morrison 1993;Cauble 1996). Numerical methods such as using the displacement-based finite element and the finite element limit analysis were also reported by numerous researchers (see e.g., Bransby and Yun 2009;Geer 1996;Gourvenec 2008;Gourvenec and Barnett 2011;Mana et al. 2013;Ukritchon et al. 2018;Keawsawasvong and Lawongkerd 2021;Keawsawasvong and Ukritchon 2016; Abstract This paper investigates the use of stability factors for estimating the ultimate bearing pressures on cylindrical caissons in cohesive-frictional soils. Rigorous upper and lower bound limit analyses with finite elements in axisymmetric (AX) condition are used for a series of numerical studies. The bearing capacity factors (N c0 , N q0 and N γ0 ) for a surface circular footing are firstly revisited. This is followed by a study on the effect of caisson's embedded depth ratio (L/D). A comprehensive set of depth factors (F cd , F qd and F γd ) is then reported as a function of caisson's embedded depth ratio (L/D) and soil internal friction angle (ϕ). The obtained results are compared with published solutions in the literature. Several examples are given to validate the principal of superposition as well as to illustrate on how to use the produced factors to estimate the ultimate bearing pressures on cylindrical caissons in cohesive-frictional soils. The study should be of great interests to practitioners. Keawsawasvong 2016;Yun and Bransby 2007). Most of the research work done was for the bearing and pull-out capacity of caissons in undrained clays. Very few works were reported in relation to soils in drained condition. In addition, the studies on the behaviours of open caissons were also considered by Lai et al. (2020) and O'Dwyer et al. (2018. The bearing capacity of strip foundations with structural skirts was presented by Al-Aghbari and Mohamedzein (2004) who performed a series of model tests on these footings embedded in dense sand. Later, Al-Aghbari and Dutta (2008) and Eid et al. (2009) conducted several experiments to investigate the behaviours of square skirted footings in sand. A numerical finite element analysis was employed by Eid (2013) to evaluate the bearing capacity and settlement of skirted shallow foundations on sand. The effect of different types of sands (e.g., medium dense and dense sands) on the bearing capacity of circular skirted footings in sand was examined by Wakil (2013). Khatri et al. (2017) later performed a series of small-scale model test to study the behaviour of rectangular and square skirted footings on sand. Using the finite element limit analysis (FELA), Khatri and Kumar (2019) proposed the depth factor for soil weight F γd for circular and strip skirted footings with small embedment depth ratios (L/D) varying from 0 to 2.
Very recently, Shiau and Al-Asadi (2020a, b) adopted a stability factor approach, that is analogous to the traditional bearing capacity problem, to study the drained stability solutions of underground tunnelling. The method has been proven to be both efficient and effective. So far, there is no comprehensive results of the depth factors (F cd , F qd and F γd ) for the caisson problems in the literature. Following the stability factor approach in Shiau and Al-Asadi (2020c), the focus of this study is to assess the ultimate bearing pressures on cylindrical caissons in cohesivefrictional soils using the advanced upper and lower bound limit analysis. Comprehensive bearing capacity factors (N c0 , N q0 and N γ0 ) and depth factors (F cd , F qd and F γd ) are presented in tables and figures to assist designers and practising engineers in calculating the critical pressure that can apply to the caissons in their preliminary stage of design. Note that the bearing capacity factors (N c0 , N q0 and N γ0 ) are used for the cases of circular footings resting on the surface of soils. To consider the circular skirted footings with an embedment depth (L/D > 0), the depth factors (F cd , F qd and F γd ) are then adopted. This will be demonstrated in the example section of the paper.

Stability Factor Approach
The problem definition of a cylindrical caisson in axisymmetry (AX) for a cohesive-frictional soil is shown in Fig. 1. The caisson has a diameter D and an embedment depth L. The ground surface is subject to a vertical surcharge (q), while a bearing pressure at the top of the caisson is (q u ). The soil medium obeys the Mohr-Coulomb failure criterion with three parameters including cohesion (c), unit weight (γ), and friction angle (ϕ).
Terzaghi's bearing capacity equation, as shown in Eq. (1), can be used to calculate the ultimate uniform pressure q u applied at the top cap of a 2D plane strain caisson when L = 0.
where c is the soil cohesion; q is the surcharge loading; γ is the soil unit weight; N c is the 2D plane strain bearing capacity factor for cohesion; N q is the 2D plane strain bearing capacity factor for surcharge Line of axial symmetry qu q Fig. 1 The axisymmetric problem of a caisson loading; N γ is the 2D plane strain bearing capacity factor for soil weight. Equation (1) is further modified to Eq.
(2) to include the axisymmetric effects in 3D using (N c0 , N q0 and N γ0 ) as well as the depth factor (L/D) effects using (F cd , F qd and F γd ).
where (N c0 , N q0 , and N q0 ) are the axisymmetric bearing capacity factors for cohesion, surcharge, and soil weight respectively of a surface circular footing (L/D = 0). The three bearing capacity factors (N c0 , N q0 and N γ0 ) are a function of only soil internal friction angle (ϕ). On the other note, the depth factors (F cd , F qd and F γd ) of cylindrical caissons are a function of both the soil internal friction angle (ϕ) and the depth ratio (L/D) of the caissons. Interestingly, the three depth factors are unity when L/D = 0, and Eq. (2) reduces to Eq. (3) for the problem of surface circular footings.
As stated in Shiau and Al-Asadi (2020b), to obtain the respective factors such as (N c0 , N q0, and N γ0 ) and (F cd , F qd, and F γd ), it is necessary to impose zero value to c (cohesion) or q (surcharge) or γ (soil unit weight). This procedure works well to obtain individual pressure q u for cohesion, surcharge, and unit weight. They are then used to calculate the respective factor of a given depth ratio (L/D). The factor N c0 can be obtained by assuming no surface pressure (q = 0) and weightless soil (γ = 0) whereas the factor N q0 be acquired by assuming no cohesion (c = 0) and weightless soil (γ = 0). In addition, the factor N γ0 can also be obtained by assuming no cohesion (c = 0) and surface pressure (q = 0). The principal of superposition has been discussed and validated by Shiau and Al-Asadi (2020a), and in this paper again, has been successfully applied to the examples in a later section.

FELA Model
Recent advances in Finite Element Limit Analysis (FELA) with both upper bound (UB) and lower bound (LB) estimates are powerful as they can provide an error indicator towards the true collapse load (Sloan 2013). The development began with linear Sloan (1988Sloan ( , 1989. Nonlinear programming formulations provide better solution accuracy with shorter cpu solution time (Lyamin and Sloan 2002a, b;Krabbenhoft et al. 2007). Recently, the FELA (Optum CE 2021) has been widely used to solve a variety of drained and undrained stability problems in geotechnical engineering (Shiau and Smith 2006, Shiau et al. 2016a, b, c, 2021aAl-Asadi 2020c, d, 2021;Keawsawasvong and Ukritchon 2017a, b, 2020Ukritchon and Keawsawasvong 2017a, b, 2019a, b, 2020aUkritchon et al. 2019Ukritchon et al. , 2020Yodsomjai et al. 2021) and it was chosen in this study to compute the bearing capacity factors (N c0 , N q0 and N γ0 ) and the depth factors (F cd , F qd and F γd ) of the cylindrical caissons in axisymmetric condition.
In the LB method, three-node triangular elements are used in the analysis. Each triangular element has the nodal stresses of σ r , σ z , σ θ , and τ rz for an axisymmetric problem. The statically admissible stress discontinuities are allowed for producing the continuity of normal and shear stresses along with the interfaces of all the elements. The conditions of stress equilibrium, stress boundary condition, and the Mohr-Coulomb failure criterion are the constraints in a typical LB analysis, in which the objective function is to maximize the collapse load of problems. On the other hand, the upper bound theorem requires a kinematically admissible velocity field where the external work is greater or equal to the plastic shear dissipation. In the UB method, three-node triangular elements are used in the formulation. At each node of the element, there are the horizontal (u) and vertical velocities (v) defined as the basic unknown variables. The setting of kinematically admissible velocity discontinuities is applied at the interfaces of all the elements. The material is set to obey the associated flow rule. These two theorems are perfectly fitted to the nonlinear programming optimization problems using the second-order cone programming (SOCP). The constraints involved in this procedure are nonlinear and non-smooth but remain convex and amenable to analysis.
Since the problem of an cylindrical caisson can be modelled under an axisymmetric (AX) condition, only half of the domain is employed in the simulation, as shown in Fig. 2. The line of axial symmetry is set to be located at the left of the domain. Noting that the underlying bound theorems assume a rigid-perfectly plastic material with associated plasticity, the soil mass is discretised as triangular elements and modelled as Mohr-Coulomb material. The caisson is modelled by using rigid plate elements. The interface condition at the contact surface between the caisson and the soil is set to be fully rough. The feature "Fan Mesh" in the program is activated at the tip of the caisson to improve the solution accuracy (Krabbenhoft et al. 2015). The bottom boundary of the model is fixed in both horizontal and vertical directions, while the left and the right boundary can move in the vertical direction. The domain size was chosen to be large enough so that the solution is not affected by the development of the overall velocity field.
An automatically adaptive mesh refinement was employed in both the UB and LB simulations to compute the tight UB and LB solutions of the ultimate pressure q u . Figure 3 presents a typical example of the mesh refinement technique for a depth ratio of (L/D) = 2. The initial FELA mesh with approximately 5,000 elements is shown in Fig. 3a, whilst the final mesh with 10,000 elements is shown in Fig. 3b. It should be also noted that all presented numerical results hereafter are the average solutions from LB and UB FELA after the five adaptive mesh refinement steps with approximately 10,000 elements. N c0 , N q0 ,  and N γ0 Numerical results of the axisymmetric (AX) cohesion factor N c0 (L/D = 0) with soil frictional angles ϕ varying from 5° to 35° are presented in Fig. 4. Noting that N c0 increases nonlinearly with the increasing ϕ, they are in excellent agreement with the slip line solutions in Martin (2004) and Kumar and Khatri (2011). As  (2013), the slip line solutions by Martin (2004) can be considered as nearly exact solutions. The published N c0 results of Terzaghi (1943), Meyerhof (1963), and the displacement-based finite element of Chavda and Dodagoudar (2019) are consistently lower than those in the present study and they are considered as conservative solutions. Figure 5 presents numerical results of the axisymmetric (AX) surcharge factor N q0 (L/D = 0) with ϕ varying from 5° to 35°. Similar trends as in the cohesion factor N c0 , the nonlinearly increasing results have compared well with published ones, except those conservative solutions of Terzaghi (1943), Meyerhof (1963) and Chavda and Dodagoudar (2019). For the axisymmetric (AX) soil unit weight factor N γ0 (L/D = 0), again, numerical results in Fig. 6 have suggested great confidence in using the present averaged results of rigorous upper and lower bounds. Interestingly, the finite difference results reported by Zhao and Wang (2008) are unsafe (or unconservative) as they are consistently greater than our solutions. The values of N c0 , N q0, and N γ0 are explicitly shown in Table 1. It should be noted that the all presented results are in good accordance with the research of Kumar and Khatri (2011) since this study and Kumar and Khatri (2011) employed the same FELA method to solve the solutions of circular foundation problems. However, Kumar and Khatri (2011) only proposed the bearing capacity factors for surface circular foundations. The solutions of caissons with a given depth ratio (L/D) have never been presented in the past. To the best of our knowledge, this study is the first work to consider the depth factors F cd , F qd, and F γd for the caisson problem.

Discussing the Depth Factors F cd , F qd , and F γd
Numerical solutions of the axisymmetric (AX) cohesion depth factor F cd with soil frictional angles ϕ varying from 5° to 35° are presented in Fig. 7 for (L/D = 0 to 6). The value of F cd starts from one at L/D = 0, and it increases linearly with the increase in L/D for all considered fiction angles ϕ = 5° to 35°. The larger the ϕ, the greater the F cd . The effect of ϕ is insignificant when ϕ > 30° since both lines of ϕ = 30° and 35° overlap perfectly. Figure 8 shows that, for all friction angles, the surcharge depth factor F qd increases as L/D increases. This increase is approximately linear. The deeper the caisson is, the larger the F qd (surcharge effect) is. The greater the ϕ, the larger the F qd . For the soil unit weight depth factor F γd , see Fig. 9, a highly nonlinear relationship between F γd and L/D is observed. An increase in L/D yields an increase in F γd . Noting that F γd is the smallest when ϕ = 35°, whilst the largest when ϕ = 5°. This is dissimilar to the other factors such as F cd and F qd . A possible explanation for this might be due to the local punching failure mechanism near the end bearing point of the deep foundation. The complete values of (F cd , F qd, and F γd ) are explicitly presented in Tables 2, 3 and 4. Using the bearing capacity factors (N c0 , N q0, and N γ0 ) and the depth factors (F cd , F qd, and F γd ), Eq. (2) can be used by  practical engineers to evaluate the ultimate capacity q u of cylindrical caissons in cohesive-frictional soils. Figure 10 compares the unit weight depth factor F γd between the present study and the footing embedment study by Lyamin et al. (2007), as well as the study of skirt foundation on the sand by Khatri and Kumar (2019). Note that the work by Khatri and Kumar (2019) was for cylindrical skirt footings with the soils inside the skirt, whereas Lyamin et al. (2007) studied rigid embedded footings without the soils inside. The comparison shows a good agreement between the present results with those in Khatri and Kumar (2019).
Though not entirely the same as our current study, the solutions provided by Lyamin et al. (2007) are conservative as they predict lower F γd factors than those in Khatri and Kumar (2019) and in this study.

Examples
Several examples are presented in this section to demonstrate how to use the produced results to evaluate the uniform bearing capacity of cylindrical caissons by using the formulation shown in Eq.
(2). The principal of superposition using the stability factor approach is also validated with the following examples.

Example 1: Cohesionless Soil Without Surcharge Loading
A cylindrical caisson has the length L = 12 m and the diameter D = 3 m. The design parameters are given as: the unit weight γ = 18 kPa and the soil internal friction angle ϕ = 35°. The soil cohesion c is zero since the soil is cohesionless. The surcharge loading q is also zero in this example. Given ϕ = 35°, the value of N γ0 = 41.298 is obtained from Table 1 and F γd = 73.189 is also obtained from Table 4 using L/D = 4. Note that the values of N c0 , N q0 , F cd and F qd are not required since c = 0 kPa and q = 0 kPa. Using Eq. 2, the ultimate uniform An actual analysis in the program using the real parameters gives q u = 95,164.657 kPa, which is 1.2% greater than the solution using stability factors and the principal of superposition.

Example 3: Cohesive-Frictional Soil Without Surcharge Loading
In this example, the cylindrical caisson has a diameter D = 4 m (i.e. L/D = 3). The design parameters are given as the unit weight γ = 16 kPa, the soil internal friction angle ϕ = 10°, and the soil cohesion c = 25 kPa. There is no surcharge loading at the ground surface so that q = 0 kPa. From Table 1 An actual computer analysis of the problem gives q u = 2,567.307, which is very close to the solution using stability factors and the principal of superposition.
The above examples have proven that the principal of superposition with the stability factors and the depth factors can be used to evaluate the unltimate bearing pressures of caissons.

Conclusions
This paper has successfully produced lower and upper bound solutions of the bearing capacity factors (N c0 , N q0, and N γ0 ) and the depth factors (F cd , F qd, and F γd ) for cylindrical caissons in cohesive-frictional soil. The following conclusions are drawn based on this study: 1. The axisymmetric bearing capacity factors (N c0 , N q0, and N γ0 ) for a surface footing are a function of only the soil friction angle (ϕ). An increase in the soil friction angle (ϕ) results in an increase in all three bearing capacity factors (N c0 , N q0, and N γ0 ). 2. The depth factors (F cd , F qd, and F γd ) are functions of both the soil friction angle (ϕ) and the depth ratio (L/D). An increase in the depth ratio (L/D) causes an increase in all depth factors (F cd , F qd, and F γd ). 3. For the unit weight depth factor, F γd , an increase in the soil friction angle (ϕ) results in a decrease in F γd . This is different from F cd and F qd, where they increase with the increasing soil friction angle (ϕ). A possible explanation for this might be due to the local punching failure near the end bearing point of the foundation. 4. The present solutions of bearing capacity factors (N c0 , N q0, and N γ0 ) are in good agreement with the previous solutions using the method of characteristic. In addition, the present solutions of the depth factor F γd are also in good agreement with those published ones. 5. The illustrated examples using the principal of superposition have proven that the use of the bearing capacity factors (N c0 , N q0, and N γ0 ) and the depth factors (F cd , F qd, and F γd ) is both convenient and accurate. 6. The proposed rigorous solutions of the bearing capacity factors (N c0 , N q0, and N γ0 ) and the depth factors (F cd , F qd, and F γd ) are useful and they can be used by pratising engineers to evaluate bearing capacity requirements of caissons in cohesivefrictional soils with great confidence.
The results in this study are applicable to skirted circular footing subjected to drained loading. For the open caissons which are generally sunken within the ground and have an opening at the bottom and top (during sinking), the current drained solutions may not be suitable. In addition, the interface between caissons and soils is only set to be fully rough, where the adhesion factor is set to be one. The influence of the adhesion factor on the bearing capacity factors may be significant, and future works should consider the full range of the adhesion factor (0 to 1) for more realistic simulations.