The Effect of the Selection of Three-Dimensional Random Numerical Soil Models on Strip Foundation Settlements

This paper delivers a probabilistic attempt to prove that the selection of a random three-dimensional finite element (FE) model of a subsoil affects the computed settlements. Parametric analysis of a random soil block is conducted, assuming a variable subsoil Young’s modulus in particular finite elements. The modulus is represented by a random field or different-sized sets of random variables; in both cases, the same truncated Gaussian model is assumed. Mean values and standard deviations of random soil settlement are estimated by a Monte Carlo simulation procedure. With regard to the adopted FE model, the estimated settlement mean values do not vary significantly, but standard deviations do strongly. Similarities also appear in the diagrams of random field correlation length versus settlement standard deviation and the diagrams displaying a total number of model random variables versus settlement standard deviation. Thus, relevant single random variable models represent the random field approach well with regard to settlement parameter estimation. This remark is verified upon a settlement analysis of a three-dimensional FE model of a hypothetical strip foundation. Following the preliminary model observations, various probabilistic geotechnical analyses may be supported, e.g., continuous footing design, slope stability analysis, and foundation reliability assessment.


Introduction
Soil medium displays an inherent material parameter variability. Numerical modelling of actual spatial soil variability is a considerable task, even when extensive in situ tests are carried out at the construction site. Hence, the only way to relevantly capture the uncertainty of subsoils leads to probabilistic methods. A strong recommendation for random methodology is now active, and the standards, e.g., Eurocodes, follow this. A direct introduction to the probabilistic approach in standard engineering computations is a present-day requirement and a future development perspective in the field [1][2][3].
Probabilistic methods are extensive due to the enormous development of computational power, introducing advanced finite element analysis (FEA) software to solve numerical problems of the subsoil media. Both the bearing capacity of foundations in their limit states (ultimate limit state-ULS) and the assessment of settlements (serviceability limit state-SLS) employ probabilistic analytical means. An extensive literature review in the field is presented in, e.g., [4], including the application of structural reliability analysis (SRA) and the stochastic finite element method (SFEM) to geotechnical problems. However, the classical Monte Carlo (MC) simulation method is the most frequent, as it is easy in numerical implementation, widely available in the software market, and simple in the interpretation of its results. Many MC variance reduction techniques are developed on a high scale, e.g., stratified sampling and Latin hypercube sampling [5]. The MC simulation method is fundamental in reliability-based robust geotechnical design (RGD) [6].

Random Spatial FE Models of the Subsoil
The issue of the numerical uncertainty definition of the subsoil is analysed in various FE models of a spatial subsoil block presented in Figure 1. The analysis is performed in ZSoil ® software [31][32][33] coupled with Python ® programming language [34] to aid the series of repetitive parametric computations. The representative random subsoil volume area is modelled as a block with dimensions of 10 × 10 × 10 m (all finite elements have dimensions equal to 1 × 1 × 1 m), distinguished in pink in Figure 1.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 3 of 20 the mechanical response variability of the subsoil, impossible to find in the deterministic analysis.

Random Spatial FE Models of the Subsoil
The issue of the numerical uncertainty definition of the subsoil is analysed in various FE models of a spatial subsoil block presented in Figure 1. The analysis is performed in ZSoil ® software [31][32][33] coupled with Python ® programming language [34] to aid the series of repetitive parametric computations. The representative random subsoil volume area is modelled as a block with dimensions of × × 10 10 10 m (all finite elements have dimensions equal to × × 1 1 1 m ), distinguished in pink in Figure 1. In the task, the Young's modulus of the subsoil E is considered the only random variable, with a Gaussian probability density function (PDF) of the mean value = 28,000 kPa and standard deviation σ = 7000 kPa E . While the coefficient of variation is , there is a 0.00317% probability that the generated i E value will be negative; The boundary conditions at the interface with the surrounding subsoil are reflected by an additional layer of elements whose Young's modulus is E bc = 560 kPa and Poisson's ratio is ν bc = 0.25, distinguished in green in Figure 1.
These parameters are established based on a preliminary numerical test. First, a simple compression test of the boundary element stratum is performed, assuming Young's modulus E bc,t = 60, 000 kPa and Poisson's ratio ν bc,t = 0, so that the boundary element soil stratum does not contribute to lateral displacements of the main subsoil block. The lateral displacement of the boundary layer δ bc = 0.000803 m is determined. Subsequently, this value is substituted as the target lateral displacement of the global subsoil block FE model, with the lateral contribution (ν bc = 0.25) of the boundary layer assumed. This allows determining the required effective stiffness parameter of the boundary layer as E bc = 560 kPa, ensuring proper subsoil deformability in all directions.
The total number of finite elements is equal to n el,tot = 12 × 12 × 11 = 1584, whereas only n el = 10 × 10 × 10 = 1000 constitute the regarded subsoil. The subsoil loading q = 170 kPa is situated on the soil block indirectly, represented by a flexible metal sheet (plate) with a width w = 6.0 m, thickness t = 0.01 m, and material parameters corresponding to S235 steel (Young's modulus E s = 210 MPa, Poisson's ratio ν s = 0.30). The subsoil is governed by the Coulomb-Mohr material model.
In the task, the Young's modulus of the subsoil E is considered the only random variable, with a Gaussian probability density function (PDF) of the mean value µ E = 28, 000 kPa and standard deviation σ E = 7000 kPa. While the coefficient of variation is v = 0.25, there is a 0.00317% probability that the generated E i value will be negative; therefore, a truncation parameter of t = 4.0 (E i ∈ µ E ± tσ E ) is introduced, and hence no negative Young's moduli appear in the subsequent generation.
Note that a relatively large Young's modulus variability of the subsoil was purposefully assumed to introduce a significant impact on the scatter of the numerical results of the settlements. Thus, the random variable distribution is not realistic. However, in most design cases in real-life engineering, the data on subsoil variability acquired from the in situ surveys are limited. Hence, suggestions in design codes and literature guidelines vicariously reflect the random parameters, adopting their high variability. As an example, variability v = 0.3 is incorporated in the reliability assessment of foundations resting on piles [35] or in stability analysis of a coastal cliff [36].
Both cohesion c = 10 kPa and internal friction angle φ = 30 • are assumed deterministic but, in fact, can be assumed random. The algorithms proposed in the paper do not limit the number of material parameters to be considered random variables. However, to compare the influence of the adopted FE models on the variability of the computed settlements, Young's modulus E is assumed the only random variable.
Various random FE models of the subsoil block are assembled in the ZSoil ® software, representing the global randomness of the FE model in different ways. Upon a given number of crude Monte Carlo simulations, all FE models examine the positive displacement of the central node of the upper (loaded) subsoil block surface (under the centre of the metal sheet transmitting the uniformly distributed load). The mean values and standard deviations of the central node displacement are investigated.

The Subsoil Material Parameters Described by Three-Dimensional Random Fields
The geotechnical parameter uncertainty of soils may be relevantly reflected by threedimensional random fields. This action requires advanced numerical tools and theoretical background in advanced probabilistic methods. Although random field generation is an advanced issue, implemented in relevant software [11,37], the random field methodology is still rare in real-life engineering applications. The only software of wide present use includes simple generators of random fields to be directly entered as input data into standard FEM software.
In this paper, such a simple computational path is also executed. Dedicated software is used to generate the FE input-ready Gaussian two-dimensional and three-dimensional random fields, and the software is described in detail in [38][39][40]. The software conducts random variable generation employing a conditional acceptance and rejection concept.
A conditional function f s (X u /X k ) is applied to directly generate an unknown vector X u , assuming that the field part defined by a vector X k is already generated and is known: where K c is a conditional correlation matrix, X c is a conditional vector of mean values, and m stands for the dimension of the correlation matrix K.
The most important aspect of random field modelling is the correlation function adoption, dedicated to a given problem. Various correlation function definitions are found in the literature; usually, the function is assumed a priori. Comprehensive studies on parameter determination based on targeted research results can be found (see, e.g., [21]). Random fields of standard engineering materials (e.g., soils, concretes, composites) make use of the exponential correlation function related to the first-order autoregressive function or a Markov process model. In a three-dimensional problem case, it has the following definition [41]: where d x , d y , and d z are damping parameters describing the correlation decay, and |x|, |y|, and |z| are the distances between the points of the random field, along respective axes. Due to simple engineering interpretation, the correlation function (2) has been used in many cases (see, e.g., [17][18][19]42]).
In the calculations, diverse random field types were generated, differing in the value of the adopted damping parameter d ≡ d    The random field applied in the FE model of the subsoil block is based on the previous truncated Gaussian variable of Young's modulus E (µ E = 28, 000 kPa, σ E = 7000 kPa, t = 4). A number of 200 three-dimensional Gaussian fields of n re = 1000 points are generated, corresponding to each of the five correlation decay parameter values, directly implemented in the FE model. The number is limited to 200 samples only, as the time required to generate the fields, to create FE models incorporating these data, to conduct the numerical calculations, and to analyse the results is extremely significant. A graphical example of a generated highly correlated random field (d = 0.01 m −1 ) is shown in Figure 3.   The assessed displacements δ i RF (i = 1, . . . , 200) of a highly correlated random field (d = 0.01 m −1 ) are shown in Figure 4, along with the convergence diagrams for both mean value µ δ RF and standard deviation σ δ RF of model response. The figure detects two positive vertical axes-standard deviation is presented on a positive upward axis, while displacements refer to a positive downward axis, the latter reflecting the settlement direction. This assumption is valid in subsequent similar figures. The assessed displacements δ i RF (i = 1, … , 200) of a highly correlated random field are shown in Figure 4, along with the convergence diagrams for both mean value δ μ RF and standard deviation δ σ RF of model response. The figure detects two positive vertical axes-standard deviation is presented on a positive upward axis, while displacements refer to a positive downward axis, the latter reflecting the settlement direction. This assumption is valid in subsequent similar figures. The data in Figure 4 are representative of all five cases of different correlation decay parameters, due to the similar mechanical response of the soil block in the case of the other four damping parameters d . The MC simulation convergence criterion is adopted: convergence holds if the mean value relative error is lower than μ Δ = 0.01% and the standard deviation relative error is lower than σ Δ = 1% in the course of a 20-sample computation. In all five cases, the adopted criterion was fulfilled. However, higher MC simulation precision requires computing a larger number of samples. Computation of the extreme central node settlements δ RF , their mean values δ μ RF , and standard deviations δ σ RF for all considered FE models are shown in Table 1. Note that Figure 4 displays peaks in the settlement diagram, slightly disturbing the response parameter convergence. In a highly correlated random field (  Table 1). The peaks stem from the FE models with a low Young's modulus taken as the basis for random field generation. Due to high The data in Figure 4 are representative of all five cases of different correlation decay parameters, due to the similar mechanical response of the soil block in the case of the other four damping parameters d. The MC simulation convergence criterion is adopted: convergence holds if the mean value relative error is lower than ∆ µ = 0.01% and the standard deviation relative error is lower than ∆ σ = 1% in the course of a 20-sample computation. In all five cases, the adopted criterion was fulfilled. However, higher MC simulation precision requires computing a larger number of samples. Computation of the extreme central node settlements δ RF , their mean values µ δ RF , and standard deviations σ δ RF for all considered FE models are shown in Table 1. Note that Figure 4 displays peaks in the settlement diagram, slightly disturbing the response parameter convergence. In a highly correlated random field (d = 0.01 m −1 ), the extreme value δ d=0.01 RF = 0.1825 m is reached ( Table 1). The peaks stem from the FE models with a low Young's modulus taken as the basis for random field generation. Due to high correlation, an initially low modulus propagates to a significant number of finite elements in the block volume, thus generating a weak subsoil, prone to excessive displacements.
The peaks vanish with correlation decay increase, as the propagation of the initial low modulus quickly fades, not allowing for an extensive material weakening of the FE subsoil model. In a slightly correlated random field (d = 5.0 m −1 ), the peaks are not observable, so the extreme displacement δ d=5.0 RF = 0.1085 m does not come too far from the mean value estimated: in this case, µ d=5.0 RF = 0.1048 m ( Table 1). The computational results (Table 1) can be summarized as a closed-form relation between the correlation decay parameter and standard deviation variability of the analysed settlement. The logarithmic approximation σ δ RF (d) = −0.00201831 ln(d) + 0.00363733 is shown in Figure 5. As the settlement mean value remains at a specified level throughout the computations of distinct sets (Table 1), it is not considered here. tween the correlation decay parameter and standard deviation variability of the analysed settlement. The logarithmic approximation δ σ = − + RF ( ) 0.00201831ln( ) 0.00363733 d d i shown in Figure 5. As the settlement mean value remains at a specified level throughou the computations of distinct sets (Table 1), it is not considered here.  Figure 5 shows that the standard deviation change of the settlement in a random field approach to FE models may be represented by appropriate damping parameters d in Equation (2). However, random field generation is difficult to perform and interpret.
The field definition also requires an advanced calibration of the FE model, while the type of correlation function has to be taken in advance with its parameters, e.g., correlation lengths or the damping parameter.

The Monte Carlo FE Material Model
There are obstacles in the application of random fields in engineering problems hence, an alternative variant is proposed. The generated Young's moduli are linked with specific finite elements of the numerical model. Merging the elements into sub-block makes it possible to produce the global subsoil response. A variation of the method re sults, e.g., foundation settlements, is reflected by a direct Monte Carlo simulation. Thus the paper proposes a Monte Carlo FE material model to directly represent the uncertain character of the subsoil volume. Parametric analysis of four FE models illustrates the com putational course, and the generated Young's moduli are linked with the variable sub block structure.  Figure 5 shows that the standard deviation change of the settlement in a random field approach to FE models may be represented by appropriate damping parameters d in Equation (2). However, random field generation is difficult to perform and interpret.

The Homogenous FE Model Described with One Random Variable
The field definition also requires an advanced calibration of the FE model, while the type of correlation function has to be taken in advance with its parameters, e.g., correlation lengths or the damping parameter.

The Monte Carlo FE Material Model
There are obstacles in the application of random fields in engineering problems; hence, an alternative variant is proposed. The generated Young's moduli are linked with specific finite elements of the numerical model. Merging the elements into sub-blocks makes it possible to produce the global subsoil response. A variation of the method results, e.g., foundation settlements, is reflected by a direct Monte Carlo simulation. Thus, the paper proposes a Monte Carlo FE material model to directly represent the uncertain character of the subsoil volume. Parametric analysis of four FE models illustrates the computational course, and the generated Young's moduli are linked with the variable sub-block structure.

The Homogenous FE Model Described with One Random Variable
The initial model assigns a single random variable n E = 1 (with the truncated normal distribution of Young's modulus of the soil: µ E = 28, 000 kPa, σ E = 7000 kPa, t = 4) to all n r = 1000 elements. Thus, the entire block stiffness is represented by a single randomly generated value. Such an FE model is straightforward in its numerical operation; hence, frequently introduced in probabilistic analysis.
The first stage of the analysis checks the convergence of the investigated displacement with 200 FE models computed for different E i random variable generations. The mean value µ 200 δ 1 = 0.1087 m and standard deviation σ 200 δ 1 = 0.0091 m are obtained here. In the second stage, the test was expanded to a total of 10,000 samples, bringing a sample mean µ 10000  Figure 6. The 200-sample test is separated in Figure 6a, and the complete 10,000-sample test is presented in Figure 6b. Appl. Sci. 2021, 11, x FOR PEER REVIEW 9 of 20  If a homogeneous FE soil model is incorporated, the dispersion is unclear. The dispersion decreases with a lower truncation parameter of Young's modulus distribution, e.g., t = 1.0; however, such an assumption is not justified by in situ soil studies. Soils in their natural state show significant variability, and the coefficient of variation may reach 0.4 (see [43][44][45]). Thus, larger values of the truncation parameter are recommended (t 1.0). Although the presented results (Figure 6a, 200 samples) are improper and only accidentally reliable (no possible peaks are observable), they serve as a reference point for further analysis.

The FE Models Described with Several Uncorrelated Random Variables
An alternative variant of the FE model is proposed to produce results of higher numerical quality. The model introduces a specified, relatively small number of random Young's moduli, assigned to smaller subsets of adjacent elements of the entire soil block.
First, a subdivision of the analysed subsoil block of n r = 1000 elements into eight equal sub-blocks was assumed (n E = 8), as presented in Figure 7. If a homogeneous FE soil model is incorporated, the dispersion is unclear. The dispersion decreases with a lower truncation parameter of Young's modulus distribution, e.g., = 1.0 t ; however, such an assumption is not justified by in situ soil studies. Soils in their natural state show significant variability, and the coefficient of variation may reach 0.4 (see [43][44][45]). Thus, larger values of the truncation parameter are recommended  ( 1.0) t . Although the presented results (Figure 6a, 200 samples) are improper and only accidentally reliable (no possible peaks are observable), they serve as a reference point for further analysis.  As a result, each sub-block consists of × × = 5 5 5 125 elements. For each sub-block, a distinct random Young's modulus E is independently generated; in turn, the FE model is described by = 1000 / 125 8 different Young's moduli. The analysis is limited to 200 FE models only. Satisfactory convergence of the solution is reached under the conditions of μ Δ , and σ Δ is set in the random field case ( Figure 8). As a result, each sub-block consists of 5 × 5 × 5 = 125 elements. For each sub-block, a distinct random Young's modulus E is independently generated; in turn, the FE model is described by 1000/125 = 8 different Young's moduli. The analysis is limited to 200 FE models only. Satisfactory convergence of the solution is reached under the conditions of ∆ µ , and ∆ σ is set in the random field case (Figure 8). Similar to the homogeneous subsoil block, the unjustified central node displacement peaks may also appear, following a large number of samples. However, contrary to the single variable model, the peaks occur only when small Young's moduli are generated for the majority of eight sub-blocks, and the probability of such a generation seems relatively low.
To further refine the analysis, another subdivision of the analysed subsoil block is introduced. Here, the 1000-element block is divided into 125 equal sub-blocks (Figure 9), each consisting of 2 × 2 × 2 = 8 elements. This means an input of 1000/8 = 125 different Young's moduli into the global FE model (n E = 125). reached. Similar to the homogeneous subsoil block, the unjustified central node displacement peaks may also appear, following a large number of samples. However, contrary to the single variable model, the peaks occur only when small Young's moduli are generated for the majority of eight sub-blocks, and the probability of such a generation seems relatively low.
To further refine the analysis, another subdivision of the analysed subsoil block is introduced. Here, the 1000-element block is divided into 125 equal sub-blocks (Figure 9   Although convergence holds for 200 MC samples, under the previous conditions, a broader scale analysis of 10,000 FE models is performed. The results and convergence of the crucial parameters of the random response are shown in Figure 11. Although convergence holds for 200 MC samples, under the previous conditions, a broader scale analysis of 10,000 FE models is performed. The results and convergence of the crucial parameters of the random response are shown in Figure 11.  Summing up, all essential results of the considered computational cases (Young's moduli modelled by sets of random variates) are presented in Table 2. With regard to FE models, the mean settlement is not affected by Young's modulus variation, but the standard deviation of the central node settlement decreases with the random variable domain rise. The extreme settlements are also listed in Table 2. It should be noted that if different sample sets are computed, the presented output is about to vary. This is due to extreme Summing up, all essential results of the considered computational cases (Young's moduli modelled by sets of random variates) are presented in Table 2. With regard to FE models, the mean settlement is not affected by Young's modulus variation, but the standard deviation of the central node settlement decreases with the random variable domain rise. The extreme settlements are also listed in Table 2. It should be noted that if different sample sets are computed, the presented output is about to vary. This is due to extreme deflections since a larger number of analysed samples increases the probability of unrealistic settlement peaks occurrence. Similar to random field modelling, an approximation of settlement standard deviation versus the number of the assumed random Young's modulus variables is performed ( Figure 13). An assumption is made on the horizontal axis-a cubic root of the number of variables c E = 3 √ n E (0 < c E < 10) is adopted, not the direct number n E (0 < n E < 1000).
The logarithmic function is approximated σ δ n E = −0.00339796 ln(c E ) + 0.00783447 (see Figure 13). Note that the diagram resembles the one in Figure 5 in the analysis of FE models with the random field modelling of Young's modulus E. The similarity between the diagrams in Figures 5 and 13 is not accidental. An apparent relationship is observable between the random field regions with similar stiffness moduli E defined by the correlation decay parameter d (Figure 3) and the number of sub-blocks E n that divide the entire FE model volume according to a given mode ( Figure   9). Although it is impossible to analytically define these relationships, they are visible in terms of extreme case analysis. When defining the total block volume by a single random value ( = 1  The similarity between the diagrams in Figures 5 and 13 is not accidental. An apparent relationship is observable between the random field regions with similar stiffness moduli E defined by the correlation decay parameter d (Figure 3) and the number of sub-blocks n E that divide the entire FE model volume according to a given mode (Figure 9). Although it is impossible to analytically define these relationships, they are visible in terms of extreme case analysis. When defining the total block volume by a single random value (n E = 1), ideal homogenization of the material emerges, and it resembles the random field with infinitely strong correlation, d → 0 . Both simulations imply a large result scatter and significant standard deviation of the investigated settlement. While the number of sub-blocks is equal to the number of finite elements (n E = 1000), the FE model resembles the white noise random field of extremely weak correlation (d 5). It is evident that the cases corresponding to the real-life soil conditions are included between these extremes. Since random field computations are difficult, they can be represented by accurate block division into sub-volumes of single Young's moduli.
Note that the determined logarithmic approximations (Figures 5 and 13) can be related to the soil model presented in the paper only. However, it is supposed that the observed relation can be estimated in any other FE model of a geotechnical structure resting on random subsoil or the FE models of slopes of a complex multi-layered geological structure with uncertain parameters, e.g., prone to stability loss.
The optimal number of singular variables to be randomly generated in a simple approach n E can be devised. The data collected in Figure 13 lead to an assumption that the refined FE model with n E = 125 variables may be regarded as the finest one. If n E = 8 is assumed, the generated FE model becomes too similar to a model defined by a homogeneous field. The results show unjustified large jumps, so to reach the convergence of the MC solution, too many samples are required, discarding the analysis from real-life engineering problems. For n E = 125, the results are practically independent of the number of FE models applied in the computational series, and the result scatter is realistic. It balances the complexity of model generation and output reliability well. Nevertheless, accepting n E = 125 is arbitrary and strongly related to the analysed case. However, the diagram in Figure 13 shows that for such a value, the variability range of FE model parameters is small enough to produce engineering-sound solutions, decisive when investigating real-life problems.
In a general case, for any FE model, the number of sub-blocks with the same material parameters n E should be specified. The easiest way to define this number is based on the percentage of the number of sub-blocks p E (0.0 < p E < 1.0) to the total number of the finite elements in the model n el : For example, in the considered FE model of the block involving n el = 1000 finite elements, the optimal number of different random values generated for the task should be n E = 0.125 · 1000 = 125 (p E = 0.125). In turn, the entire volume is divided into 125 subblocks, each including 1000/125 = 8 finite elements. In the case of an irregular mesh of finite elements, specific random variables are linked with relevant volumes to represent their contribution to the entire modelled structure. Such computations, however, require additional, pre-generated external procedures assigning a certain number of Young's moduli to precisely set sub-volumes, matching the FE model syntax. A regular mesh of finite elements, usually applicable in subsoil models, corresponds to the presented approach.
Unfortunately, no clear distinction arises on the FE model definition, and the model adoption is intuitive. Despite this, such a simplified random variable approach is dedicated to advanced geotechnical engineering design to replace the three-dimensional random field analysis.

FE Model of a Continuous Strip Foundation
The advanced computational three-dimensional FE model of a perfectly continuous foundation resting on a three-dimensional random subsoil verifies the prior FE solution of a subsoil block. A sufficiently long segment of the foundation is analysed to provide numerical results of a displacement of the central node of the foundation longitudinal cross-section free of any interference of the boundary elements. The considered FE model of a segment of the strip footing foundation of dimensions 2.0 × 0.5 × 20.0 [m], resting on a spatially random subsoil, is presented in Figure 14. × × , resting o a spatially random subsoil, is presented in Figure 14. A volume of × × 18.0 6.0 6.0 [m] i generated in the ZSoil ® software, and the outer elements are properly constrained to ac curately reflect the mechanical response of the entire subsoil under the foundation. Th subsoil under the foundation level is only modelled, and the soil layers lying above th foundation-subsoil interface are replaced by a uniformly distributed load = 18.0 kPa q situated on the upper surface of the subsoil FE model. The dead load of the foundatio wall is represented by a linear load = 1.0 MN P per FE length. The wall-foundation con nection makes the deflections compatible; the rotations do not follow the rule.  Probabilistic analysis requires multiple repetitions of similar computations; it is a significantly time-consuming process, so assuming an optimal refinement of the finite element mesh is necessary. Division of the subsoil into a given number of elements should ensure the desired result accuracy, also providing a relatively short calculation time of a single case. Therefore, the subsoil volume under the foundation is divided into several sub-volumes, variable in the degree of finite element mesh condensation ( Figure 14). Moreover, the approach to the subsoil volume definition is also diverse, as the subsoil volume (10.0 × 6.0 × 6.0 [m]) directly under the strip foundation is assumed random (marked in yellow in Figure 14), while the remaining volume is modelled deterministic (marked in pink in Figure 14). To facilitate the description of the random soil sub-volume with random fields, a regular mesh of three-dimensional continual eight-node elements (Q4 continuum 2D) with dimensions 0.5 × 0.75 × 0.75 [m] is adopted. This brings a total number of 3200 finite elements to constitute the random sub-volume.
In the random sub-volume of the soil, Young's modulus parameters are assumed the same as in the case of a soil cube (Chapter 2). In the deterministic part, the parameter was set to its mean value E = µ E = 28 MPa.
The concrete of the foundation is considered elastic, and its parameters correspond to the C30/37 concrete (with Young's modulus of E = 33 GPa).
The soil block analysis shows an optimal finite element domain of diverse Young's moduli equal to 12.5% (p E = 0.125, Equation (3)) of the total FE domain. Here, the Monte Carlo results converge fast, and their dispersion is moderate. The random domain of the strip foundation subsoil model (3200 elements, Figure 14) is split into 400 smaller subvolumes, each of them made of eight finite elements (400/3200 · 100% = 12.5%). While this generated random variable is applied to each sub-volume, 400 separate Young's modulus random variables are incorporated in a single model, n E = 400.
The computations are performed for 200 generated foundation models. To assess solution convergence, parameter variation is examined on the central point settlement of the foundation, neglecting the boundary condition impact (Figure 15). Based on the computations, variable settlement shows its mean value µ    To provide a comparison, two additional computational series are conducted: with 50 models where the entire type 1 volume is defined by a single Young's modulus variable (full homogenization) and 50 models where 3200 variables are used (white noise-like simulation). The computational results are shown in Table 3.  To provide a comparison, two additional computational series are conducted: with 50 models where the entire type 1 volume is defined by a single Young's modulus variable (full homogenization) and 50 models where 3200 variables are used (white noise-like simulation). The computational results are shown in Table 3. Table 3. Strip foundation settlement parameters; the Monte Carlo FE material model.

Number of Different Young
Moduli n E The presented results show that all the conclusions drawn for the initial subsoil block model, presented in Chapter 2, are confirmed. The extreme-value simulations for n E = 1 and n E = 3200 form an envelope of the possible task solutions. If the number of variables to capture subsoil uncertainty is chosen according to Equation (3), numerical results obtained by the FE model (computed settlements, their mean value, and standard deviation) are more reliable from the engineering viewpoint. Moreover, balance is achieved between modelling accuracy of the subsoil-foundation system, the total time cost of numerical computations, and the quality of the output data.

Conclusions
Based on the probabilistic analysis presented in the paper, the following remarks and conclusions are made:

•
The proposed Monte Carlo simulation-based material variability definition of an FE model allows for an alternative advanced analysis, competitive with the use of three-dimensional random fields; • There is a possibility to create specific FE models with a relatively regular FE mesh to directly implement the randomly generated data from simple pre-generated external procedures, supporting the presented approach; • It is optimal to sub-divide the subsoil domain into smaller sub-volumes and to sequentially represent them by unique, independently generated random geotechnical parameters; • Although the subsoil parameters are assumed identical Gaussian variables (the same µ E and σ E are used to generate their parameters in subsequent computations), the global variability description of the FE model is defined the other way: it significantly affects the standard deviation scatter of the foundation settlement, while the mean values are unchanged; • From an engineering standpoint, a high result scatter in probabilistic settlement response is unacceptable; thus, modelling the entire random subsoil by a single variable (even with its parameters supported by real-life surveys) is improper; • Unjustified result dispersion may be reduced, e.g., by introducing truncated (bounded) input variables; however, it is contrary to the natural high variability of soil parameters, displayed by real-life in situ measurements; • The proposed number of unique random geotechnical parameters is about 12.5% of the total number of FE elements of the model to balance modelling accuracy with the numerical quality of the solutions; • The method to assign material parameters to finite elements is considered general, and it allows covering subsoil parameters other than Young's modulus, e.g., yield parameters, cohesion, internal friction angles, and porosity; however, to analyse the impact of all decisive parameters of the adopted soil model, an appropriate parametric analysis should be provided.
This paper finds it advantageous to incorporate relatively simple and regular FE models, with their randomness attributed to a set of randomly generated parameters assigned to adjacent finite elements. An approximate optimal ratio of these parameters to the total number of finite elements is proposed. Hence, simplified, engineer-oriented modelling produces computationally efficient, engineering-sound solutions.
Funding: The APC was funded by the Faculty of Civil and Environmental Engineering, Gdańsk University of Technology, Narutowicza 11/12, Gdańsk 80-233, Poland.