Evolutional law and elastoplastic constitutive model of structural loess considering moisture contents

TL;DR

This study develops a new model to understand how moisture affects the structural behavior and deformation of loess soil, improving engineering assessments in loess regions.

Contribution

A novel elastoplastic constitutive model for structured loess incorporating structural parameters and moisture effects is proposed.

Findings

The model accurately reflects the impact of initial moisture content on loess structural strength and deformation.

The structural attenuation coefficient effectively captures structural degradation with axial strain.

Numerical simulations in ABAQUS confirm the model's accuracy when compared to triaxial tests.

Abstract

The deformation properties of natural undisturbed loess are intricately linked to its structural qualities. Structural characteristics can markedly improve soil stability, allowing structured loess to sustain a greater void ratio than disturbed loess under equivalent stress circumstances. This study examines the effect of initial moisture content on structural behavior by comparing the isotropic compression curves of disturbed and undisturbed loess under varying moisture conditions, thereby identifying the structural parameters of undisturbed loess and their changes during compression. A structural attenuation coefficient is developed to enhance the stress ratio-based structural parameter, and a quantitative approach for depicting structural degradation with axial strain is proposed, facilitating an effective evaluation of structural damage. Additionally, within the modified Cam–clay model framework, the yield surface equation is adjusted by integrating structural parameters, resulting in an elastoplastic constitutive model for structured loess that uniformly characterizes structural evolution under both spherical and shear stress. Finally, the suggested structural constitutive model is implemented in ABAQUS for numerical analysis, and comparisons with triaxial tests on undisturbed loess and published data confirm the model’s accuracy and applicability. The findings indicate that the established model accurately reflects the impact of initial moisture content on the structural strength, deformation properties, and degradation behavior of loess, offering a dependable theoretical foundation for the quantitative evaluation of loess structure and for engineering assessments in loess regions.

Figures

Production process of samples with different moisture contents.

One-dimensional compression curves of loess: (a) undisturbed loess; (b) disturbed loess.

Sketch map of compression curves for structural soil.

Deviatoric stress-strain curve of undisturbed loess: (a) σ3 = 50kPa; (b) σ3 = 100kPa; (c) σ3 = 150kPa; (d) σ3 = 200kPa.

Deviatoric stress-strain curve of disturbed loess: (a) σ3 = 50kPa; (b) σ3 = 100kPa; (c) σ3 = 150kPa; (d) σ3 = 200kPa.

Damage characteristics of loess samples: (a) Softening damage; (b) Stabilizing damage; (c) Hardening damage.

Evolution law of structural parameters for stress ratio: (a) σ3 = 50kPa; (b) σ3 = 100kPa; (c) σ3 = 150kPa; (d) σ3 = 200kPa.

Model fitting curves of structural parameter variation patterns: (a) σ3 = 50kPa; (b) σ3 = 100kPa; (c) σ3 = 150kPa; (d) σ3 = 200kPa.

Yield surface of structural soil.

Relationship between λμ and 𝐦μ.

Relationship between κ and 𝐦μ.

Relationship between 𝐩1 and 𝐦μ.

Critical state lines of loess with different moisture contents: (a) w = 14%; (b) w = 17%; (c) w = 20%; (d) w = 23%; (e) w = 26%; (f) w = 35.2%.

Flowchart of an improved Euler integration algorithm with error control.

Comparison of predicted and test values of the loess constitutive model: (a) w = 14%; (b) w = 20%; (c) w = 26%.

Comparison of calculated and experimental values for literature specimens: (a) Case 1; (b) Case 2; (c) Case 3.

Equations

Full text

1. Introduction

The unique climatic and historical processes in Northwest China have resulted in the formation of a wide variety of undisturbed structural loess. These loess deposits have a high void ratio and are sensitive to water. The mechanical and deformation characteristics of undisturbed loess are primarily influenced by its structural properties, which sets it apart from disturbed loess. The modified Cam-clay model, which is based on disturbed loess [1], provides a more accurate description of the deformation characteristics of this type of loess. Nevertheless, the impact of structural features on the deformation characteristics of undisturbed loess has been overlooked [2–5], resulting in significant discrepancies from actual conditions when it comes to directing engineering design in loess regions [6]. Hence, it is highly significant from a theoretical and technical perspective to develop a more comprehensive constitutive relationship for undisturbed structural loess, considering its structural properties.

In recent years, significant progress has been made in the study of constitutive relationships of structural soils. Several researchers have developed constitutive models of structural soils within the framework of critical state theory. For example, Liu et al. [7], Mendoza et al. [8], and Chowdhury et al. [9] successfully captured experimental phenomena such as the pronounced increase in compression index after yielding [10–12] and strain-softening behavior observed in triaxial shear tests [13–15]. Zhujiang Shen et al. [16], Mingjing Jiang et al. [17], and Suebsuk et al. [18] established constitutive models of structural soils by introducing a damage function based on damage mechanics theory. Dingyi Xie et al. [19] proposed the concept of comprehensive structural potential to describe the deformation and strength properties of structural soils. Rouainia et al. [20] incorporated structural parameters into the constitutive model of disturbed soils and further defined their evolution law to establish a structural soil model. Meanwhile, Enyang Zhu et al. [21,22], Park et al. [23], and Wu et al. [24] integrated the volumetric strain induced by structural degradation into a unified hardening model, thereby extending the theoretical framework for structural soils. Saberi et al. [25], Desai et al. [26,27], Liu et al. [28], and Zhang et al. [29] formulated constitutive models of structural soils based on the disturbed state concept. Collectively, these studies have provided strong theoretical support for understanding the mechanical properties and structural evolution mechanisms of structural soils. Nonetheless, the current methodologies for assessing structural characteristics continue to have considerable limits. Strain-based parameters are insufficient for characterizing the structural degradation processes caused by shear deformation. While stress-based structural parameters can encompass the cumulative effects of spherical stress and shear stress, they frequently neglect to consider the impact of pre-shear consolidation stress and the function of spherical stress during shearing, nor can they accurately depict the specimen’s initial structural state. Furthermore, many structural characteristics cannot be directly incorporated into constitutive frameworks, creating a distinct disparity between the quantitative analysis of soil structure and the constitutive models used for engineering purposes. Consequently, it is essential to formulate a structural quantification approach and a constitutive model that can integrate the complete process of initial structure, stress/strain response, and structural damage, while facilitating direct application in numerical simulations.

This study examines the disparities in isotropic compression and shear deformation between disturbed and undisturbed loess at varying initial moisture contents, undertaking four tasks to formulate an elastoplastic constitutive model that accurately represents the complete structural evolution process. (1) Coupling mechanism between structural properties and compression behavior: Isotropic compression tests are conducted to determine the initial structural yield stress and the slope of the linear segment of the e–lnp compression curve, thereby revealing the quantitative relationships among the yield surface, soil structure, and compression characteristics; (2) Evolution of structural properties and strength: Based on the stress-ratio structural parameter, a structural attenuation coefficient is introduced to establish a quantitative method for describing the progressive structural deterioration with axial strain, thus obtaining the evolution law of structural strength during structural attenuation and overcoming the limitations of existing parameters in representing shear-induced structural damage; (3) Development of an elastoplastic constitutive model for loess: Within the framework of the modified Cam–clay model, structural parameters that characterize the initial structural state and its degradation are incorporated to modify the yield surface equation, leading to the construction of a structural elastoplastic model that can uniformly describe structural evolution under combined spherical and shear stresses; (4) Numerical implementation and validation: The proposed model is implemented in ABAQUS through the UMAT interface using Fortran [30,31]. Its correctness and applicability are verified by comparing numerical predictions with triaxial test results on undisturbed loess and relevant data reported in the literature. This study aims to establish the inherent connection between the mechanical behavior of loess and its structural evolution, thereby offering an effective theoretical framework and scientific basis for the quantitative characterization of loess structure, investigations into its mechanical properties and constitutive relationships, as well as engineering applications and disaster prevention in loess regions.

2. Deformation characteristics of structural loess

2.1. Compression test of loess

Undisturbed loess has large voids in its structure and a strong cementitious mass connecting them, giving it a high structural strength in undisturbed conditions [20,32]. Once the loess is disturbed, the cement characteristics that hold the huge void together weaken, causing the structure’s strength to be damaged. As a result, the large voids transform into smaller voids that are evenly distributed. Therefore, there exists a fundamental distinction in the compressive characteristics of undisturbed loess and disturbed loess. The physical property parameters of undisturbed loess in the Luochuan Tunnel were determined using the cutting ring method, oven drying method, and liquid-plastic limit combined method. The moisture content of the loess samples ranged from 14.4 to 27.4%, the natural density from 1.63 to 2.1 g/cm3, the void ratio from 0.53 to 1.25, the plasticity index from 8.6 to 14.2, and the liquidity index from −0.27 to 0.84.

A research study was conducted on the mechanical characteristics of undisturbed loess and disturbed loess using 1-D consolidation tests and triaxial tests. The compression process can be considered as a form of structural damage, with a progressive shift from undisturbed loess to disturbed loess. The soil utilized in the experiment is the undisturbed loess obtained from the Luochuan Tunnel, with a soil extraction depth ranging from 34.6 ~ 35m. Immediately after sampling, the specimens were wrapped with multiple layers of plastic film, followed by an additional sealing layer using air bubble film. This double-sealing approach not only enhanced moisture retention but also improved impact resistance, thereby preventing structural alterations caused by water loss. The sealed samples were then placed into rigid cardboard boxes, with inflatable cushioning bags filling the gaps between the samples and the box walls. This buffering arrangement effectively absorbed external vibrations and shocks during transportation, minimizing potential disturbances to the soil specimens.

Experiments were carried out to create samples of undisturbed loess and disturbed loess with moisture contents of 14.0%, 17.0%, 20.0% (plastic limit), 23.0%, 26.0%, 35.2% (saturated moisture content), respectively. Undisturbed soil samples were obtained by the cutting ring method, and samples with different moisture content were obtained by the air-drying method and titration water infusion method. The disturbed soil samples were prepared by the undisturbed soil crushing and compression sampling method, and the samples with different moisture content could be prepared by pre-calculated blending. The testing instruments employed were the GJ-16 medium-pressure consolidation apparatus and the TSY10−2 strain triaxial instrument. The production process for test samples with varying moisture content is illustrated in Fig 1, while the results of one-dimensional consolidation tests for both undisturbed and disturbed soil samples are presented in Fig 2.

Production process of samples with different moisture contents.

One-dimensional compression curves of loess: (a) undisturbed loess; (b) disturbed loess.

Fig 2 (a) shows the one-dimensional compression curve of undisturbed loess, while Fig 2 (b) shows the one-dimensional compression curve of disturbed loess. Fig 2 (b) illustrates that the compression curves of undisturbed loess samples with varying moisture contents resemble nearly linear slopes, indicating a minimal impact of moisture content on compression characteristics. The compression curves of undisturbed loess in Fig 2(a) are significantly affected by structural properties, exhibiting non-linear behavior with a pronounced dependence on moisture content. At lower moisture levels, higher structural strength results in compression curves that can be approximated by two linear segments. With increasing moisture content, structural strength diminishes, and the curve progressively approaches that of disturbed loess. The deformation properties of structural soil are examined based on its development process, with the compression curve of structural soil illustrated schematically in Fig 3. At the initial state (point A), structural strength is minimal. The e-lgp curve is approximated as a linear segment with a minimal slope (A → B). The structural soil specimen contains soil particles that serve as a cementing material. At this point, the pressure is less than the structural strength of undisturbed loess, resulting in compressive deformation that is recoverable as elastic deformation. Upon reaching a specific pressure threshold, the absolute value of the slope of the compression curve increases, at which point the sample’s compression curve is represented by the red solid line SCL (B → C) in the figure. Consequently, the pressure surpasses the soil’s strength, resulting in compression deformation comprising both elastic and plastic components. The compression curve is typically perceived as comprising two linear segments due to its structural characteristics [33]. However, the structural strength of soil varies with moisture content; at lower moisture levels, the soil exhibits greater structural strength, resulting in an increased compressive stress necessary to attain the structural damage threshold. The compression curve of the sample in its remolded form can be shown by the NCL in the figure.

Sketch map of compression curves for structural soil.

2.2. Stress ratio structural parameters

By comparing the triaxial test results of undisturbed soil, disturbed soil, and saturated soil. This paper introduces a structural parameter $[eqn]$ , derived from comprehensive structural potential theory [19,34], which encapsulates both spherical stress and shear stress effects. It analyzes the correlation between $[eqn]$ , moisture content, and confining pressure, specifically expressed in terms of stress ratio:

$$[eqn]$$

Where $[eqn]$ is the stress ratio structural parameters; $[eqn]$ is the structural stability, and $[eqn]$ is the structural variability; where $[eqn]$ is the spherical stress and $[eqn]$ is the generalized shear stress, i.e.,:

$$[eqn]$$ $$[eqn]$$

Where $[eqn]$ , $[eqn]$ and $[eqn]$ are the stresses, respectively. Which $[eqn]$ represents the ratio of the spherical stress to the generalized shear stress, and $[eqn]$ , $[eqn]$ , and $[eqn]$ stand for the shear stress ratios for undisturbed soils, disturbed soils, and saturated soils. This structural parameter, which is the inheritance and development of stress structural parameter, can be seen from the defining equation to reflect the effect of spherical stress on soil structural properties in addition to the effect of shear stress.

The structural characteristics constitute a complicated system of stresses resulting from strong cementation and an unstable arrangement within the soil. Nonetheless, its macroscopic manifestation is that the soil exhibits great strength and limited deformation capacity. According to Equation (1), an increase in $[eqn]$ correlates with enhanced structural stabilizability and a heightened reduction in soil strength resulting from disturbance. The smaller $[eqn]$ is, the higher the structural variability and the more significant the loss of strength resulting from the degradation of soil structure due to applied stresses during immersion. Equation (1) indicates that the stress ratio structural parameter can simultaneously represent the influence of both spherical stress and shear stress on the alteration of soil structural properties. The prior stress structural parameters don’t accurately reflect the impact of spherical stress on soil structure, while the strain structural parameters inadequately depict the soil’s shear damage pattern and strength degradation, among other deficiencies. Thoroughly analyze the alteration of soil structural properties in response to external loads, ensuring enhanced rationality.

2.3. Structural loess evolutional law

2.3.1. Shear characteristics of loess.

Figs 4 and 5 present the deviatoric stress-strain curves for undisturbed and disturbed loess at varying moisture content and confining pressures, illustrating the alterations in soil strength throughout triaxial tests.

Deviatoric stress-strain curve of undisturbed loess: (a) σ3 = 50kPa; (b) σ3 = 100kPa; (c) σ3 = 150kPa; (d) σ3 = 200kPa.

Deviatoric stress-strain curve of disturbed loess: (a) σ3 = 50kPa; (b) σ3 = 100kPa; (c) σ3 = 150kPa; (d) σ3 = 200kPa.

The figure illustrates that the deviatoric stresses of the undisturbed samples exceed those of the disturbed samples under identical conditions, signifying a substantial reduction in the shear strength of the soil samples following disturbance. The moisture content and confining pressure significantly influence the shear strength and shape of the samples’ curves. The increase of moisture content diminishes the connection strength of the soil; as moisture content rises under identical confining pressure, shear strength values of the samples decrease. Conversely, for a constant moisture content, increased confining pressure results in higher shear strength values of the samples. Soil deformation damage can be divided based on stress-strain characteristics, and Fig 6 presents images illustrating the damage characteristics of typical specimens under varying moisture content and perimeter pressure circumstances. At lower moisture content, the structural strength of the sample is enhanced; conversely, under minimal confining pressure, the friction arising from the displacement of sample particles diminishes. The deviatoric stress-strain curve exhibited characteristics of the softening damage type. (e.g., w = 14.0%, w = 17.0%; σ_3_ = 50 kPa, σ3** = 100 kPa). As confining pressure rises, the friction among the sample particles increases, compensating for the reduction in cohesion caused by structural damage, resulting in the deviatoric stress-strain curve of the specimen transitioning into a stabilizing damage type (e.g., w = 17.0%; σ_3_ = 150 kPa, σ_3_ = **200kPa). Soils with a stabilizing damage type undergo plastic deformation earlier than those with a softening damage type. A relatively large deformation is required to reach the peak strength, after which the stress–strain curve remains essentially stable, with strain continuing to accumulate. The soil mass maintains good integrity after failure. This damage type is generally observed in soils under intermediate confining pressure and low moisture content. The higher moisture content and disturbing impacts compromised the structural relationships of the soils, resulting in deviatoric stress-strain curves for undisturbed soils with moisture content over the plastic limit, as well as all disturbed soils, that were predominantly of a hardening damage type.

Damage characteristics of loess samples: (a) Softening damage; (b) Stabilizing damage; (c) Hardening damage.

2.3.2. Stress ratio structural parameter evolution law.

By analyzing the triaxial test results of undisturbed soil, disturbed soil, and saturated soil, combined with Equation. (1), the evolution of the soil stress ratio structural parameter may be derived during the triaxial shear test under varying conditions, as illustrated in Fig 7.

Evolution law of structural parameters for stress ratio: (a) σ3 = 50kPa; (b) σ3 = 100kPa; (c) σ3 = 150kPa; (d) σ3 = 200kPa.

The figure illustrates that as shear displacement increases, the soil’s structural properties progressively decrease until total loss occurs, evidenced by the structural parameter approaching 1, culminating in the complete loss of soil structural properties [34]. The evolution of the stress ratio structural parameter in relation to shear deformation is more consistent; specifically, $[eqn]$ diminishes as $[eqn]$ increases, and the two exhibit a logarithmic correlation with a strong fit, articulated as follows:

$$[eqn]$$

Where A, B are parameters related to moisture content and confining pressure. To deepen understanding of the stress ratio structural parameters under complex stress path test conditions, while simultaneously facilitating numerical theoretical analysis and application. Table 1 presents the fitting results of structural parameters A and B under various test situations, thereby establishing the mathematical model for the structural parameters of soil stress ratio under complex stress path test conditions. Data analysis shows that when confining pressure is lower than the structural strength of soil, with consistent moisture content, the structural parameters of soil exhibit an exponential relationship with confining pressure; conversely, under constant confining pressure, the structural parameters of soil demonstrate a linear relationship with moisture content. This law allows for the derivation of the A and B values in the previous formula.

Table 1: The testing conditions and fitting parameters A and B.

$$[eqn]$$ $$[eqn]$$

Equation (5,6) was employed to predict the evolution model of soil structural parameters under various test conditions, as illustrated in Fig 8. The prediction model of the fitted stress ratio structural parameters is fundamentally aligned with the variations and magnitudes of the test results, regardless of moisture content, confining pressure, or shear deformation conditions.

Model fitting curves of structural parameter variation patterns: (a) σ3 = 50kPa; (b) σ3 = 100kPa; (c) σ3 = 150kPa; (d) σ3 = 200kPa.

3. Construction of the structural constitutive model

3.1. Basic theory

Referring to the modified Cam-clay model for normally consolidated soils, as shown in Fig 3. The elastic segment curve slope prior to yielding is $[eqn]$ when the consolidation pressure (p) is smaller than the initial structural strength (p1). The following equation can be used to describe the elastic segment compression curve:

$$[eqn]$$

Where $[eqn]$ is the initial void ratio; $[eqn]$ is the initial mean stress.

The soil at this point only exhibits volumetric displacement that is classified as recoverable elastic deformation:

$$[eqn]$$

In this case, $[eqn]$ represents the total volumetric strain and $[eqn]$ the elastic volumetric strain. The void ratio at which the initial structural strength is achieved is:

$$[eqn]$$

The void ratio for the soil at yielding is represented by $[eqn]$ , whereas the structurally important initial structural strength is denoted by $[eqn]$ .

Establishes a linear relationship between the void ratio of structural soil after yielding and the logarithmic value of compressive stress. At this point, the initial structural strength of the soil and the slope of the compressive linear section are both influenced by the structural parameters. The void ratio post yielding is denoted as $[eqn]$ :

$$[eqn]$$

$[eqn]$ represents the slope of the Elastoplastic deformation phase of the soil. The void ratio increment can be expressed as:

$$[eqn]$$ $$[eqn]$$

with $[eqn]$ is the void change amount in the elastically deformed portion and $[eqn]$ is the resiliency slope. Elastic volumetric strain and total volumetric strain, respectively:

$$[eqn]$$ $$[eqn]$$

Plastic compression volumetric strain:

$$[eqn]$$

The plasticizer strain $[eqn]$ as the hardening parameter is obtained by solving Equation (10):

$$[eqn]$$

3.2. Yield surface equation

The Modified Cam-clay model, derived from plastic work, demonstrates that normally consolidated soils with structural parameters held constant at a value of 1 exhibit a loading yield surface that can be described as an elliptic function. The two endpoints of this yield surface are determined by the critical state line and the isotropic compressive spherical stress, respectively. The critical state of undisturbed loess contrasts with that of disturbed loess. Research has demonstrated [18,35]that the critical state line of structural soils does not originate at the origin but intersects the P-axis at point $[eqn]$ . $[eqn]$ is a parameter associated with the structural characteristics of the soil, referred to as structural strength. Fig 9 illustrates the yield surface of structural soil. During shear deformation, as the soil’s structure damage, its structural strength will progressively diminish, resulting in a rightward shift of the soil’s yield surface.

Yield surface of structural soil.

As the structural characteristics of the soil decline the structural strength value $[eqn]$ will progressively decline from $[eqn]$ to 0. The value of the structural parameter is closely linked to the strength of the soil’s structural properties. The initial stress ratio structural parameter $[eqn]$ is maximal; as axial strain progresses, the sample experiences gradual damage, resulting in a decrease in the structural parameter’s value. The observed decline in soil structural strength aligns with the reduction in structural parameters. Consequently, the parameter $[eqn]$ can be defined as the structural attenuation coefficient using the stress ratio structural parameter $[eqn]$ and the initial stress ratio structural parameter $[eqn]$ , which fluctuates throughout the range (0,1], as seen in Equation (17). The structural strength value $[eqn]$ of the soil sample can be computed using the structural attenuation coefficient parameter, as demonstrated in Equation (18).

$$[eqn]$$ $$[eqn]$$

The isotropic compressive spherical stresses in undisturbed loess continuously fluctuate with the structural characteristics during plastic shear, and the structural critical state is likewise evolving. The critical state line of undisturbed loess under triaxial compressive stress conditions is articulated as:

$$[eqn]$$

Where: M is the slope of the critical state line for structural loess. With reference to the modified Cam-clay model derivation, the shear expansion equation for the structural constitutive model is:

$$[eqn]$$

Combining the plastic potential and yield surface, the orthogonal condition for the plastic strain increment is obtained:

$$[eqn]$$

Solving the ordinary differential equation for Equation (21), the yield function of the structural constitutive model is:

$$[eqn]$$

The yield surface equation intersects the p-axis at point $[eqn]$ , allowing for the determination of constant C, which is then substituted into the yield function to derive the yield surface equation:

$$[eqn]$$

Substitute Equation (16) into Equation (23):

$$[eqn]$$

3.3. Derivation of constitutive relationship

Tests show that following the yielding of natural structural soils, the increment of plastic strain is vertical to the plastic potential surface, allowing for the application of the associated flow law [33], which necessitates consistency between the plastic potential function and the yield function:

$$[eqn]$$

Derivative from the yield surface equation for each variable:

$$[eqn]$$ $$[eqn]$$ $$[eqn]$$

Substitute Equation (28) into Equation (25):

$$[eqn]$$

Substitute Equations (26) ~ (27) into (29). The general expression for the plastic strain of the modified Cam-clay model for structural loess is obtained:

$$[eqn]$$

According to the shear expansion Equation (20) of the structural constitutive model, the plastic shear strain increment is obtained as:

$$[eqn]$$

Where: $[eqn]$ $[eqn]$ . Organize Equations (30) ~ (31) and write them as matrix equations:

$$[eqn]$$

4. Validation of structural constitutive model

4.1. Determination of model parameters

The stress-strain equation for the final derivation of the structural properties for loess, based on the modified Cam-clay model, is of the incremental type. There are two main groups of parameters in the constitutive model. The first group is the modified Cam-clay model parameters, mainly including the critical state stress ratio M, the initial void ratio $[eqn]$ , the Poisson’s ratio $[eqn]$ , the slope of the compression curve $[eqn]$ obtained from structural loess compression tests, and the slope of the resilience curve $[eqn]$ . The second group of model parameters is related to the structural properties of loess, including the initial structural yield stress $[eqn]$ , the initial structural strength $[eqn]$ , the structural attenuation coefficient $[eqn]$ and the structural strength $[eqn]$ .

4.1.1. Determine parameters λμ(mμ), κ(mμ) and  p1(𝐦μ)

This paper incorporates the structural influence parameters of loess into the concept of stress ratio through 1-D consolidation tests. It analyzes the compression properties of undisturbed loess and disturbed loess to establish the relationship between the initial structural yield stress $[eqn]$ , the slope of the compression curve $[eqn]$ , the slope of the resilience curve $[eqn]$ , and the structural parameters. Based on the consolidation test results of loess with different moisture contents (Fig 2), Casagrande’s method was applied to determine the initial structural yield pressure P. The slope of the compression curve $[eqn]$ and the slope of the resilience curve $[eqn]$ can also be converted from the 1D consolidation test, and they are related to the compression coefficient $[eqn]$ and resilience coefficient $[eqn]$ in the 1D consolidation test as follows:

$$[eqn]$$ $$[eqn]$$

Previous studies have shown that the compressive deformation behavior of structured soils is influenced by structural changes [36, 37]. Compared with normally consolidated soils, the development of an elasto-plastic constitutive model for structured soils must account for the effects of structure on the compression curve, the slope of the resilience curve, the initial structural yield stress. Based on these considerations, and on a comparative analysis of the compressive and shear deformation characteristics of disturbed loess and structured undisturbed loess, compression tests are employed to determine the initial structural yield stress and the slopes of the linear segments of the compression and resilience curves associated with soil structure, which characterize the compressive behavior of structured soils. After conducting tests and referring to the existing literature [38–40], the parameter values for the model calculation were comprehensively determined and are shown in Table 2. The average stress ratio structural parameters corresponding to soil structural parameters at moisture contents of 14%, 17%, 20%, 23%, 26%, and 35.2% were compiled under conditions of 50kPa, 100kPa, 150kPa, and 200kPa, respectively. Figs 10–12 illustrate the slope of the compression curve $[eqn]$ , the slope of the resilience curve $[eqn]$ , the initial structural yield stress $[eqn]$ , the outcomes of the fits to the structural parameters, and a comparison of the test results.

Table 2: Relationship between model calculation parameters and structural parameters.

Relationship between λμ and 𝐦μ.

Relationship between κ and 𝐦μ.

Relationship between 𝐩1 and 𝐦μ.

According to the results in Table 2, the slope of the compression curve, the slope of the rebound curve, and the initial structural yield stress all exhibit approximately linear relationships with the structural parameter. As shown in Fig 10, the compression parameter decreases with increasing structural parameter, indicating a negative correlation. Fig 12 demonstrates that the stronger the structure of undisturbed loess, the greater the initial structural yield stress; however, an increase in moisture content significantly reduces this stress. This suggests that as moisture content rises, the bonding strength between loess particles weakens, leading to structural degradation and a reduction in structural strength, which is macroscopically reflected in the increase of the compression curve slope. Fig 11 further indicates that although the resilience parameter is also approximately linearly related to the structural parameter, the resilience deformation characteristics are not significantly affected by structure. This is mainly because the soil structure is largely destroyed during the resilience stage, resulting in nearly consistent resilience indices that can be regarded as constant.

4.1.2. Determine parameters e0 and μ.

The initial void ratio $[eqn]$ can be determined by converting the natural moisture content and relative density of the undisturbed loess, i.e.,:

$$[eqn]$$

Where $[eqn]$ is the relative density of the soil, $[eqn]$ and $[eqn]$ are the density of water and the dry density of in undisturbed loess, respectively. Poisson’s ratio of loess can be determined by Poisson’s ratio tester, take 0.25 ~ 0.4.

4.1.3. Determine parameters M and pso.

The methodology for calculating the critical state stress ratio in the loess structural constitutive model differs from that of the modified Cam-clay model. The p-q curves can be constructed based on the triaxial test results of loess samples with varying moisture content, and the equations for the critical state line of these loess samples can be derived by fitting the test data, as seen in Fig 13. The slope of the fitted curve represents the critical state stress ratio M. The initial structural strength $[eqn]$ is determined by calculating the intersection of the fitted function with the p-axis.

Critical state lines of loess with different moisture contents: (a) w = 14%; (b) w = 17%; (c) w = 20%; (d) w = 23%; (e) w = 26%; (f) w = 35.2%.

4.1.4. Determine parameters ξ and ps.

Derived from the triaxial test results of loess samples exhibiting varying moisture content and the suggested structural parameters of stress ratio. The structural attenuation coefficient $[eqn]$ is computed using Equation (17) for various confining pressures and stresses, whereas the structural strength $[eqn]$ is determined using Equation (18) and the initial structural strength $[eqn]$ . The Calculation results of model parameters are shown in Table 3.

Table 3: Calculation results of model parameters.

4.2. Secondary Development of ABAQUS and UMAT Subroutine

4.2.1. Elastoplastic matrix of the modified Cam-clay model.

ABAQUS provides a user-defined UMAT subroutine interface for the Elastoplastic constitutive model in the intrinsic model section. This study demonstrates the secondary development of the structural loess modified Cam-clay model using the Fortran language [6,41]. To prepare the user subroutine for the modified Cam-clay model of structural loess, the Elastoplastic matrix corresponding to its incremental constitutive equations must be derived. The total strain can be considered as the combination of elastic and plastic strains. The relationship between total stress and strain is as follows:

$$[eqn]$$

Elastic stress-strain relationship is used because it aligns with the modified Cam-clay model, and matrix representation is as follows:

$$[eqn]$$

According to the theory of elasticity, it is understood that there exist multiple variations of elastic constants. However, only two of these variations are considered independent, while the remaining ones can be converted into each other [42]. The elastic volumetric strain increment $[eqn]$ and the elastic shear strain increment $[eqn]$ can be expressed as respectively:

$$[eqn]$$ $$[eqn]$$

Where $[eqn]$ is the elastic matrix of the incremental constitutive relation; K is the elastic bulk modulus; and G is the elastic shear modulus. Based on the modified Cam-clay model, the elastic volumetric strain is related to the current effective stress of the soil and can be expressed as:

$$[eqn]$$

The elastic bulk modulus K is related to two other elastic parameters E (Young’s modulus) and $[eqn]$ (Poisson’s ratio) as follows:

$$[eqn]$$

Through Equations (40)~(41), E (Young’s modulus)is obtained as:

$$[eqn]$$

The following is the relationship between shear modulus G and E, $[eqn]$ :

$$[eqn]$$

The relationship with the bulk modulus is expressed as:

$$[eqn]$$

The elastic stiffness matrix is obtained by computational methods based on the generalized Hooke’s law [27]. The elastic matrix is established using the incremental constitutive relationship $[eqn]$ , as indicated in Eq. (43).

$$[eqn]$$

According to the plasticity theory, the specific equation for the Elastoplastic constitutive relation can be determined as:

$$[eqn]$$

where $[eqn]$ is the Elastoplastic matrix of the incremental constitutive relationship; A is the hardening modulus.

The structurally modified Cam-clay model incorporates the plastic volumetric strains $[eqn]$ as the hardening parameter. This model employs isotropic reinforcement and follows the associated flow law. The matrix representation of the hardening modulus A is:

$$[eqn]$$

Where $[eqn]$ ,The elasticity matrix of the structurally modified Cam-clay model is then:

$$[eqn]$$

4.2.2. Stress integration algorithm procedure.

In this study, an improved Euler integration algorithm with error control was adopted, and a UMAT subroutine for the modified Cam-clay model of structured loess, considering the effect of initial moisture content, was developed. This stress integration algorithm automatically adjusts the increment size by controlling local stress errors, effectively ensuring both computational efficiency and accuracy [43]. Fig 14 presents the flowchart of the improved Euler integration algorithm with error control, as detailed below:

Flowchart of an improved Euler integration algorithm with error control.

(1) For a given increment, a fictitious time $[eqn]$ is introduced, with a fictitious time increment $[eqn]$ . At the beginning of the increment, $[eqn]$ and $[eqn]$ are defined, and the initial stress and hardening parameters are denoted as $[eqn]$ and $[eqn]$ , respectively. The sub-incremental strain $[eqn]$ is determined through calculation by $[eqn]$ .(2) Based on the determined sub-step strain increment $[eqn]$ , the first trial stress increment is calculated using $[eqn]$ , and the corresponding increments of plastic strain and hardening parameters $[eqn]$ are obtained.(3) Based on the stress $[eqn]$ obtained from the first trial, a new elasto-plastic stiffness matrix $[eqn]$ is determined. Using the same procedure as in the first trial, the stress increment $[eqn]$ from the second trial is computed, and the average stress increment after the two trials is taken as $[eqn]$ .(4) The relative error R is calculated based on $[eqn]$ . If R exceeds the error control threshold SSTOL, the time increment is adjusted to $[eqn]$ , and $[eqn]$ is used as the new sub-step increment to return to Step 2. If $[eqn]$ , the stress components, plastic strains, and hardening parameters are updated.(5) Based on $[eqn]$ , it is checked whether the updated $[eqn]$ exceeds 1. If the updated $[eqn]$ is greater than 1, $[eqn]$ is set, and the calculation returns to Step 2. The computation stops when $[eqn]$ is satisfied.

4.3. Validation of the model

The accuracy of the modified Cam-clay model algorithm is confirmed by the utilization of ABAQUS finite element analysis software. Initially, a 3D finite element model is created for a cylinder with a diameter of 39.1 mm and a height of 80 mm. The dimensions are identical to those of the triaxial samples from the laboratory test. The laboratory test of loess samples with varying moisture content documented the progression of axial strain from 0.5% to 15%, getting a total of 30 data points. Consequently, the modeling was segmented into 30 analytical steps, each corresponding to a distinct loading process, to facilitate comparison with the laboratory results. The model’s initially and boundary conditions are established based on the actual environmental circumstances encountered by the specimen during the laboratory test. As the triaxial test just monitors axial deformation of the material, the model computes solely axial deformation while disregarding lateral deformation. The initial condition of the Abaqus model refers to the state of the specimen prior to testing, characterized by none of load and deformation, that is, a stress condition of zero and an axial deformation of zero. The sample’s loading procedure consists of two steps: the consolidation step and the shear step. In the consolidation step, the stress conditions of the model matched the confining pressure, with four groups established at 50 kPa, 100 kPa, 150 kPa, and 200 kPa, and the initial axial deformation was zero. In the shear step, the lateral stresses of the model remain constant, the axial deformation aligns with the examined step, and the axial stresses are derived values. The test data with initial moisture contents of 14%, 20%, and 26% were selected as the control group for this paper’s model, and the prediction model, modified Cam-clay model, and test results are illustrated in Fig 15.

Comparison of predicted and test values of the loess constitutive model: (a) w = 14%; (b) w = 20%; (c) w = 26%.

The figure illustrates that the sample with 14% moisture content exhibits softening damage in the deviatoric stress-strain curve, and the model’s predicted values align closely with the test values, indicating that the structural model developed in this paper accurately represents the strain softening process of the soil samples. Fig 15(a) illustrates the results of the comparison between the numerical values of the structural model and the modified Cam-clay model with the experimental values. The computed values of the structural model are more aligned with the experimental values compared to those of the modified Cam-clay model, which only represents the compressive stiffness of the soil and fails to account for the softening behavior of the soil samples. A comparative analysis of samples with 20% and 26% moisture content shows that the deviatoric stress-strain curve exhibits hardening damage. The model’s calculated values, established in this paper, closely align with the test values, effectively representing the compression hardness and shear shrinkage of the soil samples, thereby demonstrating the model’s validity.

To verify the applicability of the proposed structural constitutive model, experimental data from loess specimens reported by Hou [44] in Jingyang Malan loess (Case 1), Gao et al. [45] in Yan’an New District loess (Case 2), and Wang [46] in Xi’an Xincheng District loess (Case 3)were adopted. Numerical simulations under different confining pressures were conducted in ABAQUS and compared with the experimental results. The initial void ratio $[eqn]$ , Poisson’s ratio $[eqn]$ , compression parameter $[eqn]$ , resilience parameter $[eqn]$ , and initial structural yield stress $[eqn]$ were obtained from the literature. The structural parameters of loess were calculated using the method proposed in this study, including the critical state stress ratio M and the initial structural strength $[eqn]$ . The model parameters are summarized in Table 4, and the comparison between the calculated and experimental results is shown in Fig 16.

Table 4: Model parameters adopted from literature specimens.

Comparison of calculated and experimental values for literature specimens: (a) Case 1; (b) Case 2; (c) Case 3.

Fig 16 shows the comparison between the simulated and experimental results of specimen models from the selected literature. It can be observed that under different confining pressures, the simulated values obtained from the proposed model are in close agreement with the experimental data, accurately capturing the variation of deviator stress with axial strain. Specifically, for the specimen under a confining pressure of 100 kPa (Fig 16a), the deviatoric stress–strain curve exhibits softening behavior, and the model reproduces the stress peak, effectively describing the softening phenomenon. In contrast, for the specimens under other confining pressures (Figs 16b and 16c), the deviatoric stress–strain curves display hardening behavior, and the discrepancies between simulated and experimental results are minimal. These findings indicate that the proposed model provides a reliable description of this characteristic of loess. In summary, the model developed in this study also shows good applicability to loess specimens from other regions.

5. Conclusions

(1) This work comprehensively examines the impact of initial moisture content on the structural characteristics of loess, enhancing the isotropic consolidation parameters and the yield surface equation of the modified Cam–clay model. A constitutive model for structural loess that effectively encapsulates the complete process of structural development is constructed. In contrast to the traditional modified Cam–clay model, the suggested model incorporates four supplementary structural factors, each with clear physical significance: initial structural yield stress $[eqn]$ , the initial structural strength $[eqn]$ , the structural attenuation coefficient $[eqn]$ and the structural strength $[eqn]$ . They provide a more accurate depiction of the structural development of undisturbed loess and establish a dependable theoretical basis for enhancing the quantitative characterization of loess structure and refining engineering assessments in loess regions.(2) The relevant structural parameters are added to the modified Cam-clay model, while the consistency condition and the associated flow law are used to establish an incremental stress-strain relationship. The modified Cam-clay model was established in this study, taking into account the influence of loess structural properties. When the loess structural properties weren’t considered, the model described in this study became a modified Cam-clay model.(3) The parameters of the proposed model can be derived from standard laboratory testing or triaxial tests. Utilizing the secondary development interface of ABAQUS, the modified Cam–clay model is enhanced and expanded to provide a user-defined subroutine for the structural constitutive model of loess. Comparative examinations of the numerical findings and experimental data of loess specimens with varying starting moisture levels validate the model’s dependability. In contrast to the traditional modified Cam–clay model, the structural model formulated in this work not only encapsulates the compression hardening attributes of loess but also proficiently depicts its dilatancy behavior. It offers a more precise characterization of the impact of soil structure on the mechanical behavior of loess, exhibiting significant logic and engineering relevance.

Supporting information

S1 FileThe figures and tables have been uploaded to the Figshare repository, accessible via: 10.6084/m9.figshare.30239785.(ZIP)