Winkler Spring Stiffness Distribution for the Structural Analysis of Mat Foundations in Clay Under Drained Conditions

Winkler Spring Stiffness Distribution for the Structural Analysis of Mat Foundations in Clay Under Drained Conditions | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Research Article Winkler Spring Stiffness Distribution for the Structural Analysis of Mat Foundations in Clay Under Drained Conditions Gaby Saad, Grace Abou-Jaoude, Dimitrios Loukidis This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-6456559/v1 This work is licensed under a CC BY 4.0 License Status: Posted Version 1 posted You are reading this latest preprint version Abstract The interaction between mat foundations and the supporting soil is often simplified in design practice by assuming that the mat rests on Winkler springs. Most of previous studies aimed at determining the spatial variation of the spring constants or the modulus of subgrade reaction that is needed to obtain realistic and accurate results regarding the bending of the mat assuming that the soil behaves as a purely linear elastic material. This paper investigates the effects of plastic strains and volume changes that happen due to the consolidation of clayey soils under long-term conditions on the distribution of the equivalent spring constants and the bending moment diagrams. For this purpose, three-dimensional finite element analyses are performed for fully drained conditions with the soil simulated using the Modified Cam-Clay constitutive model. The results show that the presence of plastic yielding significantly influences the mat-soil interaction and the assumption of a linear elastic soil overpredicts the soil stiffness and the bending moments in the mat foundation. The semi-analytical Discrete Area Method could be adopted to obtain comparable results to those generated by elastoplastic finite element analysis. Soil-structure Interaction Mat foundation Finite element analysis Winkler Springs Modified Cam-Clay model Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Figure 7 Figure 8 Figure 9 Figure 10 Figure 11 Figure 12 Figure 13 Figure 14 Figure 15 Figure 16 Figure 17 Figure 18 1. INTRODUCTION The behavior of structures is affected by the properties and conditions of the foundation soil, while at the same time the structure and its loading influence the behavior of the soil. Understanding this interplay between soil and structure is fundamental in design to ensure the stability, performance, and safety of the structures. Soil-structure interaction is particularly important for the design of mat (or raft) foundations of heavy structures, where minimum reinforcement often does not govern the structural design. Consequently, engineers are forced to adopt sophisticated methods for analysis and design. This is particularly evident in geotechnical engineering practice where there is progressive shift from traditional methods to high-level and computationally expensive numerical analysis to ensure precision and reliability involving the design of foundations of high-rise buildings (Poulos 2016). On the other hand, for more common projects, simplified subgrade models are still used in the context of structural analysis software. Engineers favor mat foundations when dealing with large loads that must be transferred to weak ground. They are preferred because they can resist differential settlements, heave and bending moments more effectively than other types of shallow foundations. Mat design is mainly based on calculations of settlements, shear forces, and bending moments to determine the required mat thickness and steel reinforcement. Several subgrade models have been developed to study the interaction between the structure and the underlying soil. The soil-structure interaction of a beam foundation modeled as a structural element laying on a bed of independent linear springs representing the soil was proposed by Winkler (1867). The stiffness of the springs depends directly on the modulus of subgrade reaction k s , defined as: $$\:{k}_{s}=\frac{q}{\delta\:}$$ 1 where \(\:\delta\:\) is the settlement generated by the applied pressure q . The spring constant ( K s ) is obtained by multiplying k s by the tributary area of each node. According to Horvath (1995), k s is the most important parameter considered in the design of mat foundations. However, its calculation is not simple because it is not a fundamental property of the soil. It is a function of the stiffness parameters of the soil, as well as the geometrical characteristics of the foundation (Coduto 2001). It has long been known that the simple bed of springs with constant k s throughout the mat (uncoupled springs approach) cannot capture the actual deformation at the soil-mat contact surface (Horvath 1995). The coupled springs method represents a step up from the original Winkler spring method, where the addition of other mechanical elements (such as a stretched membrane under tension to account for the soil deformation continuity or a shear layer of defined shear modulus to account for the effects of shearing) simulate the soil response more realistically and allow the prediction of the correct settlement profile. However, its mathematical formulation is complex and requires careful considerations regarding the boundary conditions. The pseudo-coupled method is an alternative approach that relies on assuming a variable spatial distribution of the modulus of subgrade reaction k s across the mat. This method attempts to overcome the lack of coupling in the original Winkler spring method, while retaining its simplicity. Given that many structural design software already employ the Winkler spring method, accommodating the pseudo-coupled method is feasible and represents the most practical approach to date. The pseudo-coupled method was further improved by adding an iterative procedure that relies on the analytical Boussinesq solution (Ulrich 1995). This improved method is known as the Discrete Area Method (DAM). Although Ulrich (1995) advocated the use of DAM in design and highlighted its effectiveness in producing accurate results, current design codes do not enforce or make reference to this procedure. More recently, using 3D finite element modeling, Loukidis and Tamiolakis (2017) and Leonidou (2021) performed parametric analyses of mat foundations of small planar dimensions placed on linear elastic soil to back-calculate the spatial variation of the spring constants that is needed to achieve an accurate mat structural analysis. They concluded that when the soil-structure interaction is properly simulated using a variable k s distribution, the peak moments corresponding to bottom fiber in tension of the mat foundation increase and those corresponding to top fiber in tension decrease compared to those predicted using the uncoupled Winkler approach. Alzoaby et al. (2022) applied the DAM, using commercially available software SAFE and Settle3, to study the effects of the variable k s distribution on the steel reinforcement of large mat foundations underlain by soft to very stiff linearly elastic soil. This work showed that implementing a variable k s distribution in the design of mat foundations drastically affects the moment distribution of mat foundations and thus steel reinforcement leading to a more accurate and realistic design than using uncoupled springs or other pseudo-coupled approaches. The aforementioned studies have highlighted the importance of varying the modulus of subgrade reaction under mat foundations. However, they all revolved around characterizing the soil as a linear elastic material. The focus of the present study is the effect of long-term plastic deformation on the distribution of the equivalent spring constants for mat foundations resting on clay under long-term fully drained conditions, i.e., end-of-consolidation state. Eurocode EN1997-1:2004 (CEN 2004) explicitly states as a basic requirement for geotechnical design that both short-term and long-term conditions shall be considered as design situations (§ 2.2( 1 )). This concerns also the design against structural failure of the foundation (§ 2.4.7.1, § 6.8( 5 )). The nonlinear behavior of the soil is expected to cause a redistribution of the contact pressure under the foundation, which in turn will affect the bending moment distribution. This is attributed to the fact that the soil in the vicinity of the edges of the mat foundation is subjected to much more shearing than at the mat center and, as a result, the soil near the edges is prone to plastic yielding early during the mat loading process. The aim of the present study is to explore the influence of the long-term (consolidation) settlements on the bending moments developing in mat foundations resting on clay and propose ways to take into account their effect in design practice in the context of the pseudo-coupled approach. To achieve this objective, three-dimensional finite element analyses of mat-soil interaction are conducted using the software ANSYS. The soil is simulated as an elastoplastic material using the Modified Cam-Clay constitutive model (Roscoe & Burland 1965). The Modified Cam-Clay model can simulate the nonlinear hardening of normally and slightly overconsolidated clays and the softening of strongly overconsolidated clay (under drained conditions). Most importantly, it can capture accurately the development of plastic volumetric strains due to exceedance of the preconsolidation pressure. A parametric study is performed varying the soil properties, foundation characteristics, and loading magnitude. The calculated settlements at the nodes are used to back-calculate the equivalent linear spring constants K s . The bending moment and shear force diagrams obtained are compared to those generated from three common methods of assigning modulus of subgrade reaction under the mat. 2. KNOWLEDGE BACKGROUND Soils does not exhibit a purely linear elastic behavior for the entire range of loading prior to failure and, beyond a certain shear stress level, they exhibit nonlinearly hardening. Moreover, unlike structural materials, this nonlinear behavior actually initiates at early stages of shearing. Hence, when analyzing structures that induce high levels of stresses in the soil, as in the case of multistorey buildings, it is quite inaccurate to consider solely linear elastic behavior for the soil and it is reasonable to expect that this affects negatively the design of the mat and the structure. Researchers have investigated the role of nonlinear soil behavior in the soil-foundation interaction using two types of material models: nonlinear elastic models and elastoplastic models. 2.1 Nonlinear elastic soil behavior A nonlinear elastic model describes a nonlinear behavior of the soil in the “elastic zone” prior to failure by assuming that the elastic modulus is a decreasing function of the shear strains and an increasing function of the mean effective stress. Though an improvement compared to the purely linear elastic model, it still limits the representation of the soil behavior by not considering the effects of dilatancy and plastic (irreversible) deformations. Viladkar et al. (1991) demonstrated the significant influence of soil’s nonlinear behavior on foundation response and highlighted the necessity of nonlinear analysis techniques for accurately predicting settlement profiles, shear forces, and bending moments in rigid combined footings. The soil compressibility affects the behavior and design of foundations as well as the distribution of contact pressure under the footing, which influences the equivalent k s distribution. Moreover, it was observed that the shape of the contact pressure varies based on the soil properties and the magnitude of the loading. ACI 318 − 19( 22 ) states that the distribution of soil pressure under a mat foundation should be consistent with the properties of the soil and is highly dependent on the soil type. Teng (1962) emphasized that the calculation of the distribution of the modulus of subgrade is different for clays and sand. Similar outcome is observed in analyses of mat foundations, in which a redistribution of contact pressure is obtained when considering the nonlinear elastic behavior of soil (Viladkar 1994). However, a nonlinear elastic analysis produces much higher settlements and therefore it should be limited to only small loading magnitudes, since it is unreliable at high loading magnitudes (Noorzai et al. 1995). The soil could be idealized as a nonlinear elastic material provided that the loading applied does not produce extensive regions of plastic straining in the soil mass, otherwise it should be modeled using an elastoplastic constitutive mode (Dutta and Roy 2002). 2.2 Elastic - perfectly plastic behavior The elastic – perfectly plastic model assumes a linear elastic behavior prior to failure and a perfectly plastic behavior after that, which is governed by a failure criterion (usually the Mohr-Coulomb criterion). As such, it ignores post-failure strain-softening effects and pre-failure nonlinear hardening behavior, as well as plastic volumetric deformations due to consolidation and creep. Plastic yielding occurs first at the edges of the foundation and results in the redistribution of the stresses under it. Noorzai et al. (1995) illustrated that, for combined footings, elastoplastic effects become prominent in the case of heavy structures. Their results show that the settlement profiles generated using elastoplastic soil model are almost identical to those obtained assuming linear elastic soil behavior under low loading conditions and differ at high loading conditions due to the effects of plastic yielding which occurs mainly at the edges of the mat foundation. This is explained by the observation that as loads increase, plastic deformations spread in the soil towards the center of the foundation (Noorzai et al. 1995; Abdullah 2008). Moreover, the combined footing deformed in the form of a “cup” shape, having the maximum settlement below heavily loaded areas only when the soil adheres to a yield criterion. Calculating the stress distribution using the Boussinesq theory, where the soil is treated as a linear elastic homogenous isotropic half-space, underestimates the vertical stresses underneath the foundation down to a depth approximately twice the foundation width, which is an area that is contributing significantly to the settlement of the foundation (Sadek and Shahrour 2007). This discrepancy is attributed to the development of plastic shear deformations that lead to a reduced attenuation of vertical and horizontal stresses within the soil mass under loading (Sadek et al. 2010). Abdullah (2008) argued that, once soil yielding occurs at and near the foundation’s edges, the assumption of uniform contact pressure used in practice is inaccurate as it may lead to unsafe mat foundation design. Recently, Pishilis (2022) performed finite element analyses using the computer program Abaqus to derive the distribution of equivalent Winkler spring constants considering elastoplastic behavior under constant volume (undrained clay behavior). Unlike Loukidis and Tamiolakis (2017), where they explain that for elastic soil analysis the applied loading magnitude does not affect the shape of the k s distribution, the study of Pishilis (2022) showed that the shape of the distribution varies with the applied loading magnitude. Accounting for soil yielding in design increases the accuracy in design significantly (Larkela et al. 2013; Pishilis 2022). Ulrich (1991) and Banavalkar (1995) highlighted the crucial role of both immediate and long-term subgrade responses in the design of mat foundations. They pointed out that the determination of the spring constants for the design of mat foundations should reflect the time-dependent behavior of the soil, which in the case of clays is governed by the phenomenon of consolidation. Despite this assertion and the Eurocode 7 previsions, the long-term conditions for foundations in clays are largely overlooked in the structural design practice. 3. METHODOLOGY 3.1 Continuum finite element analysis Three-dimensional finite element analyses (FEA) using solid elements for the soil (continuum FEA or cFEA) were performed using ANSYS to derive the distribution of the equivalent Winkler spring constants ( K s ) across flexible mat foundations resting on clay under drained (end-of-consolidation) conditions. The settlement profiles were extracted from these analyses at each mat node and transferred into a MATLAB code that back-calculates the equivalent K s distribution (Tamiolakis 2012). The soil domain was discretized using 8-noded (first-order) hexahedral elements (type Solid185) having 3 translational degrees of freedom per node, with full Gauss integration (Fig. 1 ). The mat foundation, placed at the center of the free surface of the soil domain, was modeled with 4-noded shell elements (type Shell 181) having 6 degrees of freedom per node (3 translational and 3 rotational degrees of freedom). The use of first-order elements eliminates the complexities of deriving an equivalent spring constant distribution associated with higher order elements, which include mid-side nodes in addition to corner nodes. As it was observed that the differences in results between using reduced and full integration shell elements are insignificant, reduced Gauss integration was adopted to reduce the computational cost (Loukidis and Tamiolakis 2017). The mat was meshed with square elements 0.5m wide, a size commonly used in engineering practice. The soil elements directly underneath the foundation were cubical with the same edge length (0.5m) as the mat shell element in order to have coincidence of mat and soil nodes at their contact. The same element mesh size was extended to an area twice the foundation’s horizontal dimensions. Beyond that region the element size was allowed to increase to reduce computational time. The interaction between the soil and the foundation was defined as “shared topology”. This configuration forces the mat and soil nodes to be rigidly connected to each other with respect to the translational degrees of freedom and prevents any possibility of slippage or separation between them. The lateral boundaries of the soil domain were placed at a distance from the mat center that was 10 times the width B of the mat, and the thickness H of the soil domain was set to be 10 times the B , which is sufficient for the boundaries to have minimal effect on the numerical results (Loukidis and Tamiolakis 2017). Each finite element parametric analysis comprised three phases. First, a geostatic stress field was established while applying gravity to the soil domain. Second, the self-weight on the mat was applied. Third, the column loads were applied on the mat as concentrated vertical forces. The mat material (reinforced concrete) was modeled as linear elastic, having Young’s modulus 32 GPa, Poisson’s ratio 0.2 and unit weight 24.5 kN/m 3 (density of 2500 kg/m 3 ). The constitutive model for the clayey soil was the built-in version of Cam-Clay model in ANSYS, which can be rendered identical to the Modified Cam-Clay (MCC) model proposed by Roscoe & Burland (1968) through proper selection of the values for the model input parameters. To ensure that the constitutive model formulation and parameters have been correctly interpreted, the results of a single element analysis performed in ANSYS simulating undrained triaxial compression test were compared with MCC predictions presented in Potts and Zdravkovic (1999). To check the adequacy of the finite element mesh and boundary conditions, trial FE analyses were performed assuming a linear elastic soil and a very stiff mat and the resulting settlements were compared to a well-established analytical solution that considers the mat as a rigid foundation, namely the Steinbrenner (1934) equation coupled with the Chow (1987) rigidity factor. The FEA results are close to those obtained from the analytical solution, with an average error of 2.44%. 3.2 Parametric analysis A parametric study was conducted to investigate the effect of plasticity under long-term conditions (end-of-consolidation) on the K s distribution, as well as the developing bending moments, which control the design of the steel reinforcement. The mat foundations considered herein have planar dimension 20m × 20m. The arrangements of column loading is shown in Fig. 2 . Each mat carries 25 equally-spaced column loads, with the column spacing ( S column ) being 5m (typical for multistorey residential buildings). Regarding the column relative magnitudes, most of the finite element analyses were performed for a typical symmetric loading configuration C1, where the columns transfer a load according to their tributary area. Configurations C2 and C3 induce single axis and double axis load eccentricity respectively (i.e., the resultant of the column loads does not pass through the center of the mat), subjecting the foundation to rotation. The thickness of the mat ranges from 0.5 m to 1m. These values satisfy punching shear criteria for the load magnitudes considered in the parametric analyses. The quantify the rigidity of the mat foundations, the relative stiffness factor was determined according to the following equation (Meyerhof 1953): $$\:{R}_{\text{s}}=\frac{{E}_{c}{d}^{3}}{12E{B}^{3}}$$ 2 where E c is modulus of elasticity of the mat concrete, d is the thickness of the mat, E is the modulus of elasticity of the soil and B is the mat width. A foundation is characterized flexible if R s < 0.5. Unlike a linear elastic soil, the back-calculated spring constants for an elastoplastic soil depend on the loading magnitude Q . Loading from inner columns (i.e., those that do not lie along the perimeter of the mat) ranges from 500kN (approximately 2 storeys) to 4000kN (approximately 16 storeys). The critical friction angle φ c of the clayey soil ranges from 25° to 35°, the Poisson’s ratio ν from 0.15 to 0.35, the compression index C c from 0.1 to 0.6 (i.e. MCC parameter λ from 0.043 to 0.26), the swelling index C s from 0.0175 to 0.07 (i.e. MCC parameter κ from 0.0076 to 0.3), and the at-rest lateral earth pressure K 0 from 0.5 to 1. The most common cause of overconsolidation in nature is the removal of overburden due to erosion, in which case the preconsolidation pressure p c is expected to increase practically linearly with depth, while the initial void ratio e ini would decrease nonlinearly with depth. Details on the determination of the p c and e ini profiles as a function of the assumed thickness D of eroded overburden is given in the Appendix. These profiles were applied as initial conditions to the ANSYS finite element model by prescribing a user-defined pseudo-temperature distribution with depth and making p c and e ini dependent to it. The D values considered herein range from 5 m to 40 m, corresponding to vertical preconsolidation stress at the ground surface from 41 kPa to 326 kPa. The geostatic stresses were established in the model using the buoyant unit weight γ ′ (assuming that the water table is at the free surface). Given that the void ratio decreases with depth, γ ′ ought to be an increasing function of depth. However, because ANSYS does not allow the unit weight to be dependent on the pseudo-temperature, γ ′ was assumed constant with depth, equal to an average representative value of 8.15 kN/m 3 . 3.3 Determination of equivalent spring constant distribution To determine the equivalent Winkler spring constants from the results of the cFEA, the vertical displacements of all the nodes of the mat from each analysis were exported from ANSYS and introduced as input in the MATLAB code developed by Tamiolakis (2012). The code was modified for the purposes of the present study to account for the stiffness matrix of the shell elements used in ANSYS and to account for the foundation self-weight. This algorithm solves inversely the problem of a mat foundation resting on Winkler springs to produce as output the equivalent spring constants ( K s ) and is based on the following system of static equilibrium equations: $$\:\left(\left[{\text{K}}_{\text{m}\text{a}\text{t}}\right]+\left[{\text{K}}_{\text{s}\text{p}\text{r}\text{i}\text{n}\text{g}\text{s}}\right]\right)\left[\text{U}\right]=\left[{\text{F}}_{\text{m}\text{a}\text{t}\:\text{s}\text{e}\text{l}\text{f}-\text{w}\text{e}\text{i}\text{g}\text{h}\text{t}}+{\text{F}}_{\text{c}\text{o}\text{l}\text{u}\text{m}\text{n}\:\text{l}\text{o}\text{a}\text{d}\text{s}}\right]\:$$ 3 where [K mat ] is the assembled stiffness matrix of the mat foundation elements, [K springs ] is the assembled stiffness matrix of the Winkler springs below the mat, [U] is the matrix of nodal vertical displacements and rotations, \(\:{F}_{mat\:self-weight}\) is the matrix of the equivalent nodal loads due to the mat self-weight, and \(\:{F}_{column\:loads}\) is the matrix of the applied column loads. Note that [K springs ] is a diagonal matrix because each Winkler springs is affecting only the vertical nodal displacement of the node to which it is connected. 3.4 Determination of representative Young’s modulus To facilitate comparisons of the distribution of equivalent spring constants and bending moments from an elastoplastic analysis to those generated by a purely linear elastic analysis, a representative value of the Young’s modulus ( E * ) was calculated at a depth below the mat center equal to one-third the mat width ( B /3) at the end of the elastoplastic analysis (after full application of the mat loading). This calculation was necessary because, in the Modified Cam-Clay model, the soil Young’s modulus E is not an input model parameter and is actually a variable that depends on the current values of mean effective stress p ′ and void ratio e : $$\:E=3\left(1-2v\right)\frac{\left(1+e\right){p}^{{\prime\:}}}{\kappa\:}\approx\:6.9(1-2v)\frac{(1+e)p{\prime\:}}{{C}_{s}}$$ 4 As such, the soil Young’s modulus varies spatially, and as mat loading is applied. Moreover, the soil near the edges undergoes plastic yielding (where shear straining is dominant). Plastic deformations may develop also near the mat center if the preconsolidation pressure p c is approached during loading, otherwise the soil in that region remains largely in an elastic state. The selection of this location for the reference point for the E * determination (depth B /3 below the mat center) was made using results of trial analyses and based on the rational that, given that in most analyses, the soil under the mat center remains in an elastic state (i.e. the stress state lies inside the Modified Cam-Clay yield surface), the spring constant at the center of the foundation from a linear elastic analysis should generally align with that produced at the center of the corresponding elastoplastic analysis. The E * values among the analyses of the parametric study range from 7.6 MPa to 30.8 MPa, leading to relative stiffness factor R s in the 0.011–0.044 range. The E * is used for calculating a reference spring constant value, K r , with respect to which the K s distributions are normalized. This facilitates comparison between the K s distributions from the various parametric analyses performed and discerning the effect of each problem parameter on the shape of the K s distribution. The K r corresponds to the spring constant that would be applicable to a perfectly rigid foundation resting on linear elastic soil having Young’s modulus E equal to E * and is calculated as follows: $$\:{K}_{\text{r}}={k}_{\text{r}}{A}_{\text{i}\text{n}}$$ 5 where A in is the tributary area of the inner mat nodes (0.25 m 2 ) and k r is the modulus of subgrade reaction for rigid foundation: $$\:{k}_{\text{r}}=\frac{{E}^{*}}{{C}_{f}\left(1-{v}^{2}\right)B}$$ 6 with the C f factor proposed by Loukidis and Tamiolakis (2017), which introduces the effects of the mat aspect ratio ( L / B ) and the deformable soil layer thickness H : \(\:{C}_{f}=\frac{0.85\:{\left(\frac{L}{B}\right)}^{0.45}}{{\left[1+0.1(2+\frac{L}{B})\frac{B}{H}\right]}^{1+{e}^{5{v}^{3}}}}\) ( 7 ) 3.5 Discrete Area Method The Discrete Area Method (DAM) using a combination of the computer programs CSI SAFE and Settle3 was first applied by Alzoaby et al. (2021) for elastic soils. Herein, we expand this computational approach by performing the settlement calculations in Settle3 using the built-in nonlinear 1-D consolidation model, where parameters such swelling index C s , compression index C c and preconsolidation pressure p c are inputted directly in the software, in order to make consistent comparisons with the cFEA-based results. Nonetheless, it has to be pointed out that, although a nonlinear model is used in Settle3, the program still calculates the stress changes induced by the application of the foundation loading using the Boussinesq solution, which stems from the theory of elasticity. Hence, due to this inconsistency, the resulting settlements, which in the case of a linear elastic soil would constitute exact solutions, in the present application of DAM should be considered approximate. The DAM is an iterative approach. First the mat-on-springs problem is solved in CSI SAFE using an assumed distribution (initial guess) of modulus of subgrade reaction. The computed contact pressures are then inputted in Settle3 and new values of modulus of subgrade reaction at each mat node are obtained by dividing the contact pressures by the corresponding settlement values at each mat node. The new distribution of the modulus of subgrade reaction is inputted back in CSI SAFE to recalculate soil pressures under the mat foundation. The iterative procedure stops when the settlement estimates from CSI SAFE and Settle3 converge. Convergence is considered to be achieved when the relative difference between the settlements calculated by CSI SAFE and Settle3 does not exceed 5% for every mat node: $$\:{err}_{U}=\left|\frac{{\left(U\right)}_{CSI\:SAFE}-\:{\:\left(U\right)}_{Settle3}}{{\left(U\right)}_{CSI\:SAFE}}\right|\times\:100\:$$ 8 Usually, an average of 11 iterations is sufficient to achieve convergence. 4. RESULTS AND DISCUSSION Figure 3 compares the vertical displacement distribution across the mat foundation from cFEA with those obtained using i) a uniform k s distribution (i.e., uncoupled approach) and ii) the pseudo-coupled approach of Coduto (2001). In the latter, the mat is divided into 3 concentric zones, each with its own value of k s . The inner zone has half the mat dimensions ( B /2 × B /2), while the intermediate and outer zones have thickness equal to B /8 (Loukidis and Tamiolakis 2017). The modulus of subgrade reaction progressively increases from the central to the outermost zone, such that the outermost zone has a k s value twice as much as that in the inner zone. The uniform k s value and the average value of k s across the mat for the pseudo-coupled approach are set equal to k r in Eq. ( 6 ). It can be seen that the choice of the k s distribution strongly affects the displacement profile. The uniform k s distribution leads to less variation in the vertical displacements and, thus, smaller differential settlements. Most importantly, it results in a deformed shape of the mat that is convex (concave downwards). On the contract, the cFEA and the pseudo-coupled approach result in a concave upwards deformed shape. As will be shown later, the differences in the mat settlement profiles have a significant impact on the bending moment diagrams. The average magnitude of mat settlement from the two linear elastic approaches (uncoupled and pseudo-coupled) is comparable to that from elastoplastic cFEA, with the differences being less than 25%, indicating the effectiveness of the choices made for the calculation of E * . As the column magnitude increases, the differential settlements from cFEA and the pseudo-coupled approach become smaller (Fig. 4 ), and the mat “dishing” is less pronounced. This is because the larger applied loading induces larger shear stresses and strains in the soil near the perimeter of the mat, leading to more plastic deformations. Due to the plastic deformations, the overall (secant) stiffness of soil close or at the mat perimeter becomes smaller compared to that near the mat center. The development of plastic straining is expected to be reflected in the shape of the equivalent K s distribution, rendering it load-magnitude-dependent. 4.1 Equivalent K s distribution and bending moment diagrams Figure 5 show examples of back-calculated spring constants K s from three finite element analyses with different clay properties and loading magnitude. Similar to the case of linear elastic soil (Loukidis and Tamiolakis 2017), the K s distributions exhibit a distinctive “cup” shape, with high values at the corners and along the edges of the mat that decrease sharply towards the foundation’s central region, where an almost constant value is attained. In the remainder, comparisons of the parametric analysis results will be made based on cross-sections of the K s distribution along the mat centerline and on along a mat edge, with the values normalized with respect to the reference value K r . 4.1.1 Effect soil yielding The effect of the critical state friction angle φ c on the normalized spring constant distribution ( K s / K r ) can be seen in Fig. 6 a which shows the results of three analyses with friction angle varying from 25° to 35°. We can see that K s / K r is unaffected by the critical friction angle in the inner region of the mat since the soil at the center of the foundation remains mainly elastic. However, K s / K r along the perimeter of the mat decreases as the friction angle decreases. This is because the lower shear strength of a clay having small φ c results to earlier plastic yielding near the edges and more plastic straining. This change in edge spring constants affects the bending moments. More specifically, increasing the friction angle leads to a slight increase in the positive (bottom fiber in tension) moments throughout the central axis of the mat foundation because of the higher spring constants at the perimeter of the mat (Fig. 6 b). The dependence of the edge spring constants on plastic straining gives also rise to dependence on the loading magnitude. So, unlike the findings of Loukidis and Tamiolakis (2017) for purely linear elastic soil behavior, increasing the loading magnitude applied on the foundation and the mat thickness (increased self-weight) leads to a progressive reduction of the K s / K r especially at the edges (Figs. 7 and 8 ). Figure 8 a indicates that as the load transferred to the ground increases, K s / K r values at the edges decrease substantially (by 50% or more for every doubling of Q ). This decrease can be attributed to large shear strains near the perimeter of the mat, but also on large contact stresses in the same region, which may cause yielding also due to increased p ′ values that approach (or even exceed p c ) leading to the development of contractive plastic volumetric (consolidation) strains. On the contrary, the decrease in the central region of the mat is rather insignificant for Q up to 1000kN (total loading, including mat self-weight 25807kN), corresponding to a mobilized factor of safety roughly equal to 5.0 under drained conditions (1.84 under undrained conditions). For larger loading (smaller mobilized factor of safety values), the spring constants decrease notably also in the central region of the mat. This observation suggests that the expansion of the plastic deformations towards the mat central region is starting to affect the central K s values. However, due to the limited presence of shearing under the central region, this rate of decrease is smaller than that for the edge springs and, as a result, the edge and central K s distributions progressively converge to each other. Regarding the bending moments, as expected, the more rigid the mat foundation is (large mat thickness) or the higher the applied loading is, the larger the bending moments across the mat foundation become. However, if the moment diagrams are normalized with respect to sum of the applied loading (Σ Q ), the fact that the K s values at the edge and the center converge leads to a decrease in the normalized positive (bottom fiber in tension) moments (Fig. 8 b). This means that doubling the mat loading does not cause a doubling of the peak moments, as would be expected if the soil were assumed to be a linear elastic material. The smaller is the coefficient of lateral earth pressure at rest K 0 (= σ ′ horizontal / σ ′ vertical ) from unity (1.0), the closer is the initial soil state of an overconsolidated clay to the yield surface. A K 0 value equal to 1 means that the effective vertical and horizontal stresses are equal, and thus the initial (in-situ) deviatoric stress is equal to zero. Contrarily, a K 0 = 0.5 means that the effective vertical stresses are twice as large as the effective horizontal stresses, bringing the in-situ soil stress state closer to shear yielding. Hence, the clay under the mat (especially close to the perimeter where shearing is most dominant) is expected to develop plastic shear strains early during loading, resulting to smaller secant stiffness and, thus, comparatively smaller equivalent spring constants. The results from analyses with K 0 varying between 0.5 and 1.0 (Fig. 9 a) indicate that, as K 0 decreases, the ratio of mat edge K s values to those at the mat center becomes smaller. However, this reduction in K s contrast is small and, as a consequence, the K 0 parameter has minimal effect on the bending moment diagram (Fig. 9 b). 4.1.2 Effect of clay deformability The swelling index C s controls the elastic stiffness in the MCC model and is directly related to the Young’s modulus E in Eq. ( 4 ). According to the results shown in Fig. 10 a, the normalized spring constant ( K s / K r ) distributions are quite insensitive to the C s value. Its influence is limited to a slight increase of the edge values as the swelling index decreases. Most importantly, the shape of the distributions remains practically the same and, as a consequence, a change in C s has minimal effect on the bending moment diagram (Fig. 10 b). Figure 11 a reveals that, by decreasing the Poisson’s ratio ν , K s / K r increases at the edges, where the soil is subjected mostly to shearing. This observation can be explained by the shear modulus equation Eq. ( 9 ) for the MCC model: $$\:G=\frac{3\left(1-2v\right)}{2\left(1+v\right)}{K}_{\text{b}}=\frac{3\left(1-2v\right)}{2\left(1+v\right)}\frac{(1+e){p}^{{\prime\:}}}{\kappa\:}$$ 9 In the MCC model, the bulk modulus K b is assumed to be independent of ν ; it is only affected by the change in mean effective stress, the swelling index and the void ratio. Hence, by decreasing ν , the shear modulus G becomes larger and the reaction to shear deformation gets enhanced. Nonetheless, this increase in the Poisson’s ratio does not practically affect the bending moments (Fig. 11 b). The impact of the compression index ( C c ) on the K s / K r distribution depends heavily on the magnitude of the applied pressure and to what extent the preconsolidation pressure is exceeded at the foundation level. Figure 12 a illustrates that C c is not affecting the K s / K r distribution for a moderately overconsolidated profile soil (eroded overburden soil D = 10m, corresponding to overconsolidation ratio OCR = 2.5 at a depth B /3 below the mat) when subjected to loading of small magnitude ( Q = 500kN, corresponding to mobilized factor of safety SF = 2.31 against bearing capacity). However, decreasing the preconsolidation pressure of the soil (eroded overburden soil D = 5m, corresponding to OCR = 1.8 at a depth B /3 below the mat), the K s / K r distribution at the center and the edges of the mat foundation decreases as C c increases (Fig. 12 c), reflecting a decrease in overall (secant) shear stiffness. Notably, increasing C c by a factor of 6 (from 0.1 to 0.6) has minimal effect on the moment diagrams (Figs. 12 b and 12 d). This is attributed to the fact that despite the equivalent springs becoming substantially softer with increasing C c (Fig. 12 c), the edge-to-center ratio of K s values remains in the 1.6 to 1.8 range. 4.1.3 Effect of overconsolidation The K s / K r distribution depends heavily on the preconsolidation pressure ( p c ) of the soil, provided sufficient loading is applied to produce plastic yielding (Fig. 13 ). In this comparison, the preconsolidation pressure is represented by D , which is the thickness of the soil removed by erosion. Four D values are examined, namely 5, 10, 20 and 40m, which correspond to p c at the ground surface equal to 27, 54, 109 and 218 kPa, respectively (with OCR at depth B /3 below the mat being equal to 1.8, 2.5, 4.0, 7.1, respectively). Figure 13 indicates that an increase in p c generally leads to larger spring constants because of increased overall resistance to settlement resulting from an enhanced shear strength. However, it can be seen that, for highly overconsolidated clay ( D = 20m and 40m), the distributions are close to each other, especially in the case of the central distribution. This is because in these cases the clay remains mostly in an elastic state. For low Q values in particular (Figs. 13 a, and 13 c), the applied loading is not large enough to produce plastic yielding at the perimeter of the foundation, causing no significant change on the K s / K r distributions (both central and edge), as well as the corresponding bending moment diagrams. As a result, having stiff (high OCR ) clay generally leads to larger bending moments within the mat compared to the case of a softer, slightly overconsolidated clay. As previously discussed, an increase in the loading magnitude causes the K s / K r distribution in the D = 10m analysis to decrease at the edges. For a slightly overconsolidated clay ( D = 5m), the center and the edge K s / K r distributions decrease even more and converge, even for low loading magnitude ( Q = 500kN), due to the early development of plastic strains. As a result, soft clay (low OCR ) tends to produce a flatter K s / K r distribution subjecting the mat foundation to less dishing behavior, which results in lower peak (bottom fiber in tension) bending moments. It is interesting to note that, although the peak moments for Q equal to 500kN for slightly overconsolidated soil ( D = 5m) are the lowest and increase as OCR increases (Fig. 13 b), a reversal in the trend is observed for higher loads (e.g., Q = 1000kN) between the cases of D = 5m and D = 10m (Figs. 13 d and 13 f). This signifies that there is a second effect coming from a factor other than the enhancement of shear strength caused by an increased degree of overconsolidation. The second factor is the increase of elastic soil stiffness caused by the overconsolidation. The effect of the elastic stiffness on the bending moment distribution can be better visualized for a purely linear elastic analysis with different Young’s modulus (Fig. 14 ). The stiffer the soil (higher the Young’s modulus) is, the smaller the bending moments tend to become. Therefore, by decreasing the OCR , the soil exhibits smaller elastic stiffness, and the bending moments tend to increase, counterbalancing the reduction of bending moments caused by the early development of plastic shear strain in the case of slightly overconsolidated clay. 4.1.4 Effect of load eccentricity The effect of eccentricity on K s / K r can be seen in Fig. 15 a, which presents the results from analysis with the same mat and soil properties but with different loading type: symmetric (C1), single axis eccentricity (C2) and double axis eccentricity (C3). As observed for linear elastic soil (Loukidis and Tamiolakis 2017), the presence of load eccentricity rotates the K s / K r distribution towards the side that the eccentricity points to. This is particularly true for the distribution of the edge spring constants. Furthermore, as eccentricity rotates the Winkler spring distribution, the central bending moment diagram behaves similarly (Fig. 15 b). 4.2 Comparison with different methods of analysis In this section, we examine how much error is introduced in the calculation of the bending moments diagrams by a) assuming that the soil is purely linear elastic, i.e., ignoring the development of plastic strains and b) using simplified modulus of subgrade reaction distributions. Figure 16 compares the bending moments diagrams predicted using the K s distribution back-calculated from the elastoplastic cFEA (Modified Cam-Clay) to those produced i) using K s distribution back-calculated from linear elastic cFEA in which the soil Young’s modulus is assumed to be equal to E * , ii) assuming a uniform modulus of subgrade reaction ( k s ) across the mat, iii) using the 3-zone pseudo-coupled approach (Coduto 2001), and iv) using the Discrete Area Method (DAM) described in the methodology section (Section 3). For the cases presented in Fig. 16 , the mobilized safety factor ( SF ) against undrained bearing capacity ranges from 1.88 to 5.56. It should be noted that, in practice, shallow foundations are required to have an SF against undrained bearing capacity at least greater than 2. As expected, the analysis with the K s distribution from purely elastic cFEA predicts larger bending moments than that from the elastoplastic cFEA. For the examples shown in Fig. 16 , the elastic cFEA overestimates the peak bending moments at the column locations by 20–60%. The reason for this can be illustrated by comparing the respective K s distributions. Figure 17 compares the K s distributions for the cFEA analyses of Figs. 16 b and 16 d. The assumption of linear elasticity leads to higher edge-to-center K s ratio irrespective of the applied load magnitude and more variability in the K s distribution because of the absence of plastic yielding at the edges of the mat foundation. Regarding existing approaches used in practice, the use of a uniform k s distribution largely underestimates the peak positive (bottom fiber in tension) bending moments and overestimates the peak negative (top fiber in tension) bending moments. This is because the foundation has less tendency to bend with a concave upwards shape (dishing) due to the complete absence of higher spring constants at the mat perimeter than at its central region. Conversely, modeling the soil through the simple pseudo-coupled approach where k s decreases from the mat edges towards the center (Coduto 2001) produces a realistic dishing profile (concave upwards) and larger peak bending moments. For mobilized safety factor SF against undrained bearing capacity close to 2 (Figs. 16 a, and 16 b), the pseudo-coupled distribution predictions turn out to be very close to those from elastoplastic cFEA. This is because plastic deformations reduce the requirement for large K s at the mat edge. However, for higher mobilized SF values (greater than 3), which are much more frequently encountered in practice, the simple pseudo-coupled distribution does not produce accurate bending moments (Figs. 16 c, and 16 d), since it assumes that k s at the edges is only twice that at the center, while in fact it has to be much higher. ACI 318 − 19( 22 ) states that doubling the k s at the edges introduces a certain degree of spring coupling effects but is inadequate. On the contrary, the Discrete Area Method employed herein, which accounts for the inelastic consolidation settlements through Settle3 calculations, exhibits a superior performance. In all cases, employing the Discrete Area Method is a better alternative to be used in practice because it predicts similar moments to those produced by the elastoplastic finite element analysis. Finally, Fig. 18 plots the shear force ( V x ) diagrams along the mat centerline predicted by the different methods of analysis. Similar to the findings of Loukidis and Tamiolakis (2017) for the linear elastic case, the choice of K s (or k s ) distribution has practically very little effect on the shear force diagram. This is because static equilibrium dictates that the peak values of the diagram need to be exactly equal to the applied concentrated (column) loads, leaving little possibility for variation across a span among the different methods. Most importantly, the small discrepancies observed midway between two column loads in Fig. 18 are inconsequential with respect to the shear resistance requirements (i.e., punching failure checks). 5. CONCLUSIONS The appropriate Winkler spring constant ( K s ) distribution for mat foundation analysis is best captured using three-dimensional continuum finite element analysis (cFEA). In this study, the spring constants suitable for analysis of mats on clay for long-term (drained) conditions were back-calculated on the basis that the soil exhibits elastoplastic behavior following the Modified Cam-Clay model. The effects of various soil and mat parameters on the resulting K s distributions as well as on the bending moments were investigated. The cFEA-based predictions were compared to those from various existing modulus of subgrade reaction variation approaches encountered in foundation engineering practice. Based on the results of this study, the following main conclusions can be derived: ( 1 ) The development of plastic deformations affects the shape of the K s distribution by lowering the equivalent spring constant values at the edges of the mat foundation where most of the shearing occurs and the contact stresses are larger. Closer to the center of the foundation, the soil remains mostly elastic, unless superstructure loading is large enough to cause exceedance of the preconsolidation pressure in the clay close to the central area of the mat. ( 2 ) The shape of the K s distribution is insensitive to the clay deformability parameters, namely the compression index, the swelling index and the Poisson’s ratio. It is strongly affected by the degree of overconsolidation of the clay deposit and the magnitude of mat loading. ( 3 ) Ignoring the elastoplastic behavior of the clay leads to an overprediction of the spring constants at the mat edges, which, in turn, results in an overly conservative estimate of the peak bending moments. ( 4 ) A variant of the Discrete Area Method in which consolidation settlements are taken into account predicts bending moment distributions that are very similar to those based on elastoplastic cFEA. Hence, DAM comes out as a superior alternative for use in foundation engineering practice. Declarations Funding Declaration: No Funding. Competing Interests: The authors have no competing interests to declare that are relevant to the content of this article. Authors Contribution: The authors have contributed equally to this work. The first draft of the manuscript was written by Gaby Saad and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript. Financial/non-financial interest: The authors have no relevant financial or non-financial interests to disclose. Ethical approval: The authors declare that this manuscript have not been submitted to any other journal for simultaneous consideration. The work is original and not published elsewhere. Data Availability: Some or all data, models, or code that support the findings of this study are available from the corresponding author upon reasonable request. References Abdullah, W. S. (2008). New elastoplastic method for calculating the contact pressure distribution under rigid foundations. Jordan journal of civil engineering , 2 (1), 71-89. ACI (2019). 318-19; Building Code Requirements for Structural Concrete and Commentary. American Concrete Institute: Farmington Hills, MI, USA . Alzoaby, H., Saad, G. & Abou-Jaoude, G. Implementation of the discrete area method and its impact on the steel reinforcement of large mat foundations. Innov. Infrastruct. Solut. 10, 131 (2025). https://doi.org/10.1007/s41062-025-01926-x Ansys® Academic Research Mechanical, Release 2023 R2 Banavalkar, P.V. (1995). Mat Foundation and Its interaction with the Superstructure. ACI Special Publication, 152. CEN (2004). Eurocode 7 Geotechnical design. Part 1: General rules. EN 1997-1:2004, Brussels: European Committee for Standardization. Chow, Y. K. (1987). Vertical deformation of rigid foundations of arbitrary shape on layered soil media. International journal for numerical and analytical methods in geomechanics , 11 (1), 1-15. ‏ Coduto, D. P. (2001). Foundation design: principles and practices . Prentice Hall. Dutta, S. C., & Roy, R. (2002). A critical review on idealization and modeling for interaction among soil–foundation–structure system. Computers & structures , 80 (20-21), 1579-1594. Edgers, L., Alonge, J. L., & Sanayei, M. (2005). Modeling the effects of soil-structure interaction on a tall building bearing on a mat foundation. Civil Eng. Pract , 20 , 51-68. Estephan, R., Skeini, H., Malek, C., Mounir, M., Rahim, M.A., Abou-Jaoude, G., Yacoub, T. and Valiulin, D. (2022). A new design tool for shallow foundations offering enhanced accuracy, reduced design time, and cost savings. Proceedings of the 20 th International Conference on Soil Mechanics and Geotechnical Engineering, Sydney, Australia, 5169-5173. Horvath, J. (1995). Subgrade modeling for mat foundations: a review and critique of analytical methods. Special Publication , 152 , 117-160. Horvath, J. S., & Colasanti, R. J. (2011). Practical subgrade model for improved soil-structure interaction analysis: Model development. International Journal of Geomechanics , 11 (1), 59-64. Jaky, J. (1944). The coefficient of earth pressure at rest. Journal of the Society of Hungarian Architects and engineers, 355-358. Larkela, A., Mengelt, M., & Stapelfeldt, T. (2013). Determination of distribution of modulus of subgrade reaction. Proc., 18th Int. Conf. on Soil Mechanics and Geotechnical Engineering, French Society for Soil Mechanics and Geotechnical Engineering (CFMS), Paris, France, Leonidou, A. K. (2021). Spatial distribution of winkler spring stiffness for the structural analysis of embedded rectangular mat foundations . MSc Thesis, University of Cyprus. Cyprus. Loukidis, D., & Tamiolakis, G.-P. (2017). Spatial distribution of Winkler spring stiffness for rectangular mat foundation analysis. Engineering Structures , 153 , 443-459. . MathWorks (2022). MATLAB version: 9.13.0 (R2022b), Natick, Massachusetts: The MathWorks Inc. https://www.mathworks.com . Mayne, P. W., and Kulhawy, F. H. (1982). K 0 – OCR relationship in soil.” J. Geotech. Eng. Div., Am. Soc. Civ. Eng., 108(6), 851–872. Meyerhof, G. G. (1953). Some recent foundation research and its application to design. The Structural Engineer, 31(6), 151-167. Noorzaei, J., Viladkar, M., & Godbole, P. (1995). Elasto-plastic analysis for soil-structure interaction in framed structures. Computers & structures , 55 (5), 797-807. Pishilis, R. (2022). The influence of the elastoplastic behavior of the soil on the spatial distribution of the stiffness of the Winkler springs for mat foundations on clay. MSc Thesis, University of Cyprus. Cyprus. (in greek) Potts, D. M., Zdravković, L. (1999). Finite element analysis in geotechnical engineering: Theory (Vol. 1). Thomas Telford London. Poulos, H. G. (2016). Tall building foundations: design methods and applications. Innovative Infrastructure Solutions , 1 , 1-51. ‏ Roscoe, K. H., & Burland, J. B. (1968). On the generalized stress–strain behaviour of the ‘wet’clay. Pp. 535–609 in:(J. Heyman and FA Leckie, editors). In Cambridge University Press. Cambridge. Roscoe, K. H., Schofield, A., & Wroth, A. P. (1958). On the yielding of soils. Geotechnique , 8 (1), 22-53. Sadek, M., Mroueh, H., & Shahrour, I. (2010). Influence of nonlinearity on the stress distribution in the soil—application to road engineering problems. Journal of transportation engineering , 136 (1), 77-83. Sadek, M., & Shahrour, I. (2007). Use of the Boussinesq solution in geotechnical and road engineering: influence of plasticity. Comptes rendus. Mécanique , 335 (9-10), 516-520. Saeed, H. H. (2022). Effect of soil nonlinearity on analysis of raft foundation. Civil and Environmental Engineering , 18 (1), 292-300.[30] Salgado, R. (2008). The engineering of foundations . McGraw-Hill Europe. Schultze, E. (1961). Distribution of stress beneath a rigid foundation. In Proc., 5th Int. Conf. on Soil Mechanics and Foundation Engineering (pp. 807-813), Paris. Steinbrenner, W. (1934). Tafeln zur setzungsberechnung. Die Straße , 1 . Tamiolakis G-P. (2012). Study of the spatial distribution of the Winkler spring stiffness for the static analysis of mat foundations. MSc Thesis, University of Cyprus, Cyprus. (in greek). Teng, W.C. (1969). Foundation Design . Prentice-Hall. Ulrich, E. (1995). Introduction to the State-of-the-Art Mat Foundation Design and Construction. Special Publication , 152 , 1-12. Ulrich, E. J. (1991). Subgrade reaction in mat foundation design. Concrete International , 13 (4), 41-50. ‏ Viladkar, M., Godbole, P., & Noorzaei, J. (1991). Soil-structure interaction in plane frames using coupled finite-infinite elements. Computers & structures , 39 (5), 535-546. Viladkar, M., Noorzaei, J., & Godbole, P. (1994). Interactive analysis of a space frame-raft-soil system considering soil nonlinearity. Computers & structures , 51 (4), 343-356. Winkler, E. (1867). Die Lehre von Elastizitat und Festigkeit (The theory of elasticity and stiffness). H. Domenicus. Prague . Additional Declarations No competing interests reported. Supplementary Files APPENDIX.docx Cite Share Download PDF Status: Posted Version 1 posted You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-6456559","acceptedTermsAndConditions":true,"allowDirectSubmit":true,"archivedVersions":[],"articleType":"Research Article","associatedPublications":[],"authors":[{"id":455141605,"identity":"0339d732-d16e-4e39-8469-526a8786719c","order_by":0,"name":"Gaby Saad","email":"","orcid":"","institution":"Lebanese American University","correspondingAuthor":false,"prefix":"","firstName":"Gaby","middleName":"","lastName":"Saad","suffix":""},{"id":455141606,"identity":"23502359-1292-4ff5-a7f1-346a66fefee5","order_by":1,"name":"Grace Abou-Jaoude","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAABGklEQVRIiWNgGAWjYBCDBChmZuBnBwswMzCw8xCpRbIZpoWZsBaISoPDBLTIRx9++OEDg10ev0TCww8f26zljA8zP5NgqLBObGDmPYBNi+G5NGPJGQzJxZIzEpIlZ7alG5sdZjOTYDiTDtTCl4BVSw+DGdABzIkbbiQkSPO2HU7cdpjBTIIRyGhg5jHAroX9G/MfhnqQluTfQC31m5vZv0kw/sOtRZ6Hxwzo08MgLWkgWxIMmHmAtjTg1mLAw1Ms2WNwHEg8SLOccS7dcMZhnmKLhGPpxm04/CLfw77xw4+K6jx+9pzkGx/KrOX529s33vhQYy3bz96LNcQMwKJgB/AgmQlismFTD7KlAc5kx2rmKBgFo2AUjAIGAAJ1WOYlNF+cAAAAAElFTkSuQmCC","orcid":"","institution":"Lebanese American University","correspondingAuthor":true,"prefix":"","firstName":"Grace","middleName":"","lastName":"Abou-Jaoude","suffix":""},{"id":455141607,"identity":"af98eb72-8727-4704-ba90-def9363b094e","order_by":2,"name":"Dimitrios Loukidis","email":"","orcid":"","institution":"Lebanese American University","correspondingAuthor":false,"prefix":"","firstName":"Dimitrios","middleName":"","lastName":"Loukidis","suffix":""}],"badges":[],"createdAt":"2025-04-15 16:08:11","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-6456559/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-6456559/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":82625788,"identity":"ee714fff-ee63-49f5-9493-bb1c6d5bcae3","added_by":"auto","created_at":"2025-05-13 13:00:34","extension":"jpg","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":71189,"visible":true,"origin":"","legend":"\u003cp\u003eTypical mesh used in the finite element analyses\u003c/p\u003e","description":"","filename":"1.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/ec84cce735faea1df3f1fa45.jpg"},{"id":82624822,"identity":"0f21a221-b073-4ba9-aea3-7cd6365531b2","added_by":"auto","created_at":"2025-05-13 12:52:34","extension":"jpg","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":64419,"visible":true,"origin":"","legend":"\u003cp\u003eLoading Configurations\u003c/p\u003e","description":"","filename":"2.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/03a74e14da51e93bc40fef46.jpg"},{"id":82623678,"identity":"28828ecd-72df-4aff-96e6-337bbc6d2f99","added_by":"auto","created_at":"2025-05-13 12:44:34","extension":"jpg","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":94343,"visible":true,"origin":"","legend":"\u003cp\u003eComparison of distribution of mat vertical displacement from analyses using (a) Uniform k\u003csub\u003es\u003c/sub\u003e, (b) Pseudo-Coupled k\u003csub\u003es\u003c/sub\u003e and (c) elastoplastic cFEA K\u003csub\u003es\u003c/sub\u003e distribution\u003c/p\u003e","description":"","filename":"3.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/2f852ead2cb375b5ca6bfe8a.jpg"},{"id":82624820,"identity":"cfcfdd34-55e1-4724-bea5-0694b01cd8ca","added_by":"auto","created_at":"2025-05-13 12:52:34","extension":"jpg","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":53856,"visible":true,"origin":"","legend":"\u003cp\u003eExamples mat vertical displacement profiles for column load Q equal to (a) 500kN, and (b) 750kN\u003c/p\u003e","description":"","filename":"4.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/2ae9c1927dbb7c0d7f182213.jpg"},{"id":82625789,"identity":"7b31ebfb-7d9e-4889-89bf-dc47e027da72","added_by":"auto","created_at":"2025-05-13 13:00:34","extension":"jpg","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":153943,"visible":true,"origin":"","legend":"\u003cp\u003eBack-calculated K\u003csub\u003es\u003c/sub\u003e distribution from different analyses with C\u003csub\u003ec\u003c/sub\u003e=0.35, K\u003csub\u003e0\u003c/sub\u003e=1, v=0.25 and loading C1\u003c/p\u003e","description":"","filename":"5.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/fc466a1956a5b05898c25136.jpg"},{"id":82623679,"identity":"ac2e6fd0-cd22-452a-bf2a-ac2f1488f370","added_by":"auto","created_at":"2025-05-13 12:44:34","extension":"jpg","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":67325,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of critical friction angle on the distribution of (a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e and (b) mat centerline bending moment\u003c/p\u003e","description":"","filename":"6.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/e95fec0bea7e51112edffaea.jpg"},{"id":82623689,"identity":"e965ebd6-bc85-432b-a31a-ec492527ab09","added_by":"auto","created_at":"2025-05-13 12:44:34","extension":"jpg","order_by":7,"title":"Figure 7","display":"","copyAsset":false,"role":"figure","size":61849,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of mat thickness on the distribution of (a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e and (b) mat centerline bending moment\u003c/p\u003e","description":"","filename":"7.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/51c90835cf9aa9ed15602f1e.jpg"},{"id":82623694,"identity":"9f53c0af-521a-4be7-9362-687eda5b5b93","added_by":"auto","created_at":"2025-05-13 12:44:34","extension":"jpg","order_by":8,"title":"Figure 8","display":"","copyAsset":false,"role":"figure","size":80964,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of column load Q on the distribution of (a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e and (b) mat centerline bending moment\u003c/p\u003e","description":"","filename":"8.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/c42ef4ce31b7a197ae279cbb.jpg"},{"id":82623682,"identity":"51e272cc-a520-432c-9bb6-1669a84803b2","added_by":"auto","created_at":"2025-05-13 12:44:34","extension":"jpg","order_by":9,"title":"Figure 9","display":"","copyAsset":false,"role":"figure","size":69599,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of Κ\u003csub\u003e0\u003c/sub\u003e on the distribution of \u0026nbsp;(a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e and (b) mat centerline bending moment\u003c/p\u003e","description":"","filename":"9.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/968e2e1f44455a2087e7da59.jpg"},{"id":82623697,"identity":"93c76008-0d10-4076-b4f0-0545f82d8725","added_by":"auto","created_at":"2025-05-13 12:44:35","extension":"jpg","order_by":10,"title":"Figure 10","display":"","copyAsset":false,"role":"figure","size":55194,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of swelling index on the distribution of \u0026nbsp;(a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e and (b) mat centerline bending moment\u003c/p\u003e","description":"","filename":"10.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/20aad58ba9df811deca7a39a.jpg"},{"id":82623723,"identity":"649dc347-a5d2-46d4-81db-10e0d5d8f0bb","added_by":"auto","created_at":"2025-05-13 12:44:36","extension":"jpg","order_by":11,"title":"Figure 11","display":"","copyAsset":false,"role":"figure","size":59488,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of Poisson’s ratio on the distribution of \u0026nbsp;(a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e and (b) mat centerline bending moment\u003c/p\u003e","description":"","filename":"11.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/aecbb15dd2981c51d2915247.jpg"},{"id":82624826,"identity":"447a4d8b-23da-4b64-84b6-396af885c40f","added_by":"auto","created_at":"2025-05-13 12:52:34","extension":"jpg","order_by":12,"title":"Figure 12","display":"","copyAsset":false,"role":"figure","size":85801,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of C\u003csub\u003ec\u003c/sub\u003e on the distribution of \u0026nbsp;(a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e (left) and (b) mat centerline bending moment (right)\u0026nbsp;\u0026nbsp;\u0026nbsp;\u003c/p\u003e","description":"","filename":"12.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/28d43b2dbbcea2a2d5b5fc2d.jpg"},{"id":82623702,"identity":"01e403b1-f36e-4900-beec-ad5c9fe8699e","added_by":"auto","created_at":"2025-05-13 12:44:35","extension":"jpg","order_by":13,"title":"Figure 13","display":"","copyAsset":false,"role":"figure","size":185464,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of degree of overconsolidation on the distribution of \u0026nbsp;(a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e (left) and (b) mat centerline bending moment (right)\u003c/p\u003e","description":"","filename":"13.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/d1680039290b5ef2067a7664.jpg"},{"id":82624828,"identity":"47901f75-59f3-49d7-a7e5-e7c2fc597b38","added_by":"auto","created_at":"2025-05-13 12:52:35","extension":"jpg","order_by":14,"title":"Figure 14","display":"","copyAsset":false,"role":"figure","size":34264,"visible":true,"origin":"","legend":"\u003cp\u003eExample of influence of soil Young’s Modulus (E) on bending moments\u003c/p\u003e","description":"","filename":"14.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/dc23f29927ff9364ff3cb2b2.jpg"},{"id":82625791,"identity":"f48beb5c-fe23-43d5-a65c-1ff343a3132f","added_by":"auto","created_at":"2025-05-13 13:00:35","extension":"jpg","order_by":15,"title":"Figure 15","display":"","copyAsset":false,"role":"figure","size":63471,"visible":true,"origin":"","legend":"\u003cp\u003eEffect of load eccentricity on the distribution of \u0026nbsp;(a) K\u003csub\u003es\u003c/sub\u003e/ K\u003csub\u003er\u003c/sub\u003e and (b) mat centerline bending moment\u003c/p\u003e","description":"","filename":"15.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/b973791dc9b8c7ef44a16e83.jpg"},{"id":82623699,"identity":"8a2a5dbd-51f5-4102-84d0-e86839f0e5b6","added_by":"auto","created_at":"2025-05-13 12:44:35","extension":"jpg","order_by":16,"title":"Figure 16","display":"","copyAsset":false,"role":"figure","size":128615,"visible":true,"origin":"","legend":"\u003cp\u003eBending moment diagrams along mat centerline from various methods of subgrade analysis\u003c/p\u003e","description":"","filename":"16.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/36b4a90099ef9b80ea22ab5b.jpg"},{"id":82624830,"identity":"ad7b42af-c43f-4466-a8ba-bbe3cfac0b61","added_by":"auto","created_at":"2025-05-13 12:52:35","extension":"jpg","order_by":17,"title":"Figure 17","display":"","copyAsset":false,"role":"figure","size":66957,"visible":true,"origin":"","legend":"\u003cp\u003eComparison of K\u003csub\u003es\u003c/sub\u003e distribution back-calculated from elastoplastic cFEA and linear elastic cFEA using E=E*\u003c/p\u003e","description":"","filename":"17.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/3ad8ab76d9b9947c0923344d.jpg"},{"id":82623733,"identity":"bbac3ca6-792d-43b1-ad4a-1227a72bcce6","added_by":"auto","created_at":"2025-05-13 12:44:36","extension":"jpg","order_by":18,"title":"Figure 18","display":"","copyAsset":false,"role":"figure","size":79320,"visible":true,"origin":"","legend":"\u003cp\u003eShear force diagrams along mat centerline from various methods of analysis\u003c/p\u003e","description":"","filename":"18.jpg","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/12f613a65b33e5954c1f2a0a.jpg"},{"id":94988309,"identity":"3135d535-656b-48c2-bad4-52f340fa2a96","added_by":"auto","created_at":"2025-11-03 07:08:34","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":2483352,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/df86f2df-d482-4558-b7b6-44563534aeac.pdf"},{"id":82623677,"identity":"25974f04-59bb-4925-9922-2e9903b89794","added_by":"auto","created_at":"2025-05-13 12:44:34","extension":"docx","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":47158,"visible":true,"origin":"","legend":"","description":"","filename":"APPENDIX.docx","url":"https://assets-eu.researchsquare.com/files/rs-6456559/v1/a84153ea0e198cba71bbd61b.docx"}],"financialInterests":"No competing interests reported.","formattedTitle":"\u003cp\u003eWinkler Spring Stiffness Distribution for the Structural Analysis of Mat Foundations in Clay Under Drained Conditions\u003c/p\u003e","fulltext":[{"header":"1. INTRODUCTION","content":"\u003cp\u003eThe behavior of structures is affected by the properties and conditions of the foundation soil, while at the same time the structure and its loading influence the behavior of the soil. Understanding this interplay between soil and structure is fundamental in design to ensure the stability, performance, and safety of the structures. Soil-structure interaction is particularly important for the design of mat (or raft) foundations of heavy structures, where minimum reinforcement often does not govern the structural design. Consequently, engineers are forced to adopt sophisticated methods for analysis and design. This is particularly evident in geotechnical engineering practice where there is progressive shift from traditional methods to high-level and computationally expensive numerical analysis to ensure precision and reliability involving the design of foundations of high-rise buildings (Poulos 2016). On the other hand, for more common projects, simplified subgrade models are still used in the context of structural analysis software.\u003c/p\u003e \u003cp\u003eEngineers favor mat foundations when dealing with large loads that must be transferred to weak ground. They are preferred because they can resist differential settlements, heave and bending moments more effectively than other types of shallow foundations. Mat design is mainly based on calculations of settlements, shear forces, and bending moments to determine the required mat thickness and steel reinforcement. Several subgrade models have been developed to study the interaction between the structure and the underlying soil. The soil-structure interaction of a beam foundation modeled as a structural element laying on a bed of independent linear springs representing the soil was proposed by Winkler (1867). The stiffness of the springs depends directly on the modulus of subgrade reaction \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e, defined as:\u003cdiv id=\"Equ1\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ1\" name=\"EquationSource\"\u003e\n$$\\:{k}_{s}=\\frac{q}{\\delta\\:}$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e1\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\delta\\:\\)\u003c/span\u003e\u003c/span\u003e is the settlement generated by the applied pressure \u003cem\u003eq\u003c/em\u003e. The spring constant (\u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e) is obtained by multiplying \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e by the tributary area of each node. According to Horvath (1995), \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e is the most important parameter considered in the design of mat foundations. However, its calculation is not simple because it is not a fundamental property of the soil. It is a function of the stiffness parameters of the soil, as well as the geometrical characteristics of the foundation (Coduto 2001).\u003c/p\u003e \u003cp\u003eIt has long been known that the simple bed of springs with constant \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e throughout the mat (uncoupled springs approach) cannot capture the actual deformation at the soil-mat contact surface (Horvath 1995). The coupled springs method represents a step up from the original Winkler spring method, where the addition of other mechanical elements (such as a stretched membrane under tension to account for the soil deformation continuity or a shear layer of defined shear modulus to account for the effects of shearing) simulate the soil response more realistically and allow the prediction of the correct settlement profile. However, its mathematical formulation is complex and requires careful considerations regarding the boundary conditions. The pseudo-coupled method is an alternative approach that relies on assuming a variable spatial distribution of the modulus of subgrade reaction \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e across the mat. This method attempts to overcome the lack of coupling in the original Winkler spring method, while retaining its simplicity. Given that many structural design software already employ the Winkler spring method, accommodating the pseudo-coupled method is feasible and represents the most practical approach to date. The pseudo-coupled method was further improved by adding an iterative procedure that relies on the analytical Boussinesq solution (Ulrich 1995). This improved method is known as the Discrete Area Method (DAM). Although Ulrich (1995) advocated the use of DAM in design and highlighted its effectiveness in producing accurate results, current design codes do not enforce or make reference to this procedure.\u003c/p\u003e \u003cp\u003eMore recently, using 3D finite element modeling, Loukidis and Tamiolakis (2017) and Leonidou (2021) performed parametric analyses of mat foundations of small planar dimensions placed on linear elastic soil to back-calculate the spatial variation of the spring constants that is needed to achieve an accurate mat structural analysis. They concluded that when the soil-structure interaction is properly simulated using a variable \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution, the peak moments corresponding to bottom fiber in tension of the mat foundation increase and those corresponding to top fiber in tension decrease compared to those predicted using the uncoupled Winkler approach.\u003c/p\u003e \u003cp\u003eAlzoaby et al. (2022) applied the DAM, using commercially available software SAFE and Settle3, to study the effects of the variable \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution on the steel reinforcement of large mat foundations underlain by soft to very stiff linearly elastic soil. This work showed that implementing a variable \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution in the design of mat foundations drastically affects the moment distribution of mat foundations and thus steel reinforcement leading to a more accurate and realistic design than using uncoupled springs or other pseudo-coupled approaches.\u003c/p\u003e \u003cp\u003eThe aforementioned studies have highlighted the importance of varying the modulus of subgrade reaction under mat foundations. However, they all revolved around characterizing the soil as a linear elastic material. The focus of the present study is the effect of long-term plastic deformation on the distribution of the equivalent spring constants for mat foundations resting on clay under long-term fully drained conditions, i.e., end-of-consolidation state. Eurocode EN1997-1:2004 (CEN 2004) explicitly states as a basic requirement for geotechnical design that both short-term and long-term conditions shall be considered as design situations (\u0026sect;\u0026nbsp;2.2(\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e)). This concerns also the design against structural failure of the foundation (\u0026sect;\u0026nbsp;2.4.7.1, \u0026sect;\u0026nbsp;6.8(\u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e)). The nonlinear behavior of the soil is expected to cause a redistribution of the contact pressure under the foundation, which in turn will affect the bending moment distribution. This is attributed to the fact that the soil in the vicinity of the edges of the mat foundation is subjected to much more shearing than at the mat center and, as a result, the soil near the edges is prone to plastic yielding early during the mat loading process.\u003c/p\u003e \u003cp\u003eThe aim of the present study is to explore the influence of the long-term (consolidation) settlements on the bending moments developing in mat foundations resting on clay and propose ways to take into account their effect in design practice in the context of the pseudo-coupled approach. To achieve this objective, three-dimensional finite element analyses of mat-soil interaction are conducted using the software ANSYS. The soil is simulated as an elastoplastic material using the Modified Cam-Clay constitutive model (Roscoe \u0026amp; Burland 1965). The Modified Cam-Clay model can simulate the nonlinear hardening of normally and slightly overconsolidated clays and the softening of strongly overconsolidated clay (under drained conditions). Most importantly, it can capture accurately the development of plastic volumetric strains due to exceedance of the preconsolidation pressure. A parametric study is performed varying the soil properties, foundation characteristics, and loading magnitude. The calculated settlements at the nodes are used to back-calculate the equivalent linear spring constants \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e. The bending moment and shear force diagrams obtained are compared to those generated from three common methods of assigning modulus of subgrade reaction under the mat.\u003c/p\u003e"},{"header":"2. KNOWLEDGE BACKGROUND","content":"\u003cp\u003eSoils does not exhibit a purely linear elastic behavior for the entire range of loading prior to failure and, beyond a certain shear stress level, they exhibit nonlinearly hardening. Moreover, unlike structural materials, this nonlinear behavior actually initiates at early stages of shearing. Hence, when analyzing structures that induce high levels of stresses in the soil, as in the case of multistorey buildings, it is quite inaccurate to consider solely linear elastic behavior for the soil and it is reasonable to expect that this affects negatively the design of the mat and the structure. Researchers have investigated the role of nonlinear soil behavior in the soil-foundation interaction using two types of material models: nonlinear elastic models and elastoplastic models.\u003c/p\u003e \u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003e2.1 Nonlinear elastic soil behavior\u003c/h2\u003e \u003cp\u003eA nonlinear elastic model describes a nonlinear behavior of the soil in the \u0026ldquo;elastic zone\u0026rdquo; prior to failure by assuming that the elastic modulus is a decreasing function of the shear strains and an increasing function of the mean effective stress. Though an improvement compared to the purely linear elastic model, it still limits the representation of the soil behavior by not considering the effects of dilatancy and plastic (irreversible) deformations.\u003c/p\u003e \u003cp\u003eViladkar et al. (1991) demonstrated the significant influence of soil\u0026rsquo;s nonlinear behavior on foundation response and highlighted the necessity of nonlinear analysis techniques for accurately predicting settlement profiles, shear forces, and bending moments in rigid combined footings. The soil compressibility affects the behavior and design of foundations as well as the distribution of contact pressure under the footing, which influences the equivalent \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution. Moreover, it was observed that the shape of the contact pressure varies based on the soil properties and the magnitude of the loading. ACI 318\u0026thinsp;\u0026minus;\u0026thinsp;19(\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e) states that the distribution of soil pressure under a mat foundation should be consistent with the properties of the soil and is highly dependent on the soil type. Teng (1962) emphasized that the calculation of the distribution of the modulus of subgrade is different for clays and sand.\u003c/p\u003e \u003cp\u003eSimilar outcome is observed in analyses of mat foundations, in which a redistribution of contact pressure is obtained when considering the nonlinear elastic behavior of soil (Viladkar 1994). However, a nonlinear elastic analysis produces much higher settlements and therefore it should be limited to only small loading magnitudes, since it is unreliable at high loading magnitudes (Noorzai et al. 1995). The soil could be idealized as a nonlinear elastic material provided that the loading applied does not produce extensive regions of plastic straining in the soil mass, otherwise it should be modeled using an elastoplastic constitutive mode (Dutta and Roy 2002).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec4\" class=\"Section2\"\u003e \u003ch2\u003e2.2 Elastic - perfectly plastic behavior\u003c/h2\u003e \u003cp\u003eThe elastic \u0026ndash; perfectly plastic model assumes a linear elastic behavior prior to failure and a perfectly plastic behavior after that, which is governed by a failure criterion (usually the Mohr-Coulomb criterion). As such, it ignores post-failure strain-softening effects and pre-failure nonlinear hardening behavior, as well as plastic volumetric deformations due to consolidation and creep.\u003c/p\u003e \u003cp\u003ePlastic yielding occurs first at the edges of the foundation and results in the redistribution of the stresses under it. Noorzai et al. (1995) illustrated that, for combined footings, elastoplastic effects become prominent in the case of heavy structures. Their results show that the settlement profiles generated using elastoplastic soil model are almost identical to those obtained assuming linear elastic soil behavior under low loading conditions and differ at high loading conditions due to the effects of plastic yielding which occurs mainly at the edges of the mat foundation. This is explained by the observation that as loads increase, plastic deformations spread in the soil towards the center of the foundation (Noorzai et al. 1995; Abdullah 2008). Moreover, the combined footing deformed in the form of a \u0026ldquo;cup\u0026rdquo; shape, having the maximum settlement below heavily loaded areas only when the soil adheres to a yield criterion.\u003c/p\u003e \u003cp\u003eCalculating the stress distribution using the Boussinesq theory, where the soil is treated as a linear elastic homogenous isotropic half-space, underestimates the vertical stresses underneath the foundation down to a depth approximately twice the foundation width, which is an area that is contributing significantly to the settlement of the foundation (Sadek and Shahrour 2007). This discrepancy is attributed to the development of plastic shear deformations that lead to a reduced attenuation of vertical and horizontal stresses within the soil mass under loading (Sadek et al. 2010).\u003c/p\u003e \u003cp\u003eAbdullah (2008) argued that, once soil yielding occurs at and near the foundation\u0026rsquo;s edges, the assumption of uniform contact pressure used in practice is inaccurate as it may lead to unsafe mat foundation design. Recently, Pishilis (2022) performed finite element analyses using the computer program Abaqus to derive the distribution of equivalent Winkler spring constants considering elastoplastic behavior under constant volume (undrained clay behavior). Unlike Loukidis and Tamiolakis (2017), where they explain that for elastic soil analysis the applied loading magnitude does not affect the shape of the \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution, the study of Pishilis (2022) showed that the shape of the distribution varies with the applied loading magnitude. Accounting for soil yielding in design increases the accuracy in design significantly (Larkela et al. 2013; Pishilis 2022).\u003c/p\u003e \u003cp\u003eUlrich (1991) and Banavalkar (1995) highlighted the crucial role of both immediate and long-term subgrade responses in the design of mat foundations. They pointed out that the determination of the spring constants for the design of mat foundations should reflect the time-dependent behavior of the soil, which in the case of clays is governed by the phenomenon of consolidation. Despite this assertion and the Eurocode 7 previsions, the long-term conditions for foundations in clays are largely overlooked in the structural design practice.\u003c/p\u003e \u003c/div\u003e"},{"header":"3. METHODOLOGY","content":"\u003cdiv id=\"Sec6\" class=\"Section2\"\u003e \u003ch2\u003e3.1 Continuum finite element analysis\u003c/h2\u003e \u003cp\u003eThree-dimensional finite element analyses (FEA) using solid elements for the soil (continuum FEA or cFEA) were performed using ANSYS to derive the distribution of the equivalent Winkler spring constants (\u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e) across flexible mat foundations resting on clay under drained (end-of-consolidation) conditions. The settlement profiles were extracted from these analyses at each mat node and transferred into a MATLAB code that back-calculates the equivalent \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution (Tamiolakis 2012). The soil domain was discretized using 8-noded (first-order) hexahedral elements (type Solid185) having 3 translational degrees of freedom per node, with full Gauss integration (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003e). The mat foundation, placed at the center of the free surface of the soil domain, was modeled with 4-noded shell elements (type Shell 181) having 6 degrees of freedom per node (3 translational and 3 rotational degrees of freedom). The use of first-order elements eliminates the complexities of deriving an equivalent spring constant distribution associated with higher order elements, which include mid-side nodes in addition to corner nodes. As it was observed that the differences in results between using reduced and full integration shell elements are insignificant, reduced Gauss integration was adopted to reduce the computational cost (Loukidis and Tamiolakis 2017). The mat was meshed with square elements 0.5m wide, a size commonly used in engineering practice. The soil elements directly underneath the foundation were cubical with the same edge length (0.5m) as the mat shell element in order to have coincidence of mat and soil nodes at their contact. The same element mesh size was extended to an area twice the foundation\u0026rsquo;s horizontal dimensions. Beyond that region the element size was allowed to increase to reduce computational time. The interaction between the soil and the foundation was defined as \u0026ldquo;shared topology\u0026rdquo;. This configuration forces the mat and soil nodes to be rigidly connected to each other with respect to the translational degrees of freedom and prevents any possibility of slippage or separation between them.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe lateral boundaries of the soil domain were placed at a distance from the mat center that was 10 times the width \u003cem\u003eB\u003c/em\u003e of the mat, and the thickness \u003cem\u003eH\u003c/em\u003e of the soil domain was set to be 10 times the \u003cem\u003eB\u003c/em\u003e, which is sufficient for the boundaries to have minimal effect on the numerical results (Loukidis and Tamiolakis 2017). Each finite element parametric analysis comprised three phases. First, a geostatic stress field was established while applying gravity to the soil domain. Second, the self-weight on the mat was applied. Third, the column loads were applied on the mat as concentrated vertical forces.\u003c/p\u003e \u003cp\u003eThe mat material (reinforced concrete) was modeled as linear elastic, having Young\u0026rsquo;s modulus 32 GPa, Poisson\u0026rsquo;s ratio 0.2 and unit weight 24.5 kN/m\u003csup\u003e3\u003c/sup\u003e (density of 2500 kg/m\u003csup\u003e3\u003c/sup\u003e). The constitutive model for the clayey soil was the built-in version of Cam-Clay model in ANSYS, which can be rendered identical to the Modified Cam-Clay (MCC) model proposed by Roscoe \u0026amp; Burland (1968) through proper selection of the values for the model input parameters. To ensure that the constitutive model formulation and parameters have been correctly interpreted, the results of a single element analysis performed in ANSYS simulating undrained triaxial compression test were compared with MCC predictions presented in Potts and Zdravkovic (1999).\u003c/p\u003e \u003cp\u003eTo check the adequacy of the finite element mesh and boundary conditions, trial FE analyses were performed assuming a linear elastic soil and a very stiff mat and the resulting settlements were compared to a well-established analytical solution that considers the mat as a rigid foundation, namely the Steinbrenner (1934) equation coupled with the Chow (1987) rigidity factor. The FEA results are close to those obtained from the analytical solution, with an average error of 2.44%.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec7\" class=\"Section2\"\u003e \u003ch2\u003e3.2 Parametric analysis\u003c/h2\u003e \u003cp\u003eA parametric study was conducted to investigate the effect of plasticity under long-term conditions (end-of-consolidation) on the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution, as well as the developing bending moments, which control the design of the steel reinforcement. The mat foundations considered herein have planar dimension 20m \u0026times; 20m.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe arrangements of column loading is shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e. Each mat carries 25 equally-spaced column loads, with the column spacing (\u003cem\u003eS\u003c/em\u003e\u003csub\u003ecolumn\u003c/sub\u003e) being 5m (typical for multistorey residential buildings). Regarding the column relative magnitudes, most of the finite element analyses were performed for a typical symmetric loading configuration C1, where the columns transfer a load according to their tributary area. Configurations C2 and C3 induce single axis and double axis load eccentricity respectively (i.e., the resultant of the column loads does not pass through the center of the mat), subjecting the foundation to rotation. The thickness of the mat ranges from 0.5 m to 1m. These values satisfy punching shear criteria for the load magnitudes considered in the parametric analyses.\u003c/p\u003e \u003cp\u003eThe quantify the rigidity of the mat foundations, the relative stiffness factor was determined according to the following equation (Meyerhof 1953):\u003cdiv id=\"Equ2\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ2\" name=\"EquationSource\"\u003e\n$$\\:{R}_{\\text{s}}=\\frac{{E}_{c}{d}^{3}}{12E{B}^{3}}$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e2\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere \u003cem\u003eE\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e is modulus of elasticity of the mat concrete, \u003cem\u003ed\u003c/em\u003e is the thickness of the mat, \u003cem\u003eE\u003c/em\u003e is the modulus of elasticity of the soil and \u003cem\u003eB\u003c/em\u003e is the mat width. A foundation is characterized flexible if \u003cem\u003eR\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e \u0026lt; 0.5.\u003c/p\u003e \u003cp\u003eUnlike a linear elastic soil, the back-calculated spring constants for an elastoplastic soil depend on the loading magnitude \u003cem\u003eQ\u003c/em\u003e. Loading from inner columns (i.e., those that do not lie along the perimeter of the mat) ranges from 500kN (approximately 2 storeys) to 4000kN (approximately 16 storeys). The critical friction angle \u003cem\u003eφ\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e of the clayey soil ranges from 25\u0026deg; to 35\u0026deg;, the Poisson\u0026rsquo;s ratio \u003cem\u003eν\u003c/em\u003e from 0.15 to 0.35, the compression index \u003cem\u003eC\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e from 0.1 to 0.6 (i.e. MCC parameter \u003cem\u003eλ\u003c/em\u003e from 0.043 to 0.26), the swelling index \u003cem\u003eC\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e from 0.0175 to 0.07 (i.e. MCC parameter \u003cem\u003eκ\u003c/em\u003e from 0.0076 to 0.3), and the at-rest lateral earth pressure \u003cem\u003eK\u003c/em\u003e\u003csub\u003e0\u003c/sub\u003e from 0.5 to 1.\u003c/p\u003e \u003cp\u003eThe most common cause of overconsolidation in nature is the removal of overburden due to erosion, in which case the preconsolidation pressure \u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e is expected to increase practically linearly with depth, while the initial void ratio \u003cem\u003ee\u003c/em\u003e\u003csub\u003eini\u003c/sub\u003e would decrease nonlinearly with depth. Details on the determination of the \u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e and \u003cem\u003ee\u003c/em\u003e\u003csub\u003eini\u003c/sub\u003e profiles as a function of the assumed thickness \u003cem\u003eD\u003c/em\u003e of eroded overburden is given in the Appendix. These profiles were applied as initial conditions to the ANSYS finite element model by prescribing a user-defined pseudo-temperature distribution with depth and making \u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e and \u003cem\u003ee\u003c/em\u003e\u003csub\u003eini\u003c/sub\u003e dependent to it. The \u003cem\u003eD\u003c/em\u003e values considered herein range from 5 m to 40 m, corresponding to vertical preconsolidation stress at the ground surface from 41 kPa to 326 kPa.\u003c/p\u003e \u003cp\u003eThe geostatic stresses were established in the model using the buoyant unit weight \u003cem\u003eγ\u003c/em\u003e\u0026prime; (assuming that the water table is at the free surface). Given that the void ratio decreases with depth, \u003cem\u003eγ\u003c/em\u003e\u0026prime; ought to be an increasing function of depth. However, because ANSYS does not allow the unit weight to be dependent on the pseudo-temperature, \u003cem\u003eγ\u003c/em\u003e\u0026prime; was assumed constant with depth, equal to an average representative value of 8.15 kN/m\u003csup\u003e3\u003c/sup\u003e.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec8\" class=\"Section2\"\u003e \u003ch2\u003e3.3 Determination of equivalent spring constant distribution\u003c/h2\u003e \u003cp\u003eTo determine the equivalent Winkler spring constants from the results of the cFEA, the vertical displacements of all the nodes of the mat from each analysis were exported from ANSYS and introduced as input in the MATLAB code developed by Tamiolakis (2012). The code was modified for the purposes of the present study to account for the stiffness matrix of the shell elements used in ANSYS and to account for the foundation self-weight. This algorithm solves inversely the problem of a mat foundation resting on Winkler springs to produce as output the equivalent spring constants (\u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e) and is based on the following system of static equilibrium equations:\u003cdiv id=\"Equ3\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ3\" name=\"EquationSource\"\u003e\n$$\\:\\left(\\left[{\\text{K}}_{\\text{m}\\text{a}\\text{t}}\\right]+\\left[{\\text{K}}_{\\text{s}\\text{p}\\text{r}\\text{i}\\text{n}\\text{g}\\text{s}}\\right]\\right)\\left[\\text{U}\\right]=\\left[{\\text{F}}_{\\text{m}\\text{a}\\text{t}\\:\\text{s}\\text{e}\\text{l}\\text{f}-\\text{w}\\text{e}\\text{i}\\text{g}\\text{h}\\text{t}}+{\\text{F}}_{\\text{c}\\text{o}\\text{l}\\text{u}\\text{m}\\text{n}\\:\\text{l}\\text{o}\\text{a}\\text{d}\\text{s}}\\right]\\:$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e3\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere \u003cem\u003e[K\u003c/em\u003e\u003csub\u003e\u003cem\u003emat\u003c/em\u003e\u003c/sub\u003e\u003cem\u003e]\u003c/em\u003e is the assembled stiffness matrix of the mat foundation elements, \u003cem\u003e[K\u003c/em\u003e\u003csub\u003e\u003cem\u003esprings\u003c/em\u003e\u003c/sub\u003e\u003cem\u003e]\u003c/em\u003e is the assembled stiffness matrix of the Winkler springs below the mat, \u003cem\u003e[U]\u003c/em\u003e is the matrix of nodal vertical displacements and rotations, \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{F}_{mat\\:self-weight}\\)\u003c/span\u003e\u003c/span\u003e is the matrix of the equivalent nodal loads due to the mat self-weight, and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{F}_{column\\:loads}\\)\u003c/span\u003e\u003c/span\u003e is the matrix of the applied column loads. Note that \u003cem\u003e[K\u003c/em\u003e\u003csub\u003e\u003cem\u003esprings\u003c/em\u003e\u003c/sub\u003e\u003cem\u003e]\u003c/em\u003e is a diagonal matrix because each Winkler springs is affecting only the vertical nodal displacement of the node to which it is connected.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec9\" class=\"Section2\"\u003e \u003ch2\u003e3.4 Determination of representative Young\u0026rsquo;s modulus\u003c/h2\u003e \u003cp\u003eTo facilitate comparisons of the distribution of equivalent spring constants and bending moments from an elastoplastic analysis to those generated by a purely linear elastic analysis, a representative value of the Young\u0026rsquo;s modulus (\u003cem\u003eE\u003c/em\u003e\u003csup\u003e*\u003c/sup\u003e) was calculated at a depth below the mat center equal to one-third the mat width (\u003cem\u003eB\u003c/em\u003e/3) at the end of the elastoplastic analysis (after full application of the mat loading). This calculation was necessary because, in the Modified Cam-Clay model, the soil Young\u0026rsquo;s modulus \u003cem\u003eE\u003c/em\u003e is not an input model parameter and is actually a variable that depends on the current values of mean effective stress \u003cem\u003ep\u003c/em\u003e\u0026prime; and void ratio \u003cem\u003ee\u003c/em\u003e:\u003cdiv id=\"Equ4\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ4\" name=\"EquationSource\"\u003e\n$$\\:E=3\\left(1-2v\\right)\\frac{\\left(1+e\\right){p}^{{\\prime\\:}}}{\\kappa\\:}\\approx\\:6.9(1-2v)\\frac{(1+e)p{\\prime\\:}}{{C}_{s}}$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e4\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eAs such, the soil Young\u0026rsquo;s modulus varies spatially, and as mat loading is applied. Moreover, the soil near the edges undergoes plastic yielding (where shear straining is dominant). Plastic deformations may develop also near the mat center if the preconsolidation pressure \u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e is approached during loading, otherwise the soil in that region remains largely in an elastic state. The selection of this location for the reference point for the \u003cem\u003eE\u003c/em\u003e\u003csup\u003e*\u003c/sup\u003e determination (depth \u003cem\u003eB\u003c/em\u003e/3 below the mat center) was made using results of trial analyses and based on the rational that, given that in most analyses, the soil under the mat center remains in an elastic state (i.e. the stress state lies inside the Modified Cam-Clay yield surface), the spring constant at the center of the foundation from a linear elastic analysis should generally align with that produced at the center of the corresponding elastoplastic analysis. The \u003cem\u003eE\u003c/em\u003e\u003csup\u003e*\u003c/sup\u003e values among the analyses of the parametric study range from 7.6 MPa to 30.8 MPa, leading to relative stiffness factor \u003cem\u003eR\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e in the 0.011\u0026ndash;0.044 range.\u003c/p\u003e \u003cp\u003eThe \u003cem\u003eE\u003c/em\u003e\u003csup\u003e*\u003c/sup\u003e is used for calculating a reference spring constant value, \u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e, with respect to which the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distributions are normalized. This facilitates comparison between the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distributions from the various parametric analyses performed and discerning the effect of each problem parameter on the shape of the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution. The \u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e corresponds to the spring constant that would be applicable to a perfectly rigid foundation resting on linear elastic soil having Young\u0026rsquo;s modulus \u003cem\u003eE\u003c/em\u003e equal to \u003cem\u003eE\u003c/em\u003e\u003csup\u003e*\u003c/sup\u003e and is calculated as follows:\u003cdiv id=\"Equ5\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ5\" name=\"EquationSource\"\u003e\n$$\\:{K}_{\\text{r}}={k}_{\\text{r}}{A}_{\\text{i}\\text{n}}$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e5\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere \u003cem\u003eA\u003c/em\u003e\u003csub\u003ein\u003c/sub\u003e is the tributary area of the inner mat nodes (0.25 m\u003csup\u003e2\u003c/sup\u003e) and \u003cem\u003ek\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e is the modulus of subgrade reaction for rigid foundation:\u003cdiv id=\"Equ6\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ6\" name=\"EquationSource\"\u003e\n$$\\:{k}_{\\text{r}}=\\frac{{E}^{*}}{{C}_{f}\\left(1-{v}^{2}\\right)B}$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e6\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewith the \u003cem\u003eC\u003c/em\u003e\u003csub\u003e\u003cem\u003ef\u003c/em\u003e\u003c/sub\u003e factor proposed by Loukidis and Tamiolakis (2017), which introduces the effects of the mat aspect ratio (\u003cem\u003eL\u003c/em\u003e/\u003cem\u003eB\u003c/em\u003e) and the deformable soil layer thickness \u003cem\u003eH\u003c/em\u003e:\u003c/p\u003e \u003cp\u003e \u003cdiv class=\"gridtable\"\u003e\u003ctable float=\"No\" id=\"Taba\" border=\"1\"\u003e \u003ccolgroup cols=\"4\"\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c1\" colnum=\"1\"\u003e\u003c/div\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c2\" colnum=\"2\"\u003e\u003c/div\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c3\" colnum=\"3\"\u003e\u003c/div\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c4\" colnum=\"4\"\u003e\u003c/div\u003e \u003ctbody\u003e \u003ctr\u003e \u003ctd align=\"left\" colspan=\"2\" nameend=\"c2\" namest=\"c1\"\u003e \u003cp\u003e\u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{C}_{f}=\\frac{0.85\\:{\\left(\\frac{L}{B}\\right)}^{0.45}}{{\\left[1+0.1(2+\\frac{L}{B})\\frac{B}{H}\\right]}^{1+{e}^{5{v}^{3}}}}\\)\u003c/span\u003e\u003c/span\u003e\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colspan=\"2\" nameend=\"c4\" namest=\"c3\"\u003e \u003cp\u003e(\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e)\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003c/tbody\u003e \u003c/colgroup\u003e \u003c/table\u003e\u003c/div\u003e \u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec10\" class=\"Section2\"\u003e \u003ch2\u003e3.5 Discrete Area Method\u003c/h2\u003e \u003cp\u003eThe Discrete Area Method (DAM) using a combination of the computer programs CSI SAFE and Settle3 was first applied by Alzoaby et al. (2021) for elastic soils. Herein, we expand this computational approach by performing the settlement calculations in Settle3 using the built-in nonlinear 1-D consolidation model, where parameters such swelling index \u003cem\u003eC\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e, compression index \u003cem\u003eC\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e and preconsolidation pressure \u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e are inputted directly in the software, in order to make consistent comparisons with the cFEA-based results. Nonetheless, it has to be pointed out that, although a nonlinear model is used in Settle3, the program still calculates the stress changes induced by the application of the foundation loading using the Boussinesq solution, which stems from the theory of elasticity. Hence, due to this inconsistency, the resulting settlements, which in the case of a linear elastic soil would constitute exact solutions, in the present application of DAM should be considered approximate.\u003c/p\u003e \u003cp\u003eThe DAM is an iterative approach. First the mat-on-springs problem is solved in CSI SAFE using an assumed distribution (initial guess) of modulus of subgrade reaction. The computed contact pressures are then inputted in Settle3 and new values of modulus of subgrade reaction at each mat node are obtained by dividing the contact pressures by the corresponding settlement values at each mat node. The new distribution of the modulus of subgrade reaction is inputted back in CSI SAFE to recalculate soil pressures under the mat foundation. The iterative procedure stops when the settlement estimates from CSI SAFE and Settle3 converge. Convergence is considered to be achieved when the relative difference between the settlements calculated by CSI SAFE and Settle3 does not exceed 5% for every mat node:\u003cdiv id=\"Equ7\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ7\" name=\"EquationSource\"\u003e\n$$\\:{err}_{U}=\\left|\\frac{{\\left(U\\right)}_{CSI\\:SAFE}-\\:{\\:\\left(U\\right)}_{Settle3}}{{\\left(U\\right)}_{CSI\\:SAFE}}\\right|\\times\\:100\\:$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e8\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eUsually, an average of 11 iterations is sufficient to achieve convergence.\u003c/p\u003e \u003c/div\u003e"},{"header":"4. RESULTS AND DISCUSSION","content":"\u003cp\u003eFigure\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e compares the vertical displacement distribution across the mat foundation from cFEA with those obtained using i) a uniform \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution (i.e., uncoupled approach) and ii) the pseudo-coupled approach of Coduto (2001). In the latter, the mat is divided into 3 concentric zones, each with its own value of \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e. The inner zone has half the mat dimensions (\u003cem\u003eB\u003c/em\u003e/2 \u0026times; \u003cem\u003eB\u003c/em\u003e/2), while the intermediate and outer zones have thickness equal to \u003cem\u003eB\u003c/em\u003e/8 (Loukidis and Tamiolakis 2017). The modulus of subgrade reaction progressively increases from the central to the outermost zone, such that the outermost zone has a \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e value twice as much as that in the inner zone. The uniform \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e value and the average value of \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e across the mat for the pseudo-coupled approach are set equal to \u003cem\u003ek\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e in Eq.\u0026nbsp;(\u003cspan refid=\"Equ6\" class=\"InternalRef\"\u003e6\u003c/span\u003e).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eIt can be seen that the choice of the \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution strongly affects the displacement profile. The uniform \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution leads to less variation in the vertical displacements and, thus, smaller differential settlements. Most importantly, it results in a deformed shape of the mat that is convex (concave downwards). On the contract, the cFEA and the pseudo-coupled approach result in a concave upwards deformed shape. As will be shown later, the differences in the mat settlement profiles have a significant impact on the bending moment diagrams. The average magnitude of mat settlement from the two linear elastic approaches (uncoupled and pseudo-coupled) is comparable to that from elastoplastic cFEA, with the differences being less than 25%, indicating the effectiveness of the choices made for the calculation of \u003cem\u003eE\u003c/em\u003e\u003csup\u003e*\u003c/sup\u003e.\u003c/p\u003e \u003cp\u003eAs the column magnitude increases, the differential settlements from cFEA and the pseudo-coupled approach become smaller (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e), and the mat \u0026ldquo;dishing\u0026rdquo; is less pronounced. This is because the larger applied loading induces larger shear stresses and strains in the soil near the perimeter of the mat, leading to more plastic deformations. Due to the plastic deformations, the overall (secant) stiffness of soil close or at the mat perimeter becomes smaller compared to that near the mat center. The development of plastic straining is expected to be reflected in the shape of the equivalent \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution, rendering it load-magnitude-dependent.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cdiv id=\"Sec12\" class=\"Section2\"\u003e \u003ch2\u003e4.1 Equivalent \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution and bending moment diagrams\u003c/h2\u003e \u003cp\u003eFigure\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003e show examples of back-calculated spring constants \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e from three finite element analyses with different clay properties and loading magnitude. Similar to the case of linear elastic soil (Loukidis and Tamiolakis 2017), the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distributions exhibit a distinctive \u0026ldquo;cup\u0026rdquo; shape, with high values at the corners and along the edges of the mat that decrease sharply towards the foundation\u0026rsquo;s central region, where an almost constant value is attained. In the remainder, comparisons of the parametric analysis results will be made based on cross-sections of the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution along the mat centerline and on along a mat edge, with the values normalized with respect to the reference value \u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cdiv id=\"Sec13\" class=\"Section3\"\u003e \u003ch2\u003e4.1.1 Effect soil yielding\u003c/h2\u003e \u003cp\u003eThe effect of the critical state friction angle \u003cem\u003eφ\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e on the normalized spring constant distribution (\u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e) can be seen in Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003ea which shows the results of three analyses with friction angle varying from 25\u0026deg; to 35\u0026deg;. We can see that \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e is unaffected by the critical friction angle in the inner region of the mat since the soil at the center of the foundation remains mainly elastic. However, \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e along the perimeter of the mat decreases as the friction angle decreases. This is because the lower shear strength of a clay having small \u003cem\u003eφ\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e results to earlier plastic yielding near the edges and more plastic straining. This change in edge spring constants affects the bending moments. More specifically, increasing the friction angle leads to a slight increase in the positive (bottom fiber in tension) moments throughout the central axis of the mat foundation because of the higher spring constants at the perimeter of the mat (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003eb).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe dependence of the edge spring constants on plastic straining gives also rise to dependence on the loading magnitude. So, unlike the findings of Loukidis and Tamiolakis (2017) for purely linear elastic soil behavior, increasing the loading magnitude applied on the foundation and the mat thickness (increased self-weight) leads to a progressive reduction of the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e especially at the edges (Figs.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003e and \u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e8\u003c/span\u003e). Figure\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e8\u003c/span\u003ea indicates that as the load transferred to the ground increases, \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e values at the edges decrease substantially (by 50% or more for every doubling of \u003cem\u003eQ\u003c/em\u003e). This decrease can be attributed to large shear strains near the perimeter of the mat, but also on large contact stresses in the same region, which may cause yielding also due to increased \u003cem\u003ep\u003c/em\u003e\u0026prime; values that approach (or even exceed \u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e) leading to the development of contractive plastic volumetric (consolidation) strains. On the contrary, the decrease in the central region of the mat is rather insignificant for \u003cem\u003eQ\u003c/em\u003e up to 1000kN (total loading, including mat self-weight 25807kN), corresponding to a mobilized factor of safety roughly equal to 5.0 under drained conditions (1.84 under undrained conditions). For larger loading (smaller mobilized factor of safety values), the spring constants decrease notably also in the central region of the mat. This observation suggests that the expansion of the plastic deformations towards the mat central region is starting to affect the central \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e values. However, due to the limited presence of shearing under the central region, this rate of decrease is smaller than that for the edge springs and, as a result, the edge and central \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distributions progressively converge to each other.\u003c/p\u003e \u003cp\u003eRegarding the bending moments, as expected, the more rigid the mat foundation is (large mat thickness) or the higher the applied loading is, the larger the bending moments across the mat foundation become. However, if the moment diagrams are normalized with respect to sum of the applied loading (Σ\u003cem\u003eQ\u003c/em\u003e), the fact that the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e values at the edge and the center converge leads to a decrease in the normalized positive (bottom fiber in tension) moments (Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e8\u003c/span\u003eb). This means that doubling the mat loading does not cause a doubling of the peak moments, as would be expected if the soil were assumed to be a linear elastic material.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe smaller is the coefficient of lateral earth pressure at rest \u003cem\u003eK\u003c/em\u003e\u003csub\u003e0\u003c/sub\u003e (=\u003cem\u003eσ\u003c/em\u003e\u0026prime;\u003csub\u003ehorizontal\u003c/sub\u003e/ \u003cem\u003eσ\u003c/em\u003e\u0026prime;\u003csub\u003evertical\u003c/sub\u003e) from unity (1.0), the closer is the initial soil state of an overconsolidated clay to the yield surface. A \u003cem\u003eK\u003c/em\u003e\u003csub\u003e0\u003c/sub\u003e value equal to 1 means that the effective vertical and horizontal stresses are equal, and thus the initial (in-situ) deviatoric stress is equal to zero. Contrarily, a \u003cem\u003eK\u003c/em\u003e\u003csub\u003e0\u003c/sub\u003e\u0026thinsp;=\u0026thinsp;0.5 means that the effective vertical stresses are twice as large as the effective horizontal stresses, bringing the in-situ soil stress state closer to shear yielding. Hence, the clay under the mat (especially close to the perimeter where shearing is most dominant) is expected to develop plastic shear strains early during loading, resulting to smaller secant stiffness and, thus, comparatively smaller equivalent spring constants. The results from analyses with \u003cem\u003eK\u003c/em\u003e\u003csub\u003e0\u003c/sub\u003e varying between 0.5 and 1.0 (Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e9\u003c/span\u003ea) indicate that, as \u003cem\u003eK\u003c/em\u003e\u003csub\u003e0\u003c/sub\u003e decreases, the ratio of mat edge \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e values to those at the mat center becomes smaller. However, this reduction in \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e contrast is small and, as a consequence, the \u003cem\u003eK\u003c/em\u003e\u003csub\u003e0\u003c/sub\u003e parameter has minimal effect on the bending moment diagram (Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e9\u003c/span\u003eb).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec14\" class=\"Section3\"\u003e \u003ch2\u003e4.1.2 Effect of clay deformability\u003c/h2\u003e \u003cp\u003eThe swelling index \u003cem\u003eC\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e controls the elastic stiffness in the MCC model and is directly related to the Young\u0026rsquo;s modulus \u003cem\u003eE\u003c/em\u003e in Eq.\u0026nbsp;(\u003cspan refid=\"Equ4\" class=\"InternalRef\"\u003e4\u003c/span\u003e). According to the results shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e10\u003c/span\u003ea, the normalized spring constant (\u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e) distributions are quite insensitive to the \u003cem\u003eC\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e value. Its influence is limited to a slight increase of the edge values as the swelling index decreases. Most importantly, the shape of the distributions remains practically the same and, as a consequence, a change in \u003cem\u003eC\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e has minimal effect on the bending moment diagram (Fig.\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e10\u003c/span\u003eb).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eFigure\u0026nbsp;\u003cspan refid=\"Fig11\" class=\"InternalRef\"\u003e11\u003c/span\u003ea reveals that, by decreasing the Poisson\u0026rsquo;s ratio \u003cem\u003eν\u003c/em\u003e, \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e increases at the edges, where the soil is subjected mostly to shearing. This observation can be explained by the shear modulus equation Eq.\u0026nbsp;(\u003cspan refid=\"Equ8\" class=\"InternalRef\"\u003e9\u003c/span\u003e) for the MCC model:\u003cdiv id=\"Equ8\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ8\" name=\"EquationSource\"\u003e\n$$\\:G=\\frac{3\\left(1-2v\\right)}{2\\left(1+v\\right)}{K}_{\\text{b}}=\\frac{3\\left(1-2v\\right)}{2\\left(1+v\\right)}\\frac{(1+e){p}^{{\\prime\\:}}}{\\kappa\\:}$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e9\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eIn the MCC model, the bulk modulus \u003cem\u003eK\u003c/em\u003e\u003csub\u003eb\u003c/sub\u003e is assumed to be independent of \u003cem\u003eν\u003c/em\u003e; it is only affected by the change in mean effective stress, the swelling index and the void ratio. Hence, by decreasing \u003cem\u003eν\u003c/em\u003e, the shear modulus \u003cem\u003eG\u003c/em\u003e becomes larger and the reaction to shear deformation gets enhanced. Nonetheless, this increase in the Poisson\u0026rsquo;s ratio does not practically affect the bending moments (Fig.\u0026nbsp;\u003cspan refid=\"Fig11\" class=\"InternalRef\"\u003e11\u003c/span\u003eb).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe impact of the compression index (\u003cem\u003eC\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e) on the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distribution depends heavily on the magnitude of the applied pressure and to what extent the preconsolidation pressure is exceeded at the foundation level. Figure\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e12\u003c/span\u003ea illustrates that \u003cem\u003eC\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e is not affecting the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distribution for a moderately overconsolidated profile soil (eroded overburden soil \u003cem\u003eD\u003c/em\u003e\u0026thinsp;=\u0026thinsp;10m, corresponding to overconsolidation ratio \u003cem\u003eOCR\u003c/em\u003e\u0026thinsp;=\u0026thinsp;2.5 at a depth \u003cem\u003eB\u003c/em\u003e/3 below the mat) when subjected to loading of small magnitude (\u003cem\u003eQ\u003c/em\u003e\u0026thinsp;=\u0026thinsp;500kN, corresponding to mobilized factor of safety SF\u0026thinsp;=\u0026thinsp;2.31 against bearing capacity). However, decreasing the preconsolidation pressure of the soil (eroded overburden soil \u003cem\u003eD\u003c/em\u003e\u0026thinsp;=\u0026thinsp;5m, corresponding to \u003cem\u003eOCR\u003c/em\u003e\u0026thinsp;=\u0026thinsp;1.8 at a depth \u003cem\u003eB\u003c/em\u003e/3 below the mat), the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distribution at the center and the edges of the mat foundation decreases as \u003cem\u003eC\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e increases (Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e12\u003c/span\u003ec), reflecting a decrease in overall (secant) shear stiffness. Notably, increasing \u003cem\u003eC\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e by a factor of 6 (from 0.1 to 0.6) has minimal effect on the moment diagrams (Figs.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e12\u003c/span\u003eb and \u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e12\u003c/span\u003ed). This is attributed to the fact that despite the equivalent springs becoming substantially softer with increasing \u003cem\u003eC\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e12\u003c/span\u003ec), the edge-to-center ratio of \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e values remains in the 1.6 to 1.8 range.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec15\" class=\"Section3\"\u003e \u003ch2\u003e4.1.3 Effect of overconsolidation\u003c/h2\u003e \u003cp\u003eThe \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distribution depends heavily on the preconsolidation pressure (\u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e) of the soil, provided sufficient loading is applied to produce plastic yielding (Fig.\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e13\u003c/span\u003e). In this comparison, the preconsolidation pressure is represented by \u003cem\u003eD\u003c/em\u003e, which is the thickness of the soil removed by erosion. Four \u003cem\u003eD\u003c/em\u003e values are examined, namely 5, 10, 20 and 40m, which correspond to \u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e at the ground surface equal to 27, 54, 109 and 218 kPa, respectively (with \u003cem\u003eOCR\u003c/em\u003e at depth \u003cem\u003eB\u003c/em\u003e/3 below the mat being equal to 1.8, 2.5, 4.0, 7.1, respectively). Figure\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e13\u003c/span\u003e indicates that an increase in \u003cem\u003ep\u003c/em\u003e\u003csub\u003ec\u003c/sub\u003e generally leads to larger spring constants because of increased overall resistance to settlement resulting from an enhanced shear strength. However, it can be seen that, for highly overconsolidated clay (\u003cem\u003eD\u003c/em\u003e\u0026thinsp;=\u0026thinsp;20m and 40m), the distributions are close to each other, especially in the case of the central distribution. This is because in these cases the clay remains mostly in an elastic state. For low \u003cem\u003eQ\u003c/em\u003e values in particular (Figs.\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e13\u003c/span\u003ea, and \u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e13\u003c/span\u003ec), the applied loading is not large enough to produce plastic yielding at the perimeter of the foundation, causing no significant change on the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distributions (both central and edge), as well as the corresponding bending moment diagrams. As a result, having stiff (high \u003cem\u003eOCR\u003c/em\u003e) clay generally leads to larger bending moments within the mat compared to the case of a softer, slightly overconsolidated clay. As previously discussed, an increase in the loading magnitude causes the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distribution in the \u003cem\u003eD\u003c/em\u003e\u0026thinsp;=\u0026thinsp;10m analysis to decrease at the edges. For a slightly overconsolidated clay (\u003cem\u003eD\u003c/em\u003e\u0026thinsp;=\u0026thinsp;5m), the center and the edge \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distributions decrease even more and converge, even for low loading magnitude (\u003cem\u003eQ\u003c/em\u003e\u0026thinsp;=\u0026thinsp;500kN), due to the early development of plastic strains. As a result, soft clay (low \u003cem\u003eOCR\u003c/em\u003e) tends to produce a flatter \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distribution subjecting the mat foundation to less dishing behavior, which results in lower peak (bottom fiber in tension) bending moments.\u003c/p\u003e \u003cp\u003eIt is interesting to note that, although the peak moments for \u003cem\u003eQ\u003c/em\u003e equal to 500kN for slightly overconsolidated soil (\u003cem\u003eD\u003c/em\u003e\u0026thinsp;=\u0026thinsp;5m) are the lowest and increase as \u003cem\u003eOCR\u003c/em\u003e increases (Fig.\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e13\u003c/span\u003eb), a reversal in the trend is observed for higher loads (e.g., \u003cem\u003eQ\u003c/em\u003e\u0026thinsp;=\u0026thinsp;1000kN) between the cases of \u003cem\u003eD\u003c/em\u003e\u0026thinsp;=\u0026thinsp;5m and \u003cem\u003eD\u003c/em\u003e\u0026thinsp;=\u0026thinsp;10m (Figs.\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e13\u003c/span\u003ed and \u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e13\u003c/span\u003ef). This signifies that there is a second effect coming from a factor other than the enhancement of shear strength caused by an increased degree of overconsolidation. The second factor is the increase of elastic soil stiffness caused by the overconsolidation. The effect of the elastic stiffness on the bending moment distribution can be better visualized for a purely linear elastic analysis with different Young\u0026rsquo;s modulus (Fig.\u0026nbsp;\u003cspan refid=\"Fig14\" class=\"InternalRef\"\u003e14\u003c/span\u003e). The stiffer the soil (higher the Young\u0026rsquo;s modulus) is, the smaller the bending moments tend to become. Therefore, by decreasing the \u003cem\u003eOCR\u003c/em\u003e, the soil exhibits smaller elastic stiffness, and the bending moments tend to increase, counterbalancing the reduction of bending moments caused by the early development of plastic shear strain in the case of slightly overconsolidated clay.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec16\" class=\"Section3\"\u003e \u003ch2\u003e4.1.4 Effect of load eccentricity\u003c/h2\u003e \u003cp\u003eThe effect of eccentricity on \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e can be seen in Fig.\u0026nbsp;\u003cspan refid=\"Fig15\" class=\"InternalRef\"\u003e15\u003c/span\u003ea, which presents the results from analysis with the same mat and soil properties but with different loading type: symmetric (C1), single axis eccentricity (C2) and double axis eccentricity (C3). As observed for linear elastic soil (Loukidis and Tamiolakis 2017), the presence of load eccentricity rotates the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e/\u003cem\u003eK\u003c/em\u003e\u003csub\u003er\u003c/sub\u003e distribution towards the side that the eccentricity points to. This is particularly true for the distribution of the edge spring constants. Furthermore, as eccentricity rotates the Winkler spring distribution, the central bending moment diagram behaves similarly (Fig.\u0026nbsp;\u003cspan refid=\"Fig15\" class=\"InternalRef\"\u003e15\u003c/span\u003eb).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e \u003c/div\u003e \u003cdiv id=\"Sec17\" class=\"Section2\"\u003e \u003ch2\u003e4.2 Comparison with different methods of analysis\u003c/h2\u003e \u003cp\u003eIn this section, we examine how much error is introduced in the calculation of the bending moments diagrams by a) assuming that the soil is purely linear elastic, i.e., ignoring the development of plastic strains and b) using simplified modulus of subgrade reaction distributions. Figure\u0026nbsp;\u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003e compares the bending moments diagrams predicted using the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution back-calculated from the elastoplastic cFEA (Modified Cam-Clay) to those produced i) using \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution back-calculated from linear elastic cFEA in which the soil Young\u0026rsquo;s modulus is assumed to be equal to \u003cem\u003eE\u003c/em\u003e\u003csup\u003e*\u003c/sup\u003e, ii) assuming a uniform modulus of subgrade reaction (\u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e) across the mat, iii) using the 3-zone pseudo-coupled approach (Coduto 2001), and iv) using the Discrete Area Method (DAM) described in the methodology section (Section 3). For the cases presented in Fig.\u0026nbsp;\u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003e, the mobilized safety factor (\u003cem\u003eSF\u003c/em\u003e) against undrained bearing capacity ranges from 1.88 to 5.56. It should be noted that, in practice, shallow foundations are required to have an \u003cem\u003eSF\u003c/em\u003e against undrained bearing capacity at least greater than 2.\u003c/p\u003e \u003cp\u003eAs expected, the analysis with the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution from purely elastic cFEA predicts larger bending moments than that from the elastoplastic cFEA. For the examples shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003e, the elastic cFEA overestimates the peak bending moments at the column locations by 20\u0026ndash;60%.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe reason for this can be illustrated by comparing the respective \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distributions. Figure\u0026nbsp;\u003cspan refid=\"Fig17\" class=\"InternalRef\"\u003e17\u003c/span\u003e compares the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distributions for the cFEA analyses of Figs.\u0026nbsp;\u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003eb and \u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003ed. The assumption of linear elasticity leads to higher edge-to-center \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e ratio irrespective of the applied load magnitude and more variability in the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution because of the absence of plastic yielding at the edges of the mat foundation.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eRegarding existing approaches used in practice, the use of a uniform \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution largely underestimates the peak positive (bottom fiber in tension) bending moments and overestimates the peak negative (top fiber in tension) bending moments. This is because the foundation has less tendency to bend with a concave upwards shape (dishing) due to the complete absence of higher spring constants at the mat perimeter than at its central region. Conversely, modeling the soil through the simple pseudo-coupled approach where \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e decreases from the mat edges towards the center (Coduto 2001) produces a realistic dishing profile (concave upwards) and larger peak bending moments. For mobilized safety factor \u003cem\u003eSF\u003c/em\u003e against undrained bearing capacity close to 2 (Figs.\u0026nbsp;\u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003ea, and \u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003eb), the pseudo-coupled distribution predictions turn out to be very close to those from elastoplastic cFEA. This is because plastic deformations reduce the requirement for large \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e at the mat edge. However, for higher mobilized \u003cem\u003eSF\u003c/em\u003e values (greater than 3), which are much more frequently encountered in practice, the simple pseudo-coupled distribution does not produce accurate bending moments (Figs.\u0026nbsp;\u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003ec, and \u003cspan refid=\"Fig16\" class=\"InternalRef\"\u003e16\u003c/span\u003ed), since it assumes that \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e at the edges is only twice that at the center, while in fact it has to be much higher. ACI 318\u0026thinsp;\u0026minus;\u0026thinsp;19(\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e) states that doubling the \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e at the edges introduces a certain degree of spring coupling effects but is inadequate. On the contrary, the Discrete Area Method employed herein, which accounts for the inelastic consolidation settlements through Settle3 calculations, exhibits a superior performance. In all cases, employing the Discrete Area Method is a better alternative to be used in practice because it predicts similar moments to those produced by the elastoplastic finite element analysis.\u003c/p\u003e \u003cp\u003eFinally, Fig.\u0026nbsp;\u003cspan refid=\"Fig18\" class=\"InternalRef\"\u003e18\u003c/span\u003e plots the shear force (\u003cem\u003eV\u003c/em\u003e\u003csub\u003ex\u003c/sub\u003e) diagrams along the mat centerline predicted by the different methods of analysis. Similar to the findings of Loukidis and Tamiolakis (2017) for the linear elastic case, the choice of \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e (or \u003cem\u003ek\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e) distribution has practically very little effect on the shear force diagram. This is because static equilibrium dictates that the peak values of the diagram need to be exactly equal to the applied concentrated (column) loads, leaving little possibility for variation across a span among the different methods. Most importantly, the small discrepancies observed midway between two column loads in Fig.\u0026nbsp;\u003cspan refid=\"Fig18\" class=\"InternalRef\"\u003e18\u003c/span\u003e are inconsequential with respect to the shear resistance requirements (i.e., punching failure checks).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e"},{"header":"5. CONCLUSIONS","content":"\u003cp\u003eThe appropriate Winkler spring constant (\u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e) distribution for mat foundation analysis is best captured using three-dimensional continuum finite element analysis (cFEA). In this study, the spring constants suitable for analysis of mats on clay for long-term (drained) conditions were back-calculated on the basis that the soil exhibits elastoplastic behavior following the Modified Cam-Clay model. The effects of various soil and mat parameters on the resulting \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distributions as well as on the bending moments were investigated. The cFEA-based predictions were compared to those from various existing modulus of subgrade reaction variation approaches encountered in foundation engineering practice. Based on the results of this study, the following main conclusions can be derived:\u003c/p\u003e \u003cp\u003e(\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e) The development of plastic deformations affects the shape of the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution by lowering the equivalent spring constant values at the edges of the mat foundation where most of the shearing occurs and the contact stresses are larger. Closer to the center of the foundation, the soil remains mostly elastic, unless superstructure loading is large enough to cause exceedance of the preconsolidation pressure in the clay close to the central area of the mat.\u003c/p\u003e \u003cp\u003e(\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e) The shape of the \u003cem\u003eK\u003c/em\u003e\u003csub\u003es\u003c/sub\u003e distribution is insensitive to the clay deformability parameters, namely the compression index, the swelling index and the Poisson\u0026rsquo;s ratio. It is strongly affected by the degree of overconsolidation of the clay deposit and the magnitude of mat loading.\u003c/p\u003e \u003cp\u003e(\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e) Ignoring the elastoplastic behavior of the clay leads to an overprediction of the spring constants at the mat edges, which, in turn, results in an overly conservative estimate of the peak bending moments.\u003c/p\u003e \u003cp\u003e(\u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e) A variant of the Discrete Area Method in which consolidation settlements are taken into account predicts bending moment distributions that are very similar to those based on elastoplastic cFEA. Hence, DAM comes out as a superior alternative for use in foundation engineering practice.\u003c/p\u003e"},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eFunding Declaration:\u0026nbsp;\u003c/strong\u003eNo Funding.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCompeting Interests:\u003c/strong\u003e The authors have no competing interests to declare that are relevant to the content of this article.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAuthors Contribution:\u0026nbsp;\u003c/strong\u003eThe authors have contributed equally to this work.\u003cstrong\u003e\u0026nbsp;\u003c/strong\u003eThe first draft of the manuscript was written by Gaby Saad and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eFinancial/non-financial interest:\u0026nbsp;\u003c/strong\u003eThe authors have no relevant financial or non-financial interests to disclose.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eEthical approval:\u0026nbsp;\u003c/strong\u003eThe authors declare that this manuscript have not been submitted to any other journal for simultaneous consideration. The work is original and not published elsewhere.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eData Availability:\u0026nbsp;\u003c/strong\u003eSome or all data, models, or code that support the findings of this study are available from the corresponding author upon reasonable request.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\n\u003cli\u003eAbdullah, W. S. (2008). New elastoplastic method for calculating the contact pressure distribution under rigid foundations. \u003cem\u003eJordan journal of civil engineering\u003c/em\u003e,\u003cem\u003e 2\u003c/em\u003e(1), 71-89.\u003c/li\u003e\n\u003cli\u003eACI (2019). 318-19; Building Code Requirements for Structural Concrete and Commentary. \u003cem\u003eAmerican Concrete Institute: Farmington Hills, MI, USA\u003c/em\u003e.\u003c/li\u003e\n\u003cli\u003eAlzoaby, H., Saad, G. \u0026amp; Abou-Jaoude, G. Implementation of the discrete area method and its impact on the steel reinforcement of large mat foundations. Innov. Infrastruct. Solut. 10, 131 (2025). https://doi.org/10.1007/s41062-025-01926-x \u003c/li\u003e\n\u003cli\u003e\u003cem\u003eAnsys\u0026reg; Academic Research Mechanical, Release 2023 R2\u003c/em\u003e\u003c/li\u003e\n\u003cli\u003eBanavalkar, P.V. (1995). Mat Foundation and Its interaction with the Superstructure. ACI Special Publication, 152.\u003c/li\u003e\n\u003cli\u003eCEN (2004). Eurocode 7 Geotechnical design. Part 1: General rules. EN 1997-1:2004, Brussels: European Committee for Standardization.\u003c/li\u003e\n\u003cli\u003eChow, Y. K. (1987). Vertical deformation of rigid foundations of arbitrary shape on layered soil media. \u003cem\u003eInternational journal for numerical and analytical methods in geomechanics\u003c/em\u003e, \u003cem\u003e11\u003c/em\u003e(1), 1-15.\u003cspan dir=\"RTL\"\u003e\u0026rlm;\u003c/span\u003e\u003c/li\u003e\n\u003cli\u003eCoduto, D. P. (2001). \u003cem\u003eFoundation design: principles and practices\u003c/em\u003e. Prentice Hall.\u003c/li\u003e\n\u003cli\u003eDutta, S. C., \u0026amp; Roy, R. (2002). A critical review on idealization and modeling for interaction among soil\u0026ndash;foundation\u0026ndash;structure system. \u003cem\u003eComputers \u0026amp; structures\u003c/em\u003e,\u003cem\u003e 80\u003c/em\u003e(20-21), 1579-1594.\u003c/li\u003e\n\u003cli\u003eEdgers, L., Alonge, J. L., \u0026amp; Sanayei, M. (2005). Modeling the effects of soil-structure interaction on a tall building bearing on a mat foundation. \u003cem\u003eCivil Eng. Pract\u003c/em\u003e,\u003cem\u003e 20\u003c/em\u003e, 51-68.\u003c/li\u003e\n\u003cli\u003eEstephan, R., Skeini, H., Malek, C., Mounir, M., Rahim, M.A., Abou-Jaoude, G., Yacoub, T. and Valiulin, D. (2022). A new design tool for shallow foundations offering enhanced accuracy, reduced design time, and cost savings. Proceedings of the 20\u003csup\u003eth\u003c/sup\u003e International Conference on Soil Mechanics and Geotechnical Engineering, Sydney, Australia, 5169-5173.\u003c/li\u003e\n\u003cli\u003eHorvath, J. (1995). Subgrade modeling for mat foundations: a review and critique of analytical methods. \u003cem\u003eSpecial Publication\u003c/em\u003e,\u003cem\u003e 152\u003c/em\u003e, 117-160.\u003c/li\u003e\n\u003cli\u003eHorvath, J. S., \u0026amp; Colasanti, R. J. (2011). Practical subgrade model for improved soil-structure interaction analysis: Model development. \u003cem\u003eInternational Journal of Geomechanics\u003c/em\u003e,\u003cem\u003e 11\u003c/em\u003e(1), 59-64.\u003c/li\u003e\n\u003cli\u003eJaky, J. (1944). The coefficient of earth pressure at rest. Journal of the Society of Hungarian Architects and engineers, 355-358.\u003c/li\u003e\n\u003cli\u003eLarkela, A., Mengelt, M., \u0026amp; Stapelfeldt, T. (2013). Determination of distribution of modulus of subgrade reaction. Proc., 18th Int. Conf. on Soil Mechanics and Geotechnical Engineering, French Society for Soil Mechanics and Geotechnical Engineering (CFMS), Paris, France,\u003c/li\u003e\n\u003cli\u003eLeonidou, A. K. (2021). \u003cem\u003eSpatial distribution of winkler spring stiffness for the structural analysis of embedded rectangular mat foundations\u003c/em\u003e. MSc Thesis, University of Cyprus. Cyprus.\u003c/li\u003e\n\u003cli\u003eLoukidis, D., \u0026amp; Tamiolakis, G.-P. (2017). Spatial distribution of Winkler spring stiffness for rectangular mat foundation analysis. \u003cem\u003eEngineering Structures\u003c/em\u003e,\u003cem\u003e 153\u003c/em\u003e, 443-459.\u003c/li\u003e\n\u003cli\u003e\u003cem\u003e. MathWorks (2022). MATLAB version: 9.13.0 (R2022b), Natick, Massachusetts: The MathWorks Inc. https://www.mathworks.com\u003c/em\u003e.\u003cem\u003e \u003c/em\u003e\u003c/li\u003e\n\u003cli\u003eMayne, P. W., and Kulhawy, F. H. (1982). K\u003csub\u003e0\u003c/sub\u003e\u0026ndash; OCR relationship in soil.\u0026rdquo; J. Geotech. Eng. Div., Am. Soc. Civ. Eng., 108(6), 851\u0026ndash;872.\u003c/li\u003e\n\u003cli\u003eMeyerhof, G. G. (1953). Some recent foundation research and its application to design. The Structural Engineer, 31(6), 151-167.\u003c/li\u003e\n\u003cli\u003eNoorzaei, J., Viladkar, M., \u0026amp; Godbole, P. (1995). Elasto-plastic analysis for soil-structure interaction in framed structures. \u003cem\u003eComputers \u0026amp; structures\u003c/em\u003e,\u003cem\u003e 55\u003c/em\u003e(5), 797-807.\u003c/li\u003e\n\u003cli\u003ePishilis, R. (2022). \u003cem\u003eThe influence of the elastoplastic behavior of the soil on the spatial distribution of the stiffness of the Winkler springs for mat foundations on clay.\u003c/em\u003e MSc Thesis, University of Cyprus. Cyprus. (in greek)\u003c/li\u003e\n\u003cli\u003ePotts, D. M., Zdravković, L. (1999). \u003cem\u003eFinite element analysis in geotechnical engineering: Theory\u003c/em\u003e (Vol. 1). Thomas Telford London.\u003c/li\u003e\n\u003cli\u003ePoulos, H. G. (2016). Tall building foundations: design methods and applications. \u003cem\u003eInnovative Infrastructure Solutions\u003c/em\u003e, \u003cem\u003e1\u003c/em\u003e, 1-51.\u003cspan dir=\"RTL\"\u003e\u0026rlm;\u003c/span\u003e\u003c/li\u003e\n\u003cli\u003eRoscoe, K. H., \u0026amp; Burland, J. B. (1968). On the generalized stress\u0026ndash;strain behaviour of the \u0026lsquo;wet\u0026rsquo;clay. Pp. 535\u0026ndash;609 in:(J. Heyman and FA Leckie, editors). In Cambridge University Press. Cambridge.\u003c/li\u003e\n\u003cli\u003eRoscoe, K. H., Schofield, A., \u0026amp; Wroth, A. P. (1958). On the yielding of soils. \u003cem\u003eGeotechnique\u003c/em\u003e, \u003cem\u003e8\u003c/em\u003e(1), 22-53.\u003c/li\u003e\n\u003cli\u003eSadek, M., Mroueh, H., \u0026amp; Shahrour, I. (2010). Influence of nonlinearity on the stress distribution in the soil\u0026mdash;application to road engineering problems. \u003cem\u003eJournal of transportation engineering\u003c/em\u003e,\u003cem\u003e 136\u003c/em\u003e(1), 77-83.\u003c/li\u003e\n\u003cli\u003eSadek, M., \u0026amp; Shahrour, I. (2007). Use of the Boussinesq solution in geotechnical and road engineering: influence of plasticity. \u003cem\u003eComptes rendus. M\u0026eacute;canique\u003c/em\u003e,\u003cem\u003e 335\u003c/em\u003e(9-10), 516-520.\u003c/li\u003e\n\u003cli\u003eSaeed, H. H. (2022). Effect of soil nonlinearity on analysis of raft foundation. \u003cem\u003eCivil and Environmental Engineering\u003c/em\u003e,\u003cem\u003e 18\u003c/em\u003e(1), 292-300.[30] Salgado, R. (2008). \u003cem\u003eThe engineering of foundations\u003c/em\u003e. McGraw-Hill Europe.\u003c/li\u003e\n\u003cli\u003eSchultze, E. (1961). Distribution of stress beneath a rigid foundation. In \u003cem\u003eProc., 5th Int. Conf. on Soil Mechanics and Foundation Engineering\u003c/em\u003e (pp. 807-813), Paris.\u003c/li\u003e\n\u003cli\u003eSteinbrenner, W. (1934). Tafeln zur setzungsberechnung. \u003cem\u003eDie Stra\u0026szlig;e\u003c/em\u003e,\u003cem\u003e 1\u003c/em\u003e.\u003c/li\u003e\n\u003cli\u003eTamiolakis G-P. (2012). \u003cem\u003eStudy of the spatial distribution of the Winkler spring stiffness for the static analysis of mat foundations.\u003c/em\u003e MSc Thesis, University of Cyprus, Cyprus. (in greek).\u003c/li\u003e\n\u003cli\u003eTeng, W.C. (1969). \u003cem\u003eFoundation Design\u003c/em\u003e. Prentice-Hall.\u003c/li\u003e\n\u003cli\u003eUlrich, E. (1995). Introduction to the State-of-the-Art Mat Foundation Design and Construction. \u003cem\u003eSpecial Publication\u003c/em\u003e,\u003cem\u003e 152\u003c/em\u003e, 1-12.\u003c/li\u003e\n\u003cli\u003eUlrich, E. J. (1991). Subgrade reaction in mat foundation design. \u003cem\u003eConcrete International\u003c/em\u003e, \u003cem\u003e13\u003c/em\u003e(4), 41-50.\u003cspan dir=\"RTL\"\u003e\u0026rlm;\u003c/span\u003e\u003c/li\u003e\n\u003cli\u003eViladkar, M., Godbole, P., \u0026amp; Noorzaei, J. (1991). Soil-structure interaction in plane frames using coupled finite-infinite elements. \u003cem\u003eComputers \u0026amp; structures\u003c/em\u003e,\u003cem\u003e 39\u003c/em\u003e(5), 535-546.\u003c/li\u003e\n\u003cli\u003eViladkar, M., Noorzaei, J., \u0026amp; Godbole, P. (1994). Interactive analysis of a space frame-raft-soil system considering soil nonlinearity. \u003cem\u003eComputers \u0026amp; structures\u003c/em\u003e,\u003cem\u003e 51\u003c/em\u003e(4), 343-356.\u003c/li\u003e\n\u003cli\u003eWinkler, E. (1867). Die Lehre von Elastizitat und Festigkeit (The theory of elasticity and stiffness). \u003cem\u003eH. Domenicus. Prague\u003c/em\u003e.\u003c/li\u003e\n\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":true,"highlight":"","institution":"","isAcceptedByJournal":false,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true},"keywords":"Soil-structure Interaction, Mat foundation, Finite element analysis, Winkler Springs, Modified Cam-Clay model","lastPublishedDoi":"10.21203/rs.3.rs-6456559/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-6456559/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003eThe interaction between mat foundations and the supporting soil is often simplified in design practice by assuming that the mat rests on Winkler springs. Most of previous studies aimed at determining the spatial variation of the spring constants or the modulus of subgrade reaction that is needed to obtain realistic and accurate results regarding the bending of the mat assuming that the soil behaves as a purely linear elastic material. This paper investigates the effects of plastic strains and volume changes that happen due to the consolidation of clayey soils under long-term conditions on the distribution of the equivalent spring constants and the bending moment diagrams. For this purpose, three-dimensional finite element analyses are performed for fully drained conditions with the soil simulated using the Modified Cam-Clay constitutive model. The results show that the presence of plastic yielding significantly influences the mat-soil interaction and the assumption of a linear elastic soil overpredicts the soil stiffness and the bending moments in the mat foundation. The semi-analytical Discrete Area Method could be adopted to obtain comparable results to those generated by elastoplastic finite element analysis.\u003c/p\u003e","manuscriptTitle":"Winkler Spring Stiffness Distribution for the Structural Analysis of Mat Foundations in Clay Under Drained Conditions","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2025-05-13 12:44:29","doi":"10.21203/rs.3.rs-6456559/v1","editorialEvents":[{"type":"communityComments","content":0}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"e984000e-a9a3-4ccc-97d7-d2b158d3f8ea","owner":[],"postedDate":"May 13th, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"posted","subjectAreas":[],"tags":[],"updatedAt":"2025-11-02T03:38:27+00:00","versionOfRecord":[],"versionCreatedAt":"2025-05-13 12:44:29","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-6456559","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-6456559","identity":"rs-6456559","version":["v1"]},"buildId":"8U1c8b4HqxoKbykW_rLl7","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}