#### Abstract

The objective of this paper is to investigate the effect of soil variability on bearing capacity of an embedded foundation in the presence of nonstationary undrained shear strength. The nonstationary undrained shear strength is simulated by a nonstationary random field generator based on the spectral representation method. An embedded foundation buried into the soil to two times of width is presented to investigate the influence of spatially variable undrained shear strength on bearing capacity. Firstly, Monte Carlo simulations are carried out to discuss the effect of distribution type, nonstationary gradient parameter, and horizontal autocorrelation length on the bearing capacity from the standpoint of mean value and standard deviation. Then, the influence of the distribution type on the failure probability of nonstationary random soil is also investigated, with the failure probability for the Beta distribution being demonstrated to be always larger than that for the Lognormal and Gamma distribution.

#### 1. Introduction

It is widely recognized that soil exhibits spatial variability in both horizontal and vertical direction due to the effects of chemical and environmental changes and soil depositions. This inherent spatial variability of soil becomes the dominating source of uncertainty for the estimation of soil properties from in situ tests [1–3]. Based on this nature, random field theory and random finite element method (RFEM) were developed to consider the inherent random variation of soil properties in geotechnical problems, such as the bearing capacity of foundations, the stability of slopes, and soil liquefactions. The past studies have clearly shown the importance of considering the spatial variability of soil property for a rational geotechnical design.

The bearing capacity of a shallow foundation is a classical geotechnical problem and has been further studied by probabilistic analysis methods that deal with the spatial variability of soil properties. For example, Popescu et al. [4] and Griffiths and Fenton [5] employed Monte Carlo simulation through RFEM of a shallow foundation and found a substantial reduction in the bearing capacity of the shallow foundation when considering the spatial variability of soil properties. In other words, the overestimation of the bearing capacity of shallow foundations indeed exists when not accounting for the spatial variability of soils. With the development of the offshore engineering in recent years, the embedded foundations have been widely applied in offshore engineering, especially in soft clay seabed soil. Taking the large circular footings of mobile jack-up platforms as an example, these footings can be buried into soft clay seabed to depths of several times of diameters [6–8]. The failure mechanism of these deeply embedded foundations and thus the bearing capacity differ from those of shallow surface footings studied previously [9, 10]. The failure mechanism of soil transforms from a general shear failure for a shallow foundation [11] to a localized failure for a deeply buried foundation [12, 13], with the ultimate bearing capacity for a deeply buried foundation being demonstrated to be much larger than that of a shallow foundation [14].

In the existing literature, the stationary random field has been widely applied to describe the spatial variability of soil parameters. However, it is well known that soil parameters often exhibit nonstationary characteristics; i.e., it was confirmed by some in situ tests. Lumb [15] pointed out two typical nonstationary types of soil properties: (a) mean value of soil parameter linearly increases with depth while with its standard deviation it keeps constant and (b) both mean value and standard deviation of soil parameter linearly increase with depth while its coefficient of variation (COV) keeps constant. The undrained shear strength was reported to accord with type (b) because it depends on the effective stress of soil [9]; nevertheless, only a few recent kinds of research have considered the nonstationary characteristics of soil properties in geotechnical probabilistic analysis. For instance, Li et al. [16] investigated the reliability of a shallow strip footing in the presence of spatially variable undrained shear strength that linearly increases with depth and found that ignoring the nonstationary characteristic of undrained shear strength will overestimate the failure probability of the footing. Wu et al. [17] adopted different nonstationary extents of random fields of the undrained shear strength associated with a shallow foundation. They found that the mean value and standard deviation of the bearing capacity increase linearly with strength gradient parameter *M.* Wu et al. [18] have systematically studied how the degree of nonstationarity of the undrained shear strength and the embedded depth cross-affect the statistics of bearing capacity and failure mechanism of soil. In this paper, due to the lack of relevant researches and more for the accuracy of the estimated bearing capacity, the nonstationary random field model, according to type (b), and the stationary random field model will be both adopted to simulate the undrained shear strength of the soil and investigate the resulting effects of spatial variability on the bearing capacity of a deeply buried foundation.

Theoretically, soil properties follow non-Gaussian probability distributions in the probabilistic analysis because physical material properties cannot be assumed as negative values [4]. Hence, the Lognormal distribution is often selected to fit the undrained shear strength of the soil. However, according to several studies reported in the literature, undrained shear strength can follow different probability distribution functions (PDFs) for different types of soils and sites [19]. Beta and Gamma distribution, which also belong to a non-Gaussian probability distribution, will considerably impact and make the estimated statistics of bearing capacity different from that of the Lognormal distribution when employed to define the soil’s undrained shear strength [19]. As a result, this study will consider three different probability distribution functions of the undrained shear strength, i.e., Lognormal, Gamma, and Beta distribution.

This paper aims to investigate the bearing capacity of a deeply embedded foundation in the presence of spatially variable undrained shear strength with different probability distribution functions. A nonlinear finite element model combined with nonstationary random fields, discretized by a random field generator (i.e., spectral representation method), is employed to achieve this aim through Monte Carlo simulation. For comparison with the results obtained from the nonstationary random fields, the results obtained from the corresponding stationary random fields under different PDFs are also provided. With the methodology above, the influences of the distribution type of the undrained shear strength under the condition of both stationarity and nonstationarity on the reliability and bearing capacity of an embedded foundation are illustrated in this paper.

#### 2. Finite Element Model

In the current study, a finite element model was performed to compute the bearing capacity of an embedded foundation under the condition of two-dimensional plane strain in this study. As shown in Figure 1, this study has considered an embedded foundation of width *B* and height *h* embedded at the depth of *D* = 2*B*, as analyzed by Li et al. [14]. According to the typical size of a spudcan foundation, the foundation width (*B*) is 20 m and its height (*h*) is 4 m, respectively. The domain of the model has a width of 6.4*B* and a height of 6*B* to ensure no obvious boundary effects. A mesh scheme consisting of 22240 elements was adopted for the bearing capacity analysis. The mesh scheme discretized the analysis domain into elements, with sizes of 2 m in the horizontal direction and 0.5 m in the vertical direction, and further discretized the region of 3*B* × 2*B* close to the footing into smaller finite elements of 0.5 m × 0.5 m. For boundary conditions, the horizontal displacement of vertical sides was constrained and the displacement of bottom nodes was fixed in all directions.

The embedded footing was considered to be rigid and has a rough interface with the soil. The soil was modeled by a linear-elastic perfectly plastic constitutive law. The failure is according to the Tresca yield criterion, with the maximum shear stress in any plane being limited to the undrained shear strength.

The parameters of Young’s modulus (*E*) and undrained shear strength (*c*_{u}) for soil model were both considered to be stochastic parameters while only the undrained shear strength was simulated by the random field generator. The undrained shear strength could be assumed to be perfectly correlated to Young’s modulus based on the relationship, i.e., *E*/*c*_{u} = 500 [20]; hence, Young’s modulus was also a random field. Poisson’s ratio for soil model was set to hold constant, i.e., 0.49, to model the undrained conditions of no volume change as well as to ensure numerical stability. The basic parameters describing the stochastic characteristics of the undrained shear strength and the simulating procedure of the random fields are presented as follows.

#### 3. Simulating Random Fields of Soil Property

In the process of Monte Carlo simulation, because of the discrete nature of numerical simulations, a continuous-parameter random filed must be discretized into random variables by a numerical method. Several discretization methods have been developed to attain this aim [21–24], and among them, the spectral representation method (SRM) is widely applied due to its good features in convergences. Therefore, the SRM is adopted to discretize anisotropic random fields of soil property in this study. The non-Gaussian distribution is often used to model the nonnegative data of soil properties, and soil parameters often exhibit nonstationary characteristics (i.e., the mean value increases with depth). The SRM method proposed by Wu et al. [25–28] can be easily and high-efficiently extended to generate the non-Gaussian stationary random fields. For more details of the general steps for simulating the realizations of a non-Gaussian stationary field with zero mean value and unit variance, the paper given by Wu et al. [25] can be referred. Once the standard non-Gaussian stationary field is generated, the desired non-Gaussian nonstationary random field *y* (*x*, *z*), whose mean value and standard deviation increase linearly with depth, but with a constant COV, can be obtained bywhere and are the target mean value and standard deviation of the desired random field at depth *z*, respectively.

In this study, to investigate the effect of distribution type of undrained shear strength *c*_{u} on the bearing capacity of the embedded foundation, the Lognormal, Gamma, and Beta distributed nonstationary and stationary random fields with the mean value , constant COV, horizontal correlation length *l*_{h}, and vertical correlation length were simulated with the simulating procedure mentioned above. The PDFs of the three standard non-Gaussian distribution with zero mean value and unit standard deviation are plotted in Figure 2.

To describe the nonstationary extent of the mean undrained shear strength , a dimensionless gradient strength parameter *M* is given aswhere is the mean value of *c*_{u} at the bottom of the model, where the depth *z* is equal to the height of foundation soil *H* (i.e., *H* = 6 *B*). Similarly, represents the mean value of *c*_{u} at the surface level (i.e., *z* = 0 m). The is set to be 10 kPa for both stationary and nonstationary random fields. Then, the corresponding dimensionless gradient strength parameter *M* for stationary and nonstationary random fields is set to be equal to 0 and 12, respectively. In fact, at the depth *z* can be expressed by , in which *k* denotes the increasing gradient with the depth and has a unit of kN/m^{3}. In this study, the value of *k* was set to be 0 and 1.

A squared exponential autocorrelation 2D autocorrelation function was adopted to define the autocorrelation structure of *c*_{u}:where and are lag distance in horizontal and vertical directions, respectively. The horizontal and vertical correlation distances, i.e., *l*_{h} and , were set to be 60 m and 6 m, respectively. The parameters involved in the random simulation of shear strength are listed in Table 1.

Based on the coordinate of elements in the finite element model, a series of realizations of the stationary and nonstationary random fields of *c*_{u} are generated by the simulation method and involved parameters mentioned above. In this study, 1000 simulations, which have been enough to ensure the estimated error to be allowed, have been performed for each random case.

#### 4. Analysis of Simulation Results

The MCS was employed to conduct probability analysis and evaluate the statistical properties of the bearing capacity of the embedded foundation in both stationary and nonstationary soils. In the framework of MSC, the 1000 simulations of *c*_{u} for each type of random field, coupled with the fixed finite element model, were recalled by a batch program to compute the bearing capacity on each of them. It is worth noting that, according to the predecessors’ researches [29], the normalized bearing capacity is usually discussed rather than the actual bearing capacity. For an embedded foundation, the actual bearing capacity should be nondimensionalized by the mean shear strength at the embedded depth of foundation because the mean bearing capacity directly relates to the soil around the foundation [17]. Hence, in the statistical analysis of this paper, the bearing capacity factor for the *i*th realization can be calculated usingwhere is the computed bearing capacity for the *i*th realization and represents the mean strength at embedded depth 2*B*.

It is worth noting that, in the statistical analysis, the normalized bearing capacity is usually discussed rather than the actual bearing capacity. For a stationary case, the bearing capacity normalized by the mean undrained shear strength is widely accepted as done by Griffiths and Fenton [29]. However, for a nonstationary case, there are different ways to normalize the bearing capacity, i.e., normalizing by the constant surface shear strength, normalizing by the mean shear strength of the whole domain, and normalizing by the mean shear strength at the embedded depth, which may lead to different trends relating to the influence of nonstationary characteristics on the statistical results. Since the mean bearing capacity was directly related to the soil strength around the foundation, the computed bearing capacity will be normalized by the mean shear strength at the embedded depth of foundations.

##### 4.1. The First Two Statistical Moments of Bearing Capacity

A series of computed bearing capacity of the embedded foundation is illustrated in Figures 3 and 4, which lead to the following analysis concerning the effect of the spatial variability of *c*_{u} of the foundation soil on the bearing capacity from the view of mean value and standard variation. It can be seen from Figures 3(a) and 3(b) and Figures 4(a) and 4(b) that the estimated mean bearing capacity considering the spatial variability of *c*_{u} is smaller than the bearing capacity for deterministic analysis. Note that *c*_{u} of the deterministic analysis corresponding to the nonstationary random fields linearly increases with depth according to *M* = 12.

**(a)**

**(b)**

**(c)**

**(d)**

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

###### 4.1.1. Effect of Probability Distribution Function

It is obvious that the distribution type of *c*_{u} has a significant influence on the bearing capacity, although the different distribution types are assigned identical mean value and standard variation of *c*_{u} in the process of random field simulation. Particularly, from Figures 3(a) and 3(c), it is easily observed that Beta distribution shows a great difference of the estimated mean bearing capacity and the corresponding standard variation from that of Lognormal and Gamma distribution. Furthermore, for the random fields with *l*_{h} = 60 m, the statistics of bearing capacity for the Beta distributed random soil exhibits the smallest mean value and largest standard deviation among the three distribution cases. A similar observation about the difference between the estimated mean bearing capacity and corresponding standard variation of Beta distribution type and that of Lognormal and Gamma distribution was also found for shallow foundation [11]. This may be due to the reason that both Gamma distribution and Lognormal distribution have richer right tails and more values distributed around the mean *c*_{u} than Beta distribution. Hence, *c*_{u} of Lognormal and Gamma distribution produces larger and more concentrative bearing capacity than *c*_{u} of the Beta distribution. However, for the random fields with *l*_{h} = 6 m, the standard deviation of bearing capacity is demonstrated to not conform to the above rule.

###### 4.1.2. Effect of M

It can be seen from Figures 4(a) and 4(b) that the mean value of bearing capacity for the nonstationary case is smaller than that of the stationary case, while the corresponding standard deviation for the nonstationary case is larger than that of the stationary case. Figure 5 shows the development of failure planes for a stationary random soil and a nonstationary random soil with the same seed. From Figures 5(a)–5(d), it is found that for the stationary case, the largest shear strain develops from the soil surrounding the foundation in a circular manner and finally forms a complete localized plane around both the left side and right side of the foundation. However, in the nonstationary case (shown in Figures 5(e)–5(h)), the largest strain appears at the four corners of foundation at first and then the shear planes at the upper two corners develop upward rather than around the foundation. Finally, it is seen that a complete shear plane throughout the left side of the foundation to the soil surface is formed, mobilizing a smaller area of soil at the left side of foundation comparable to the mobilized area at both sides in the stationary soil. The reason lies in that the shear plane in the nonstationary case, where the soil strength above foundation becomes smaller with decreasing the depth, extends along the relatively weaker soil; thus, it extends upward to the soil surface easily. Therefore, the soil with nonstationarity exhibits smaller resistance and bearing capacity than the corresponding stationary case.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(g)**

**(h)**

###### 4.1.3. Effect of the Horizontal Correlation Length

According to Li et al. [16], the autocorrelation length can bring about the changes of the mean and standard deviation of bearing capacity. Hence, the bearing capacity of the embedded foundation with different horizontal correlation lengths was provided in this paper, as Figure 4 shows. Both the mean value and standard deviation of bearing capacity for the stationary and nonstationary cases increase with the increasing of the horizontal correlation length. The similar conclusions with respect to a shallow foundation were drawn by Wu et al. [17] for nonstationary and stationary random field cases and Ye et al. [10] for nonstationary cases.

##### 4.2. Statistical Distribution Characteristics of Bearing Capacity

In this section, the statistical distribution characteristics of bearing capacity are studied. Figure 6 shows the cumulative distributions (i.e., CDFs) of the bearing capacity for different values of *M*. From Figure 6, for the three PDFs considered in this paper, especially when *l*_{h} = 60 m, the CDFs of the bearing capacity for a foundation embedded at 2*B* depth all move to the left side when *M* increases from 0 to 12, which indicates that the bearing capacity will be overestimated when the significant nonstationarity is completely neglected. This conclusion seems not to be obviously observed in Wu et al. [18], perhaps because the degree of nonstationarity for *c*_{u} in this paper is much larger than that in Wu et al. [18]. Besides this, the CDFs of the bearing capacity for the Beta distribution are found to be always on the left of the CDFs for Gamma and Lognormal distribution for both the stationary and nonstationary cases, which demonstrates that the Beta distribution produces a large amount of bearing capacity smaller than the Lognormal and Gamma distribution (as shown in Figure 7).

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(a)**

**(b)**

**(c)**

**(d)**

To better investigate the influence of the distributions of *c*_{u} on the foundation’s bearing capacity, the differential bearing capacity between the stationary and nonstationary cases will be further discussed as follows:where and are the bearing capacity computed for the *i*th realization of the stationary (i.e., *M* = 0) and nonstationary random fields (i.e., *M* = 12), respectively. It should be noted that, to make the comparisons reasonable, the *i*th realization of the stationary random field and nonstationary random field is mapped from the same “seed” random field with zero mean and unit deviation.

The CDFs of differential bearing capacity for all involved cases can be obtained and illustrated in Figure 8. It is found that when the vertical correlation length *l*_{h} is set to 60 m, in most realizations, the bearing capacity for stationary cases is larger than that for nonstationary cases and there are quite a few realizations, in which the bearing capacity for stationary cases is smaller than that for nonstationary cases. That is to say, of most realizations are larger than 0. However, the cases of smaller than 0 do exist more than larger than 0 when *l*_{h} = 6 m. In addition, *c*_{u} with Beta distribution produces a larger amount of higher for nonstationary cases than that for stationary cases and more discrete than *c*_{u} with Lognormal and Gamma distribution, which is more obviously observed in the cases of *l*_{h} = 60 m. Therefore, if we simply use the value of soil strength at depth of 2*B* where the foundation is buried to perform the design, more conservative results will be given for *c*_{u} with Beta distribution than *c*_{u} with Lognormal and Gamma distribution.

**(a)**

**(b)**

##### 4.3. Probability of Failure for Foundation

In practice, it is essential to know the potential failure probability of soil foundation from a design standpoint. According to Griffiths and Fenton [29], if the computed bearing capacity for a certain random realization is less than the bearing capacity for the deterministic case, it can be assumed that a failure of the foundation occurs. Then, the failure probability can be deemed as the probability of is less than 1. In other words, the failure probability can be easily obtained from the cumulative probabilities of normalized bearing capacity of the foundation , as shown in Figure 9. Note that, for the nonstationary case, was computed deterministically from the case with *c*_{u} linearly increasing with depth according to *M* = 12.

**(a)**

**(b)**

**(c)**

**(d)**

For a normally distributed bearing capacity, the probability of failure can be estimated by the following formula:where is the normal cumulative distribution function and and are the mean value and the standard deviation of the bearing capacity, respectively.

Note that the probability of failure utilized in this study is based on normalized bearing capacity rather than the actual bearing capacity. This defined probability of failure only reflects the comparison between bearing capacity from probabilistic analyses and those from the deterministic analysis. It does not reflect the real safety of the foundation.

In Figure 9, taking the cases of *l*_{h} = 60 m as an example, the failure probabilities are 86.4, 86.3, and 87.7% for the stationary *c*_{u} satisfied with Lognormal, Gamma, and Beta distribution, respectively, and 85.9, 85.8, and 86.6% for the nonstationary *c*_{u} satisfied with Lognormal, Gamma, and Beta distribution, respectively. It is obvious that the failure probability of the Beta distribution is the largest one among the probabilities of three distributions and the failure probability decreases when the nonstationary factor *M* turns from 0 to 12. When the vertical correlation length is 6 m, the failure probabilities are nearly 100% and thus unable to observe the influence of the PDFs with *c*_{u} and the nature of nonstationarity on the failure probability.

However, the above failure probabilities are too large to accept when directly using the bearing capacity of the corresponding uniform model to evaluate the failure probability. The probability that the bearing capacity is less than a certain level of applied load is widely utilized by introducing the factor of safety (FS). Therefore, the probability of failure in the form of with FS > 1 can be greatly reduced.

In Figure 10, the failure probabilities with FS = 1.3 for involved cases with *l*_{h} = 60 m are illustrated. The failure probabilities are 2.5, 1.85, and 4% for the stationary *c*_{u} satisfied with Lognormal, Gamma, and Beta distribution, respectively, and 6.4, 5.4, and 9.2% for the nonstationary *c*_{u} satisfied with Lognormal, Gamma, and Beta distribution, respectively. It was found that the failure probability for *c*_{u} with Beta distribution is still larger than that with Lognormal and Gamma distribution when FS is 1.3. However, it is interesting to find that, in contrast to the failure probability obtained from Figure 9, the failure probability is smaller for the bearing capacity of stationary cases than that of nonstationary cases when FS is 1.3.

**(a)**

**(b)**

To investigate the variation of failure probability for different distribution types with the value of FS, the curves of relations between the failure probability of nonstationary random soil and the value of FS are shown in Figure 11. For the random soil with *l*_{h} = 60 m, it is found that the failure probability decreases from about 10^{0} to 10^{−3} when the factor of safety decreases from 1.0 to 1.5; meanwhile, *c*_{u} with Beta distribution always exhibits higher failure probability than *c*_{u} with Lognormal and Gamma distribution. In addition to the expected trend showing the influence of the PDFs on the failure probability when *l*_{h} = 6 m, Figure 11(b) also indicates that the influence of PDFs on the failure probability is more obvious than that with *l*_{h} = 60 m. Hence, it can be concluded that from the view of probability analysis, the Beta distributed shear strength may result in an overconservative failure probability, i.e., higher failure probability, of soil foundation than the Lognormal and Gamma distributed shear strength. Although under a certain factor of safety can the failure probability of soil for Lognormal and Gamma distribution satisfy the target probability range of 10^{−2} to 10^{−3}, the failure probability for Beta distribution under the same FS may fail to meet the target probability.

**(a)**

**(b)**

#### 5. Conclusions

This paper focuses on investigating the stochastic results of bearing capacity for a deeply embedded foundation, buried at depth of 2*B*, in the presence of the nonstationary random field of undrained shear strength. The nonstationary random fields of *c*_{u} with Lognormal, Gamma, and Beta distribution and different vertical correlation lengths were generated by a generator based on the spectral representation method. Parametric studies were studied to evaluate the influence of the probability of function (PDF) of *c*_{u} in both stationary and nonstationary stochastic soils. According to this study, the following conclusions were drawn:(1)The PDF of *c*_{u} has a significant effect on the estimated mean value and standard deviation of the computed bearing capacity of the 2*B* embedded foundation. Particularly, the mean value for the Beta distribution is much smaller than that for Lognormal and Gamma distribution, while the standard deviation for Beta distribution is much larger than that for other cases when the vertical correlation length is 6 m.(2)The bearing capacity in nonstationary soils exhibits greater differences than that in stationary soils. It was found that, compared to a stationary spatial variable soil foundation, the failure plane of nonstationary soils extends more easily to soil surface rather than around the foundation and thus results in a smaller mean bearing capacity.(3)Compared with the Lognormal and Gamma distributed random fields of undrained shear strength, the Beta distributed random field may result in an overconservative failure probability in the bearing capacity. If we simply use the value of soil strength at depth of 2*B* where the foundation is buried to perform the design, more conservative results will be given for *c*_{u} with Beta distribution than *c*_{u} with Lognormal and Gamma distribution.

#### Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.

#### Acknowledgments

The support by the Foundation of Inner Mongolia Power (Group) Co., Ltd. (Grant no. 2019097) is greatly acknowledged.