Numerical Analysis of Shallow Foundations with Varying Loading and Soil Conditions

: The load–deformation relationship under the footing is essential for foundation design. Shallow foundations are subjected to changes in hydrological conditions such as rainfall and drought, affecting their saturation level and conditions. The actual load–settlement response for design and reconstructions is determined experimentally, numerically, or utilizing both approaches. Ssettlement computation is performed through large-scale physical modeling or extensive laboratory testing. It is expensive, labor intensive, and time consuming. This study is carried out to determine the effect of different saturation degrees and loading conditions on settlement shallow foundations using numerical modeling in Plaxis 2D, Bentley Systems, Exton, Pennsylvania, US. Plastic was used for dry soil calculation, while fully coupled ﬂow deformation was used for partially saturated soil. Pore pressure and deformation changes were computed in fully coupled deformation. The Mohr–Columb model was used in the simulation, and model parameters were calculated from experimental results. The study results show that the degree of saturation is more critical to soil settlement than loading conditions. When a 200 KPa load was applied at the center of the footing, settlement was recored as 28.81 mm, which was less than 42.96 mm in the case of the full-depth shale layer; therefore, settlement was reduced by 30% in the underlying limestone rock layer. Regarding settlement under various degrees of saturation (DOS), settlment is increased by an increased degree of saturation, which increases pore pressure and decreases the shear strength of the soil. Settlement was observed as 0.69 mm at 0% saturation, 1.93 mm at 40% saturation, 2.21 mm at 50% saturation, 2.77 mm at 70% saturation, and 2.84 mm at 90% saturation of soil.


Introduction
Foundations are typically built to achieve standards of strength and serviceability to support structure and equipment. Under serviceability conditions, the foundation must perform so that the structure or equipment it supports may fulfil its design purpose under typical operating loads. Settlement or other motion limits are commonly used to explain these serviceability limitations. Strength requirements are to ensure that the foundation has enough reserve strength to withstand enormous loads that may occur due to extreme environmental conditions or other sources. Serviceability or settlement and strength criteria can be addressed as separate design jobs in most circumstances. Serviceability is a longterm issue for the foundation that might be influenced by time-dependent consolidation characteristics. Foundation strength, also known as bearing capacity, may be a short-term issue, such as the construction of an embankment on an undrained clay foundation, or a long-term one in which the maximum foundation load may arise at an undetermined period in the future [1]. In designing shallow foundations, settlement is one of the main parameters [2,3]. Settlement assessment of shallow foundations and carrying capacity calculation are significant and typical geotechnical problems [4], and have been widely investigated as deterministic problems. The most prominent scientific publications emphasize the importance of the interaction between soil and footing settlement and portray this as of the impact on footing shape [5][6][7][8][9][10][11][12][13][14]. According to Das et al. [15], factoring in settlement is more important than considering the bearing capacity in shallow foundation design, particularly for foundation widths greater than 1.5 m, which is more common in engineering practices.
There are two kinds of shallow foundation settlements-immediate and secondary settlements of compression. Immediate settlement is encountered as load application after the structure is constructed [16,17]. Settlement of footings depends on many variables, including the shape and size of the footing, embedding depth, layering, soil mass nonhomogeneity, type of loading conditions, and saturation degree [18].
Loads are transferred to near-surface unsaturated soils, which change with hydrological events at shallow footings. Recent advances in unsaturated soil mechanics demonstrate that matric suction has a large influence on the strength and settlement of soils. Shallow footings have been located and built on near-surface unsaturated soils, ignoring the influence of matric suction on soil shear strength [2,[19][20][21][22][23][24]. At present, it is possible to analyze foundations by finite element methods, and limit equilibrium [25][26][27][28] and the finite difference method [29] have been widely used in recent years to determine the bearing capacity and settlement of footings.
The hydrological process in heterogeneous porous media has received a lot of interest and study [30][31][32][33], using analytical or numerical methods [34,35], since layered soil is significantly more widespread than homogeneous soil. The numerical modeling of shallow foundations in unsaturated layered soil using variably saturated conditions under varying loading and stiffness of soil has hardly been reported.
In recent years, the geotechnical problem of shallow foundation response was investigated by its stochastic nature. Scientific publications have investigated the uncertainty quantification of the material uncertainty in cohesive and non-cohesive soil materials. Uncertainty quantification analyses have led to probability density functions regarding many aspects of soil response such as the non-linear response of sand or the porous consolidation and failure of clays. Finally, footing settlement response under spatial variability of soil has also been investigated [11,[36][37][38][39][40][41][42][43][44]. Sivakugan and Johnson [45] developed a probabilistic system based on settlement records in the literature to calculate the risk associated with settlement prediction methods. Settlement and bearing capacity of foundation models with various vertical cross-sectional shapes under the vertically applied load action are presented on non-cohesive subsoil bases. Models of foundations of rectangular, wedge, and T vertical cross-sectional forms were experimentally tested and verified, with a study generally showing foundations with higher bearing capacity and lower settlement with rectangular vertical cross-sectional shapes rather than with wedge and T shapes, from which lower bearing capacity and higher settlement were reported [46].
The mechanical behavior of unsaturated soils, depending on the form of the soil and various pore-water and pore-air conditions, can be interpreted using either modified total Buildings 2022, 12, 693 3 of 15 stress or a modified efficient stress system. The technique proposed is tested in unsaturated cohesive soils with model base test results [22].
Due to challenges in measuring strength and deformation parameters as functions of matric suction and/or degree of saturation, only a few studies have been conducted to examine the effect of rainfall infiltration on the stability of shallow footings (DoS) [47][48][49][50][51][52][53][54][55][56]. The spatial and temporal change in the degree of saturation (DoS) of soil is impacted directly by numerous hydrological parameters, including water table depth, infiltration, flood, and drought [57,58].
Changes in moisture content and groundwater level highly affect the strength and the deformation properties of soil, thereby on the whole overlaying construction. However, geotechnical engineers neglect this topic in most cases, assuming that soil conditions will always remain unchanged. Thus, it is considered one of the foremost causes of foundation settlement, leading to various adverse effects on the overlaying constructions. It appears necessary to use various ground modification methods to stabilize the soil, and strengthen and restore foundations [59][60][61][62].
Floods, excessive rainfall, seasonal changes, and drought substantially impact foundation settlement behavior, which may exceed limits [63][64][65][66][67][68][69][70][71][72][73][74][75][76][77]. For locations where the near-surface soil is partially saturated during the structure's design life, the present design approach can be either conservative or unconservative, depending on the type of hydrological event. This process can cause settlements to exceed acceptable limits, jeopardizing the structure's stability. As a result, it is vital to estimate the additional settlements that may occur due to changes in water conditions to offer an adequate margin of safety [78,79].
The current study is dedicated to observing footing settlement under various degrees of saturation and matric suction in Plaxis 2D FEM Software, Bentley Systems, Exton, Pennsylvania, US. The numerical modeling of Plaxis 2D is to be conducted on soil from the Jamshoro area to assess settlement of shallow foundations under different soil saturation and loading conditions. This research consists of three phases: compilation of all the in situ and laboratory data available and extraction of soil profiles for each plot in the area; the use of well-known correlation to determine model parameters; and numerical modeling.

Experimental Methods
Undisturbed soil samples, namely shale and weathered limestone, were obtained by rotatory drilling. After sample collection, samples were transported to the laboratory of Mehran University Jamshoro for testing. Several tests were performed on the soil samples including sieve analysis and calculation of liquid limit, plastic limit, shear strength, and unconfined compression strength.
Extensive soil investigation and laboratory work have been conducted in the Jamshoro study area. The collected soil sample profile has a characteristic two-layer soil structure followed by a stiff layer of rocks, as shown in Figure 1. The soil was classified as A-7-5 according to AASHTO (8th edition), while CH was based on the USCS classification system. After soil extraction, the Casagrande testing method was used to calculate the Atterberg limits. The liquid limit was determined as 70%, the plastic limit as 30%, the shrinkage limit as 15%, and the specific gravity of the shale as 2.60. The modified proctor test found a maximum dry density of 1.9 g/cm 3 . Equations (1) and (2) were used to calculate the undrained cohesion (Cu) and the modulus of elasticity (E), which are represented in Table 1. The angle of internal friction and cohesion was determined using the shear box test as 11 • and 22 (kN/m 2 ).
The modulus of elasticity is determined using the following relationship: where Cu is the undrained cohesion, and qu is the ultimate load applied.

Parameters Values
The modulus of elasticity is determined using the following relationship: where Cu is the undrained cohesion, and qu is the ultimate load applied. The modulus of the welasticity of soil was determined with the correlation mentioned above and cited from the work of Sivrikaya et al. [80,81]. Thus, E is determined as 24,711 kN/m 2 , as mentioned in Table 1.
The Poisson's ratio was calculated as suggested by Pusadkar et al. [82] for CH and CL soil, determined as 0 from the relationship (Ψ = ϕ − 30).

Numerical and Boundary Conditions
The finite element analysis was conducted using Plaxis 2D ultimate v.21, Bentley Systems, Exton, Pennsylvania, US. The plain strain condition using the elastoplastic Mohr-Columb model was selected to simulate the behavior of the soil sample (Jamshoro Shale), with the parameters selected as drained. The overlaying foundation is modeled as a typical linear elastic with non-porous media. The symmetry axis and the proper vertical boundaries are constrained laterally. In both vertical and horizontal directions, the bottom boundary was restrained. Under load, soil and rock exhibit highly non-linear behaviour. The well-known Mohr-Coulomb model can be thought of as a first-order approximation of real soil behaviour. Mohr-Coulomb model parameters such as the modulus of the elasticity of soil, shear strength parameters such as cohesion, and angle of internal friction are inserted as input parameters. The dilatancy angle was zero. As a first-order (linear) simulation of actual soil efficiency, the known Mohr-Coulomb constitutive law was used.

Parameters Values
Angle of internal friction φ (degrees) The modulus of the welasticity of soil was determined with the correlation mentioned above and cited from the work of Sivrikaya et al. [80,81]. Thus, E is determined as 24,711 kN/m 2 , as mentioned in Table 1.
The Poisson's ratio was calculated as suggested by Pusadkar et al. [82] for CH and CL soil, determined as 0 from the relationship (Ψ = φ − 30).

Numerical and Boundary Conditions
The finite element analysis was conducted using Plaxis 2D ultimate v.21, Bentley Systems, Exton, PA, USA. The plain strain condition using the elastoplastic Mohr-Columb model was selected to simulate the behavior of the soil sample (Jamshoro Shale), with the parameters selected as drained. The overlaying foundation is modeled as a typical linear elastic with non-porous media. The symmetry axis and the proper vertical boundaries are constrained laterally. In both vertical and horizontal directions, the bottom boundary was restrained. Under load, soil and rock exhibit highly non-linear behaviour. The wellknown Mohr-Coulomb model can be thought of as a first-order approximation of real soil behaviour. Mohr-Coulomb model parameters such as the modulus of the elasticity of soil, shear strength parameters such as cohesion, and angle of internal friction are inserted as input parameters. The dilatancy angle was zero. As a first-order (linear) simulation of actual soil efficiency, the known Mohr-Coulomb constitutive law was used. Hook's law of isotropic elasticity underpins the linear elastic portion of the MCM. Based on the Mohr-Coulomb failure criterion, the component is perfectly plastic. The MCM does not incorporate either the stress nor the stress-path dependencies of the stiffness. In general, failure stress states can be adequately represented using the Mohr-Coulomb failure criterion and effective strength parameters. There was no assumption of material hardening or softening. Because there was no subterranean water in the profile at testing, no pore pressure was assumed during the study [83].
The distributed load was applied on the footing as a line load in Plaxis; and through the staged construction option, the footing was initially activated after the line load. After the geometry of the model was established and the material properties were assigned to all clusters, the next step was to divide the geometry into elements in mesh generation. A medium-sized mesh was selected for more accuracy and reduced program processing time. The Plaxis 2D software allows the 'Robust Triangulation Scheme' to automatically generate finite element meshes. A relatively coarse mesh may fail to capture the domain's significant responses, whereas probability accumulates numerical errors beyond the optimally fine mesh. Additionally, very fine meshing should be avoided because calculations take excessive time. With further provision of local refinements, as required by the merit of the problem and the position of the answer points in the numerical simulation, any simple meshing scheme can be adopted.
The plane strain model and 15 nodes were selected to simulate the soil medium. Plane strain models are used for shapes with a (more or less) uniform cross-section, stress state, and loading scheme along a certain length perpendicular to the cross-section (z direction). It is assumed that there are no displacements or strains in the z direction. Plane strain assumes that the problem being analyzed has an infinite length normal to the segment of the plane being analyzed. In a plane strain analysis, the out-of-plane displacement (strain) is zero by definition. The axisymmetric analysis is commonly applied to circular tunnels. For circular structures with a (more or less) uniform radial cross-section and loading scheme around the central axis, an axisymmetric model is utilized, in which the deformation and stress state are assumed to be the same in any radial direction. Compared to the 6 noded triangular components, it offers more nodes, and Gauss points to assist in the comparatively precise determination of displacements and stresses. The model dimensions have been selected so that the deformation in the soil does not intersect the model's boundaries. The width of the soil model was set as 20 m and the depth as 10 m. There is a need to determine the initial stresses in Plaxis simulations. For the specification of these stresses, two possibilities are available in the software: 'K 0 procedure' and 'gravity loading'. As a guideline, in the case of a horizontal surface and for any soil layer and phreatic lines parallel to the surface, the 'K 0 procedure' should be used. The "standard fixity" condition was employed in the numerical model. On the vertical edges, horizontal fixity was added, while the model's bottom edge is considered to be non-yielding and constrained from both vertical and horizontal movements. The numerical modeling was performed using an official Plaxis 2D v.21 Ultimate Licence at Saint Petersburg polytechnic university, Russian Federation, subscription type: CNTI-SPbPU (1006650066). Figure 2 shows the geometry of the model, in which the foundation modeled with a plate element and positive and negative interfaces was applied. The positive and negative interfaces are applied to soil and foundation.

Results and Discussion
Foundation settlement under load model conditions was simulated similar to under field conditions using Plaxis 2D. The initial stress state of the soil was generated before construction. The initial stage is often called the K0 procedure. K0 and Rinter are 0.67 and 1.0, respecively. The displacements after this initial stage are set to zero. The model's geometry from the second phase is inherited, and the vertical load is applied. The footing has a uniform load of 12 or 24 kPa, roughly corresponding to the equivalent load of a typical one-story and two-story building [84]. Settlement of the footing under a load of 12,

Results and Discussion
Foundation settlement under load model conditions was simulated similar to under field conditions using Plaxis 2D. The initial stress state of the soil was generated before construction. The initial stage is often called the K 0 procedure. K 0 and R inter are 0.67 and 1.0, respecively. The displacements after this initial stage are set to zero. The model's geometry from the second phase is inherited, and the vertical load is applied. The footing has a uniform load of 12 or 24 KPa, roughly corresponding to the equivalent load of a typical one-story and two-story building [84]. Settlement of the footing under a load of 12, 24, and 36 KPa is depicted in Figure 3, and it is noted that with the increase in load, settlement increased. The corresponding bending moment of the footing is shown in Figure 4.

Results and Discussion
Foundation settlement under load model conditions was simulated similar to under field conditions using Plaxis 2D. The initial stress state of the soil was generated before construction. The initial stage is often called the K0 procedure. K0 and Rinter are 0.67 and 1.0, respecively. The displacements after this initial stage are set to zero. The model's geometry from the second phase is inherited, and the vertical load is applied. The footing has a uniform load of 12 or 24 kPa, roughly corresponding to the equivalent load of a typical one-story and two-story building [84]. Settlement of the footing under a load of 12, 24, and 36 kPa is depicted in Figure 3, and it is noted that with the increase in load, settlement increased. The corresponding bending moment of the footing is shown in Figure 4.  With the intensity of load increased to 200 kPa, the load-deformation curve under a load of 50, 100, and 200 kPa is shown in Figure 5, and the crorresponding observed settlements are 3.7, 9.8, and 42.2 mm, respectively. Numerical modeling was performed by Altaweel et al. [85]. The purpose was to offer a numerical study that uses 3D Plaxis appli- With the intensity of load increased to 200 KPa, the load-deformation curve under a load of 50, 100, and 200 KPa is shown in Figure 5, and the crorresponding observed settlements are 3.7, 9.8, and 42.2 mm, respectively. Numerical modeling was performed by Altaweel et al. [85]. The purpose was to offer a numerical study that uses 3D Plaxis application's finite element analysis to evaluate the impact of clay soils on foundation settlement. This effect is explored using a 1 m wide strip foundation under different situations. With the intensity of load increased to 200 kPa, the load-deformation curve under a load of 50, 100, and 200 kPa is shown in Figure 5, and the crorresponding observed settlements are 3.7, 9.8, and 42.2 mm, respectively. Numerical modeling was performed by Altaweel et al. [85]. The purpose was to offer a numerical study that uses 3D Plaxis application's finite element analysis to evaluate the impact of clay soils on foundation settlement. This effect is explored using a 1 m wide strip foundation under different situations. In Figure 6, numerical modeling was performed on two-layered soil, with a difference in stiffness. Settlement decreased with the inclusion of a higher stiffness layer. In turn, the bending moment of the footing decreased, as shown in Figure 7. The thickness of the shale layer (top layer) was 1 m and that of the lower layer was 9 m. Similarly, settlement decreased considerably in the limestone layer (rock layer) below the shale layer. The effect of the deformation modulus on the deformation properties is shown in Figure   Figure 5. Load-deformation curve.
In Figure 6, numerical modeling was performed on two-layered soil, with a difference in stiffness. Settlement decreased with the inclusion of a higher stiffness layer. In turn, the bending moment of the footing decreased, as shown in Figure 7. The thickness of the shale layer (top layer) was 1 m and that of the lower layer was 9 m. Similarly, settlement decreased considerably in the limestone layer (rock layer) below the shale layer. The effect of the deformation modulus on the deformation properties is shown in Figure 8. The thickness of the footing was 0.75 m. The inteface between the shale and limestone layers influenced the deformation properties, as shown in Figure 9.         In Figure 10, settlement decreased due to unsaturated conditions. This is due to m tric suction dominating settlement. The unsaturated conditions were modeled with fu coupled flow deformation. The research presented by Liu et al. [86] is a novel investigat In Figure 10, settlement decreased due to unsaturated conditions. This is due to matric suction dominating settlement. The unsaturated conditions were modeled with fully coupled flow deformation. The research presented by Liu et al. [86] is a novel investigation into the numerical modelling of rainfall-induced shallow landslides in unsaturated layered soil using the variably saturated flow equation. A one-dimensional, transient, unsaturated groundwater flow problem in two-layered soil was investigated. The permeability coefficient of the layers was different. This demonstrates that the lowest FOS occurs during a rainstorm event at the interface between two consecutive soil layers.  In Figure 11, the effect of the degree of saturation on settlement of footing is depicted; the degree of saturation varied from 0.4 to 0.9, and as the degree of saturation increased, settlement increased. Changing the matric suction and saturation level affects the soil shear modulus, which directly impacts shallow foundation elastic settling. In general, Figure 10. Effect of dry and unsaturated conditions. In Figure 11, the effect of the degree of saturation on settlement of footing is depicted; the degree of saturation varied from 0.4 to 0.9, and as the degree of saturation increased, settlement increased. Changing the matric suction and saturation level affects the soil shear modulus, which directly impacts shallow foundation elastic settling. In general, raising the matric suction (or lowering the saturation level) significantly impacts foundation settlement reduction [87]. The change in the degree of numerical saturation modeling was performed in FORTON with the cam clay model by Mehnedritta and Sawant [88] to observe the effect of saturation on settlement. In that study, the degree of saturation varied from 85% to 95% and 100%, and they stated that the instantaneous displacement is reported to be significantly lower at 100 percent saturation than at lower degrees of saturation. At 100 percent saturation, the load is transmitted to the soil particles and water present in the voids, resulting in a decrease in soil volume due to pore water expulsion alone. However, in a partially saturated situation, the volume can be reduced by compressing the air voids and the water and soil particles. ship between negative pore pressure and degree of saturation. The matric suction depends on the soil type and permeability. For simplicity, the same modulus of elasticity for unsaturated soil is assumed as saturated soil. Figure 11. Effect of the degree of saturation (DOS). Figure 11. Effect of the degree of saturation (DOS). Figure 12 shows the variation in deformation with time from settlement to heave of the footing due to the development of matric suction in the soil under a load of 12 KPa. Figure 13 shows the contours of settlement and heave. In the first phase of dry conditions and the second phase, arbitrary small rainfall is applied to soil to activate negative pore pressure (matric suction). In dry conditions, under the load, settlement is depicted in Figure 13 as 0.00080 m; and in unsaturated conditions, under the same load, the footing shows a heave of 0.000485 m. The generation of negative pore pressure is the main reason behind this phenomenon. The soil water characteristic curve parameters were taken from the database of Plaxis 2D for clay type of soil. The soil water characteristic curve governs the behavior of matric suction of soil, more often termed as SWCC, which is the relationship between negative pore pressure and degree of saturation. The matric suction depends on the soil type and permeability. For simplicity, the same modulus of elasticity for unsaturated soil is assumed as saturated soil.

Conclusions
Finite element research was undertaken using the PLAXIS 2D to predict the settlement behavior of shallow foundations. Analyses were carried out utilizing the Mohr-Columb model under different soil conditions, loading, and saturation degree (saturated and unsaturated). Based on the findings of this study, the following conclusions can be

Conclusions
Finite element research was undertaken using the PLAXIS 2D to predict the settlement behavior of shallow foundations. Analyses were carried out utilizing the Mohr-Columb model under different soil conditions, loading, and saturation degree (saturated and unsaturated). Based on the findings of this study, the following conclusions can be drawn: • Settlement is depicted with load-settlement and time-dependent settlement variation. Settlement decreased with increased soil stiffness (the modulus of elasticity) and increased with loading intensity. • Settlement at the center of footing at a load of 200 KPa is 28.81 mm, which is less than 42.96 mm in the case of the full-depth shale layer. Thus, settlement reduced by 33 percent when the underlying limestone rock layer was present. These results show that the presence of an underlying limestone layer decreasse settlement of the shallow foundation to significant level.

•
The footing settlement under various degrees of saturation (DOS) was observed. It was found that settlment is increasing by increasing the degree of saturation, which increases pore pressure and decreases the shear strength of the soil. Settlement was observed as 0.69 mm at 0% saturation 1.93 mm at 40% saturation, 2.21 mm at 50% saturation, 2.77 mm at 70% saturation, and 2.84 mm at 90% saturation of soil. The increase in the degree of the saturation of soil settlement increased.

•
The heave of the footing was observed with an increase in matric suction followed by settlement.

•
The changing of matric suction due to water conditions with time should be considered while designing of shallow foundations.

•
Heave of the foundation was observed when soil conditions changed from saturated to unsaturated. • ThAe proper drainage system should be provided for the foundation, so that there is a minimum effect of changing water conditions of soil.