Spatially variable soils affecting geotechnical strip foundation design

Joanna Pieczy′nska-Koz?owska, Giovanna Vessia

a Department of Civil Engineering, Wroc?aw University of Science and Technology, Wroc?aw, 50-370, Poland

b Department of Engineering and Geology, Università degli Studi G. d’Annunzio Chieti-Pescara, Chieti, 66100, Italy

Keywords:Bearing capacity Shallow foundation Random finite element method (RFEM)Ordinary kriging (OK)Soil property variability structure Reliability-based design

A B S T R A C T

1. Introduction

Soil characterisation for geotechnical design commonly deals with large uncertainties due to inherent variability(e.g.Phoon and Kulhawy, 1999a; Vessia et al., 2017), spatial variability of natural sediments(Phoon and Kulhawy,1999b;Cherubini et al.,2007),and uncertainty related to design models (Lesny et al., 2017; Mo et al.,2021; Tang and Phoon, 2021). To consider most of the aforementioned uncertainties,the random finite element method(RFEM)is often used(Griffiths and Fenton,1993).This method falls within the fully probabilistic approaches according to the current international standards for reliability design in civil engineering(ISO 2394,2015).However,a critical point is the proper estimation of the true probability distributions of soil parameters. Considering the stochastic geotechnical design of shallow foundation bearing capacity,the RFEM is a suitable method(Griffiths and Fenton,2001;Fenton and Griffiths,2003;Vessia et al.,2009;Pieczy′nska-Koz?owska et al.,2015) to handle random variables as random fields (Vanmarcke,1983). This assumption implies that the following features must be estimated for every design parameter or related measured variable:the scale of fluctuation(δ),the autocorrelation function,and the global variance. Several applications of RFEM have been devoted to analysing the worst-case in bearing capacity design(that is,the minimum value of the bearing capacity) related to the value of the scale of fluctuation. For several cases, when isotropic scales of fluctuation(δH=δv)are considered,the worst-case for the bearing capacity value of shallow foundation equals to the footing width (B) (e.g. Vessia et al., 2009). When actual anisotropic cases are considered (δH≠δv) (e.g. Fenton and Griffiths, 2003; Soubra et al., 2008), the worst-case varies case by case. It is mostly affected by both the coefficient of variation (COV) and the δH/δvratio of the design soil parameters.For this reason,hereinafter,the anisotropic conditions are considered, and the spatial variability structure of soil property random fields is effectively described by geostatistical tools, such as the Kriging interpolator (Matheron,1973). The approach presented in this paper is a complete answer to how to estimate strength parameters from large databases and design foundations based on soil spatial variability. The paper shows the relationship between the foundation width and the spatial variability of the soil. A similar problem, that is to find out the relationship between the soil spatial variability and the dimension of the foundation,was also addressed by Chwa?a(2021).

In the present study, the spatial variability structure of soil properties obtained from 182 cone penetration tests (CPTs) performed in Bologna Province, Italy, has been derived through the ordinary Kriging(OK)approach(Pieczy′nska-Koz?owska et al.,2017;Vessia et al.,2020a,b).Hereinafter,only the first 10 m depth of the a three-dimensional (3D) model based on the cone tip (qc) and the shaft resistance (fs) of CPT readings has been used. Then, RFEM analyzes calculated mean values and standard deviations of the bearing capacity under undrained and drained conditions according to the random field characters of the mechanical properties of fine soils,such as cohesion(c′)and friction angle(φ′).The scales of fluctuation have been calculated by the OK method and the other two methods, whereas two simple models have been considered for the autocorrelation functions:Markov and squared exponential(Gaussian). Subsequently, in Sections 2-4, methods and material used in this study have been introduced, in Section 5, the method used for the calculation of the bearing capacity of shallow foundations is illustrated,and the input variables are described through the random field theory.Section 6 shows the results and discusses the mean and standard deviation of the bearing capacity for shallow foundations with different widths.Finally,conclusions are drawn, and some key points are highlighted.

2. RFEM

The RFEM consists of three components: (1) the random field theory used to model the soil parameters, the first application of which was conceived and proposed by Vanmarcke (1977,1983) to soil modeling; (2) the classical finite element method (FEM) for calculation; and (3) the Monte Carlo method to generate a set of realisations of the random parameters and to estimate the statistical moments of the calculated variables(i.e.the bearing capacity).As the basis of the soil modeling used in RFEM, the random field theory assumes that each variable is defined by an autocorrelation function and a probability distribution.The fields might have one-dimensioan(1D),two-dimensional(2D),and 3D characteristics.The current task adopts 2D fields and the associated FEM analysis in the plane strain state. Combining the random field model with the FEM mesh, the local average subdivision(LAS)method is implemented(Fenton and Vanmarce, 1990). This method discretises the random field into a given FEM mesh considering the number of elements and their sizes in different directions(Lh-the size in the horizontal direction;Lvthe size in the vertical direction).In the meshed field,the mean value remains constant while the variance is reduced using the reduction function(γvar),considering the size of the scale of fluctuation:

where τhand τvare distances between paired points in horizontal and vertical directions, respectively. The variance value of a single element changes with the size of the element. Within a single realisation, feature values are randomly assigned. FEM analysis is carried out for the discretised mesh. Repeating this process many times and using the classical Monte Carlo Method, convergence of the estimated values to the mean solution can be reached.

3. Case study of the alluvial sediment in Po River plain

The 3D mechanical model of the subsoil consists of cone resistance(qc)and sleeve resistance(fs)values which have been obtained from 182 CPT soundings performed in Bologna Province, Italy. It covers about 900 km2, located in the southern portion of Padania Plain in the Emilia Romagna Region.It is the largest alluvial plain in Italy.Details about the dataset and the original database ofqcandfsprofiles can be found in Vessia et al.(2020b)and Di Curzio and Vessia(2021). These deposits are made up of undifferentiated mixtures of clay,sand,and silt with embedded gravelly bodies of alluvial fans in the south of the investigated area (Fig.1). Moving northward, narrower sandy paleo-channels (i.e. fluvial deposits) and silty-clayey lenses (i.e. lacustrine deposits) can be found. From this large area and dataset,11 cone and shaft resistance readings were selected.The selected CPTs are reported in Fig.1,limited to a study area of 57 km2in the northwestern portion of the investigated area.

Most of the profiles start from 1 m depth under the ground level and reach 20 m depth.A lithological selection of the soundings was made based on the average value of the soil behaviour type index(ISBT)(Robertson,2009)(red lines in Fig.2).The first 10 m of these profiles are considered in the present study.Along with this depth,the layering has not been considered according to the lithological structure of the soil mixtures (Pieczy′nska-Koz?owska et al., 2017).The averageISBTvalues in the first 10 m fall in zone 3 (Fig. 2a).According to Robertson (2009)’s classification, these alluvial deposits are made of fine soils, such as clay and silty clay.

Fig. 2. (a) ISBT for selected CPT soundings from the Emilia Romagna; and (b) ISBT at a depth of 10 m below ground level located in zone 3.

Therefore, the variability structure of a unique layer of 10 m depth has been estimated to describe the random input fields of the RFEM code(Griffiths and Fenton,1993;Fenton and Griffiths,2008).The fields have been used to investigate the bearing capacity of shallow foundation under both undrained and drained conditions.Besides, the stochastic input variables used in RFEM analyses include the cohesion(c′),the friction angle(φ′),and the undrained shear strength (su.).

Fig.1. The study area with the CPTs selected for numerical analyses.

4. RFEM numerical simulations

CPT soundings are widely used to display soil mechanical properties due to their quasi-continuous characteristics. They are also particularly useful when the spatial variability of soil properties is investigated. They record the mechanical variations with depth of the deterministic trend due to the increasing the confining stresses and the inherent variability of soil properties (named fluctuations). The spatial variability structure of measured properties can also be characterized through kriging methods based on the regionalised variable theory (Matheron, 1973). The Ordinary Kringing has been used in the present study to scale of fluctuation by variogram functions to be used in the RFEM analyses has been used in the present study to identify the variation range of δ in the RFEM analyses. However, δ has also been assessed by fitting the empirical correlation coefficient (ρe) (the method commonly known as the Vanmarcke graphical method) with a theoretical correlation coefficient(ρt).The empirical method involves a single δ estimation for CPT sounding. However, along the length of the profile,the soil may change(with differentISBTvalues).As a result,the scale values could be misvalued.Then,the primary objective of the current study was to check the variability inISBTclassified as zone 3. Therefore, soundings were selected from the base, where profiles overlapped the investigated zone,and δ was estimated for average values ofqcandfsalong the depth. The values of δ are significantly different when they are derived from the semivariogram range and the empirical method. Nonetheless, Vanmarcke’s approach,based on the random field theory,has been widely used in geotechnical design through fully probabilistic approaches.Thus,the semivariogram method is used to validate the other two methods and calculate the variation range of δ.

4.1. Estimation of the scale of fluctuation

Vanmarcke (1983) proposed to fit the theoretical correlation coefficient ρtto the empirical correlation coefficient ρeresulting from soundings.The ρeof the autocorrelation function was adopted as suggested in Cami et al. (2020):

whereXiis the point values;jandkdescribe the number of pairs points at a distance τj; and μeand σ2eare the mean value and the variance resulting from the soundings, respectively.

The analyses demonstrated that depending on the theoretical autocorrelation function, the scale of fluctuation values were measured from the first 10 m of the 11 selected CPT profiles. It is equal to δ=0.81 m for the Markov function(ρMt),and δ=0.65 m for the Gauss (square exponential) function.

where τ is the distance between two points in space.

The fitting results are shown in Fig. 3.

4.2. Design variables derived from the CPT profiles

Besides the scale of fluctuation,soil strength parameters are also drawn from CPT profiles. One of them is the undrained shear strength (su): The method to calculatesufromqcis described in Lunne and Kleven (1981) based on the corrected total cone resistance (qt)， that is:

whereu2is the pore pressure measured behind the cone;anis the net area ratio;σzis the total overburden pressure at the elevation of the cone; andNktis the tip factor, which is an empirical quantity varying between 14 and 16 (Robertson and Cabal, 2015).Nkt= 15 is assumed in the current study.

Variable values ofsuwith depth were determined for the 11 soundings by assuming that the soil unit weight (γ) is also a variable quantity given by the following formula proposed by Bagi′nska(2016):

As can be seen in Fig.4a,sushows a clear trend with depth.The histogram of values (Fig. 4b) shows that it fits well to a lognormal distribution with a mean valuesu，mean= 68:07 kPa, and a standard deviationsu，st:dev= 14:42 kPa: Then, the COV ofsu，COV= 0:21 has been determined.

Fig. 3. Fitting theoretical autocorrelation coefficients to the empirical value derived from CPT soundings.

To fully illustrate the effect of the spatial variability of soil properties on the bearing capacity of the shallow foundation, the effective strength parameters (c′and φ′)were also estimated from CPT profiles. Several methods can be found in the literature for estimating the effective strength parameters of fine soils. The Technical University of Norway method (NTH method) has been widely used,especially in soft soils(Senneset and Janbu,1985).This method assumes that effective strength parameters of soil are determined mainly through the φ′,and the relationship betweenc′and φ′is inverse. The preceding variables will be assumed to be independent when modeling the soil in the RFEM code. For this reason, the equation method suggested in (Eslami et al.,2019) has been used for estimating the effective strength parameters. However, the authors determinedc′using the procedure proposed in Sorensen and Okkels(2013),issued in the Danish Standard(DS.415,1998). The relationship between effective cohesion and undrained shear strength is based on a comparison of the ultimate bearing capacity(qult)under undrained and drained conditions for the plate loading test on clay-based on Jacobsen (1970)’s work:

However, to obtain the φ′, a method derived from the assumption that the static cone is treated as a deep micro foundation(Eslami and Fellenious, 1997) can be used. Then φ′can be determined by approximating the Terzaghi equation according to the proposition shown in Eslami et al. (2019):

whereNc，Nq，Nγ are the strength coefficients;andB=0.0375 m is the size of the cone. The load (q) is taken as the soil weight (γ)multiplied by the cone depth (z).

The strength coefficients are calculated as follows:

where β is the angle of plastification. In this study, β = 0 is assumed.

The values of the estimatedc′and φ′are shown in Fig.5a and b.Bothparametersvary withdepthininversetrend.Asmaller valueofc′corresponds to a larger φ′,and asc′increases withdepth,φ′decreases.

Fig. 4. (a) Variability of su parameter along with depth, and (b) histogram of su.

Accordingly, the strength parameters obtained from the soundings were fitted to the theoretical distribution functions, i.e.lognormal for the cohesion and bounded firor the internal friction angle, as widely adopted in the previous studies (Fenton and Griffiths, 2008; Pieczy′nska-Koz?owska et al., 2015). The strength parameter values are presented in Table 1.

The random field model assumes a mean value and a standard deviation estimated for a certain layer thickness (in the current study, a 10 m layer is assumed). A natural phenomenon in soils is the changes of parameters with depth.This effect is related,among other phenomena,to sedimentation.In addition to the changes,the soil properties show residual changes, which are used to estimate the spatial variability by the stationary random fields in RFEM models. The formulae transforming the CPT measurement profiles into design variables (su(Eq:(5),c′(Eq:(7) and φ′(Eq:(8)) do not reduce the trend. This must be determined (as shown in Figs. 4a and 5a) and then removed to estimate the scale of fluctuation from the residuals.

Fig.5. The variation of effective strength parameters with depth:(a)Cohesion(c′)and(b) Internal friction angle (φ′).

5. Numerical model

The RFEM was used for calculating the bearing capacity of a shallow foundation, with the width varying from 0.6 to 6 m. By estimating the strength parameters and their spatial variability,the LAS can be applied, and the RFEM simulations can be carried out.The random field theory could be applied to model the variability structure of the soil. The individual parameters described by probability distributions were modelled as independent random fields. In undrained conditions, onlysuwas modelled as a random variable. In drained conditions, two independent variables were used to calculate the strip foundation bearing capacity:c′and φ′.The classical 2D FEM based on rbear2d algorithm by Fenton and Griffiths (2008) modified by Pieczy′nska (2012) was used to perform the calculations. In the following study, a mesh was adopted (see Fig. 6) with 88 × 36 elements of 0.25 m × 0.25 m,corresponding to an overall dimension of 22 m×9 m.The size was determined using a single grid model to analyse the spatial variation of the bearing capacity under different foundation dimensions.The investigated width dimensions areB=0.5 m,1 m,1.5 m,2 m,2.5 m, 3 m, 4 m, 5 m, and 6 m. The foundation depth (D) and the unit weight γ of soil under the foundation are assumed deterministic. All parameters used in these numerical analyses are summarised in Table 1.

Table 1The strength and other design parameters in this study.

Fig. 6. FEM model used in these analyses.

Applying the adopted mesh size in FEM,a comparison between the deterministic FEM values and the analytical bearing capacity estimated through numerical simulations under conditions was performed.

6. Results and discussion

6.1. Deterministic vs. RFEM analyses

At the beginning of the calculation,analytical calculations of the bearing capacity for different foundation dimensions were conducted based on the mean values of the soil strength parameters.Subsequently, the deterministic FEM calculations were calibrated to represent the analytical case as closely as possible. The results are shown in Fig.7.The coloured columns(yellow,blue,and green)illustrate the quantitative contribution of each part of Terzaghi’s formula for the boundary resistance of the soil.The yellow columns are related to the cohesion (133.33 kPa), the blue columns are related to the foundation embedment (201 kPa), and the green columns are related to the soil unit weight (from 40.54 kPa forB=1 m to 486.43 kPa forB=6 m depending on foundation width).The red lines correspond to the summed analytical values under drained and undrained conditions.The black crosses are the values estimated deterministically by the FEM.

It is noted that the values are not dependent on the foundation dimensions for bearing capacity calculation under undrained conditions. Conversely, the values related to drained conditions increase linearly with the foundation’s width.

Performing calculations for 5000 realisations ensured satisfactory convergence of mean values and standard deviation of bearing capacity for comparison. Increasing the number of realisations in Monte Carlo simulations improved the stability of the obtained results.Assuming 5000 realisations allows studying the changes at different levels. It is possible to analyse mean values, standard deviations and distributions. The optimal number of realisation in RFEM analyses has been discussed in the previous papers(Pieczy′nska, 2012; Pieczy′nska-Koz?owska et al., 2015; among others).

Comparing the two calculation conditions for the case of Gaussian autocorrelation function, Fig. 8a and b shows that the stochastic mean values are equal to the deterministic ones in both cases. An apparent difference is shown at extreme values. Under undrained conditions, when onlysuis considered as the design parameter, the random values of the bearing capacity decrease when the foundation width increases. Conversely, when the two variables (c′and φ′) are considered, the random values of the bearing capacity increase as the foundation width increases.Furthermore, as the vertical scale of fluctuation increases (0.65 m,2 m, and 6 m), the random values of the bearing capacity increase accordingly (Fig. 8, from left to right, from the top to the bottom).The black lines in Fig. 8 show that the ultimate bearing capacity(qult) varies slightly with δv.

Fig. 7. Comparison of analytical and deterministic values of the bearing capacities in strip foundation for different foundation width under drained and undrained conditions.

Comparing the mean values generated for each scale of fluctuation (Fig. 9), it can be noted that a larger scale of fluctuation generates a greater mean value of the bearing capacity.This effect is visible in both drained and undrained conditions. In undrained conditions, this effect changes when foundation width is greater than 5 m. Interestingly, in perfectly cohesive soil with increasing foundation dimension, the random value of the bearing capacity decreases.This effect is not observed in the case of soil described by effective strength parameters.

Furthermore,a complex increasing trend can be observed under drained conditions although all the δvvalues give almost the same mean value of bearing capacity.Therefore,whenBranges between 2.5 m and 4.5 m, the stochastic values are higher than the deterministic ones, the contrary whenBis larger than 4.5 m. This occurrence can be seen in Figs. 8-10 and 13. These small changes around the deterministic trend are about 3%-4%: they are numerical noise and not meaningful for geotechnical design. The same results can be observed for the undrained conditions.

Additionally, the analysis of the horizontal scale effect on the bearing capacity was performed. Due to the large distance among soil investigations in the Po River Plain area, the horizontal scale values were not obtained through the semivariogram method. It was set to δH= 1400 m in this study.

Concerning the adopted FEM computational model (Fig. 6),the horizontal scale of fluctuation was modelled as infinite. To have a complete overview of the horizontal and vertical variability effect, the variability of the bearing capacity has been investigated for the horizontal scale falling within the theoretical value δH= 10 m. A comparison of the results for both scales is presented in Fig. 10.

Fig.10. Comparison of horizontal scale effect on the bearing capacity at δv = 0:65 m under (a) undrained and (b) drained conditions.

The horizontal scale of fluctuation only affects the minimum values estimated for the perfectly cohesive case. The scale δh=10 m slightly narrows the spectrum of random bearing capacity values without any visible effect on its mean value.

When studying the soil parameters variation, the standard deviation (σqult) shall be considered. Comparing the two cases of drained and undrained, the dimensionless COV depends on both the mean value (μqult) and the mentioned standard deviation:

The relationship between the COV and the foundation size has been investigated, and the results are shown in Fig.11.

Fig.11. Comparison of the bearing capacity COV for different foundation widths using Gaussian correlation function with the horizontal scales δh = 1400 m under drained(black filled square)and undrained (red open square)conditions: (a) A constant vertical scale of fluctuation (δv = 0.65 m) is considered; and (b) Different vertiacal scales of fluctuation.

Fig. 11a and b can be drawn an inverse relationship between COV and foundation width that is not linear but exponential.Nonetheless, the two equations show different determination coefficients(R2),and the points for the drained case are sparser than the undrained case. Thus, based on the trend line, a high fit to the power function fading with increasing foundation dimension is seen. The drained case shows more variability for foundations smaller than 1.5 m compared to width greater than 1.5 m.Here the greater variability is also maintained in the undrained case. By comparing the maximum and minimum values of COV(Fig.11a),it can be found that the difference for drained case (⊿COV(drained)) is less than the undrained case (⊿COV(undrained)), about 8%-11%.

The difference between the maximum value of COV (forB= 0.5 m) and its minimum value (forB= 6 m) is related to the scale of fluctuation.As the vertical scale of fluctuation increases,the COV increases. This effect is observed for each foundation size(Fig. 11b). Interestingly, the difference between maximum and minimum values is almost the same for each δv. In the undrained case,the values are 8%for δv=0.65 m,8.6%for δv=2 m,and 6.6%for δv= 6 m. The difference between COV values is greater in the drained case: 11% for the δv= 0.65 m and δv= 6 m, and 12% for δv=2 m.The δv=2 m shows the highest variability in both cases.

In the undrained case,whenB=0.5 m,the COV for δv=0.65 m is 70% less than that for δv= 6 m. WhenB= 6 m, the difference between the scales reaches 43%.The deviation between the values is slightly smaller in the drained case:it varies from 68%forB=0.5 m to 36% forB= 6 m, respectively.

The horizontal scale of fluctuation does not influence the mean value of the bearing capacity. However, the analysis of the COV allows observing an influence as shown in Fig.12.

In Fig. 12, divergence can be observed as the foundation dimension increases. The influence of two shapes of the autocorrelation functions is illustrated in Fig.13. The bearing capacity results for the Gauss and Markov autocorrelation functions show no significant differences between the two shapes.

The results for δv= 2 m in drained and undrained cases using two autocorrelation functions are shown in Fig. 13. The results show no significant differences. The Markov function generates slightly larger differences in mean values of bearing capacity depending on the horizontal scale of fluctuation. In the perfectly cohesive case, the bearing capacity values are slightly smaller for δh= 10 m,while in thec-φ case for δh= 10 m,a slightly larger mean is calculated.

Considering the effect of the autocorrelation function,it is most evident in the perfectly cohesive case.Here,the use of the Markov function generates a rapid reduction in the range of maximum values while the observed minimum values are increased.

Comparing the values of the COV for both autocorrelation functions in Fig.14, larger differences among curves can be noted related to different δhin Gaussian function case. For the case of δh= 10 m under undrained conditions, theCOVqultvalues are almost the same regardless of the autocorrelation function. An interesting phenomenon can be observed when the Gaussian function is used, i.e. COV values for small foundation widths are almost the same under drained conditions. A similar effect is observed for large foundation widths in the case of Markov function.

Fig. 8. The values of random bearing capacity for foundations with different dimensions, at δh = 1400 m under (a) undrained and (b) drained conditions.

Fig. 9. Comparison of the mean values of bearing capacity at δh = 1400 m for foundations with different dimensions under (a) undrained and (b) drained conditions.

Fig. 12. Comparison of COVqult under drained and undrained conditions (δv =0:65 m).

Fig. 13. Comparison of RFEM bearing capacity values with the applied Gaussian (a, b) and Markov (c, d) autocorrelation functions under drained (b, d) and undrained (a, c)conditions, with the assumed vertical scale of fluctuation δv = 2 m. Each figure considers two horizontal scales of fluctuation: δh = 10 m (dots and crosses) and δh = 1400 m(solid and dotted lines). The RFEM mean values (blank dots and black solid lines), the RFEM minimum values (black dots and dotted lines),and the RFEM maximum values (black crosses and grey solid lines) are compared with the deterministic FEM calculations (yellow solid lines).

7. Conclusions

In this paper, the strip foundation design in terms of bearing capacity has been investigated by considering the subsoil spatial variability structure from alluvial mixtures of silt, clay, and sand located in the Po River Plain. Based on the popular calculation method of RFEM, the paper shows a comprehensive method for obtaining soil modeling parameters using stationary random fields directly from CPT profiles. The analyses show that it is possible to estimate all necessary parameters such as strength parameters(under drained and undrained conditions) and parameters describing the related autocorrelation functions with a database of CPT measurements. The most relevant results are listed below:

(1) In drained and undrained conditions, the ultimate resistances resemble each other, i.e. the range of resistances between the minimum and maximum values increase with the scale of fluctuation.

(2) As the foundation dimension increases,the value of the scale of fluctuation has a more pronounced effect on the mean value of the bearing capacity. In the drained condition, a small scale of fluctuation (0.65 m) significantly reduces the mean value.In the undrained case,the largest differences of the mean values concerning the deterministic values were observed for the scale of 6 m.

(3) The horizontal value of the scale of fluctuation showed no significant effect in the current task.

(4) The use of different autocorrelation functions pointed out that they have a low impact on the mean value of the bearing capacity. However, the Gaussian function generated a maximum value of bearing capacity higher than the Markov function in undrained conditions of a few tens of kilopascals.The results from RFEM analyses gave additional information on the COV and then on the probability of failure of the bearing capacity:

(5) The shape of the autocorrelation function induces some differences in final values of COV.

Fig.14. The COV of bearing capacity versus the foundation size for the two shapes of the autocorrelation functions: (a) Gaussian and (b) Marko (δv = 2 m).

(6) The bearing capacities calculated under drained and undrained conditions show a small difference in COV values,which decrease as the foundation size increases.

The present study suggested that it is a severe simplification to disregard the spatial variability of soils considering only the mean values of the design variables instead of considering the parameters such as the scale of fluctuation and the COV. Comparing the mean valueof the stochasticbearing capacity withthe deterministic one,no significant differences can be recognised. The COV analysis shows differences due to the scale of fluctuation.It causes an increase in the probability of failure, which reduces with the increasing scale of fluctuation and rises up with increasing the foundation width. The worst design case in this study is the foundation size ofB=6 m and the scale of fluctuation equals 0.65 m.None of the results presented herein can be drawn through a deterministic analysis but the mean values.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.