Resilient modulus and cumulative plastic strain of frozen silty clay under dynamic aircraft loading

This paper describes an investigation into the factors influencing the resilient modulus and cumulative plastic strain of frozen silty clay. A series of dynamic triaxial tests are conducted to analyze the influence of the temperature, confining pressure, frequency, and compaction degree on the resilient modulus and cumulative plastic strain of frozen silty clay samples. The results show that when the temperature is below − 5 °C, the resilient modulus decreases linearly, whereas when the temperature is above − 5 °C, the resilient modulus decreases according to a power function. The resilient modulus increases logarithmically when the frequency is less than 2 Hz and increases linearly once the frequency exceeds 2 Hz. The resilient modulus increases as the confining pressure and compaction degree increase. The cumulative plastic strain decreases as the temperature decreases and as the confining pressure, frequency, and compaction degree increase. The research findings provide valuable information for the design, construction, operation, maintenance, safety, and management of airport engineering in frozen soil regions. The empirical equation for resilient modulus and temperature is a segmented function; the segmentation point is − 5 °C. The empirical equation for resilient modulus and frequency is a segmented function; the segmentation point is 2 Hz. The cumulative plastic strain is positively correlated with temperature and negatively correlated with frequency. The empirical equation for resilient modulus and temperature is a segmented function; the segmentation point is − 5 °C. The empirical equation for resilient modulus and frequency is a segmented function; the segmentation point is 2 Hz. The cumulative plastic strain is positively correlated with temperature and negatively correlated with frequency.


Introduction
The pavement stability of airport engineering is highly dependent on the resilient modulus and cumulative plastic strain of subgrade soils [1,2]. The resilient modulus and cumulative plastic strain are important parameters in reflecting the mechanical properties of such soils [3][4][5][6]. Many scholars have studied the resilient modulus and cumulative plastic strain of unfrozen soils. For instance, Xenaki and Athanasopoulos [7] investigated the influence of the confining pressure, moisture content, and cyclic strain amplitude on the resilient modulus. Guisasola et al. [8] studied the effect of soil category, moisture content, and stress history on the resilient modulus, and found that increasing the moisture content helped to weaken the resilient modulus. Chen et al. [9] conducted dynamic triaxial tests to analyze the resilient modulus of sandy pebble soil. Their results indicate that the growth of dynamic strain results in reduced resilient modulus, whereas the addition of consolidation pressure increases the resilient modulus. Wang et al. [10] stated that the resilient modulus increased with increases in frequency, consolidation ratio, and confining pressure and a decrease in moisture content. Guo et al. [11] researched the resilient modulus of intact soft clay in Wenzhou, China, and concluded that the confining pressure had a significant effect on the resilient modulus. Bao and Mohajerani [12] believed that the resilient modulus increased as the deviatoric stress decreased, and showed that the deviatoric stress and moisture content greatly affected the resilient modulus. Han and Vanapalli [13] established three equations for predicting the resilient modulus of unsaturated subgrade soils. Wang et al. [14] divided the cumulative plastic strain into reversible and irreversible components, and established a prediction model for saturated soft clay using cyclic loading tests. Chen et al. [15] used a cyclic triaxial apparatus to investigate the effects of stress, cycle number, and the initial deviatoric stress on the cumulative plastic strain of soft clay from Shanghai. Xiao and Liu [16] analyzed the influence of the compaction degree, moisture content, and dynamic stress level on the cumulative plastic strain of silt soils. Zhang et al. [17] used long-term cyclic triaxial tests to study the influence of the moisture content, number of loading cycles, deviatoric stress, and confining pressure on the cumulative plastic strain, and used the results to construct a prediction model for subgrade soils.
Several scholars have also studied the resilient modulus and cumulative plastic strain of frozen soils. Ling et al. [18] investigated the effect of temperature, confining pressure, and moisture content on the resilient modulus of frozen soil, and stated that the resilient modulus increased as the temperature decreased and the confining pressure increased. The moisture content had little influence on the resilient modulus when it was greater than 0.21%. Zhu et al. [19] used subgrade soil from the Qinghai-Tibet Railway to explore the effect of temperature, moisture content, frequency, and confining pressure on the cumulative plastic strain. Zhao et al. [20] found that the cumulative plastic strain varied little with different dynamic stress amplitudes based on dynamic creep tests. Tang et al. [21] analyzed the resilient modulus and cumulative plastic strain of silty clay soils with different freezing temperatures, dynamic stress amplitudes, and loading frequencies. They found that the freezing temperature had a small influence on the resilient modulus, and that freeze-thaw action was beneficial to weakening the resilient modulus. Zhang et al. [22] studied the influence of temperature, confining pressure, frequency, and dynamic deviatoric stress on the resilient modulus and cumulative plastic strain of frozen aeolian soils.
There have been many studies investigating the influence of frequency, confining pressure, and dynamic deviatoric stress on the resilient modulus and cumulative plastic strain of unfrozen soils. Relatively few scholars have researched the effect of temperature, frequency, confining pressure, and dynamic deviatoric stress on the resilient modulus and cumulative plastic strain of frozen soils. Some studies have investigated the frequency, confining pressure, and dynamic deviatoric stress of subgrade soils with respect to an aircraft taxiing on the pavement, and there has been some consideration of the effect of the compaction degree on the subgrade in airport engineering.
In this paper, based on the heat balance control equation and dynamic aircraft load equation, an aircraft-pavement model is established to analyze the range of temperature, frequency, and stress of subgrade soils in permafrost regions. This model provides the testing parameters for triaxial tests of frozen soil, which are conducted to analyze the influence of temperature, confining pressure, frequency, and compaction degree on the resilient modulus and cumulative plastic strain. Finally, the relationship between temperature, confining pressure, frequency, compaction degree, and resilient modulus is established and discussed. The research findings provide valuable information for the design, construction, operation, maintenance, safety, and management of airport engineering in frozen soil regions, especially frozen silty clay regions. following practical assumptions: (1) the pavement structure layer, graded gravel layer, and subgrade are uniform and isotropic; (2) moisture migration can be ignored at the boundary of each structure layer; (3) the subgrade follows the energy conservation law and only includes heat transfer effects [23].
where ρ is the density, C is the specific heat, T is the temperature, t is the time, λ is the thermal conductivity coefficient, and x, y are the rectangular coordinates.
Based on different freezing and thawing states, λ and C can be written as follows: (1) where T m is the temperature of the phase transformation, ΔT is the temperature variation of the phase transformation, L is the latent heat, λ f is the thermal conductivity coefficient of the frozen soil, C f is the specific heat of the frozen soil, λ u is the thermal conductivity coefficient of the unfrozen soil, and C u is the specific heat of the unfrozen soil.
There are three kinds of boundary conditions according to heat-conduction theory [24]. The first boundary condition is established based on the monitoring temperature of the permafrost region in Northeast China: The second boundary condition is the adiabatic boundary condition, which is applied on both sides of the numerical model: The third boundary condition is the constant temperature gradient of 0.03℃/m, which is applied at the bottom of the numerical model:   where T 0 is the boundary layer temperature of 4.5 °C, A is the amplitude of the upper boundary temperature, which is 23 °C for natural surfaces and 25 °C for concrete pavement, t is the time, b is the initial phase of π/2, and a is the annual average temperature rise of 0.03 °C. The numerical model of the pavement structure in airport engineering is illustrated in Fig. 1 [3,4]. The boundary conditions of the numerical model are given by Eqs. (4) and (5). The bottom of the numerical model is assigned a constant temperature gradient of 0.03 °C/m, as shown in Eq. (6). The apparent heat capacity method is used to simulate the phase change occurring in a certain temperature range. The thermal parameters of the base layer, cushion layer, graded gravel layer, silty clay, and weathered rock are listed in Table 1 [24].
The temperature field in Northeast China was first calculated for 20 years. The temperature field over these 20 years was then taken as the initial temperature to calculate the 30-year temperature field, as shown in Fig. 2.
As shown in Fig. 2, for a given subgrade depth, the temperature in March increases with the thickness of the graded gravel layer. For a given graded gravel layer thickness, the temperature in March increases with the depth of the subgrade. The temperature of the 0-m graded gravel layer with a 0.2-m cushion layer is higher than that with a 0-m cushion layer. The temperature of the 0-m graded gravel layer with no cushion layer varies from 0 to − 7 °C. Because the critical freezing temperature of the silty clay is about − 1 °C [25], the triaxial tests considered temperatures from − 1 to − 9 °C.

Frequency parameters
The dynamic equation of aircraft loading are as follows [26]: where P v is the load of the aircraft when taxiing, M is the mass of the aircraft, g is the acceleration due to gravity, and c 0 is a coefficient, set to 10 -3 m −0.5 s 0.5 . IRI denotes the international roughness index, v is the taxiing speed, and Y v is the elevating force.
When the aircraft has left the ground, P v = 0.
where v 0 and Y v0 denote the speed and elevating force when the aircraft is away from the ground, respectively. The dynamic load coefficient can be obtained by Eqs. (7) and (8) as: where k is the dynamic load coefficient. Based on the parameters of a B737-800 [26], the aircraft is constructed using beam elements (MPC184). The suspension and non-suspension systems are represented by mass elements (Mass21), as are the inertia moments of pitching and rolling rotations. The spring and damper are represented by a spring-damper element (Combine14). Transient dynamic analysis was then used to determine the aircraft dynamic load, as shown in Fig. 3. Figure 3 shows that the aircraft is static and k = 1 when the taxiing speed is 0 m/s. The aircraft leaves the ground and k = 0 when the taxiing speed reaches 69.2 m/s. k increases and then decreases as the taxiing speed rises from 0 to 69.2 m/s. For a given taxiing speed, k increases as IRI increases. The maximum value of k occurs at a taxiing speed of 10 m/s. As IRI increases from 1 to 6, the variation of k is less than 0.1 for a taxiing speed of 10 m/s. Based on the technical specifications of aerodrome pavement evaluation and management, the roughness level is good when IRI is less than 2 and poor when IRI is greater than 4. Because most airports have a medium roughness level (IRI between 2.0 and 4.0), the dynamic aircraft load was The elastic-plastic constitutive model and Mohr-Coulomb yield criteria are used to establish the numerical model of the subgrade. The elastic modulus, Poisson's ratio, cohesive force, and internal friction angle of frozen soil can be calculated from the relationship between temperature and the mechanical parameters, as shown in Eqs. (10)-(13) [27]. The pavement structure parameters are listed in Table 2 [27].

Fig. 3 Dynamic load coefficient influenced by IRI and v
where E is the elastic modulus, ν is Poisson's ratio, c is the cohesive force, ϕ is the internal friction angle, and m is a nonlinear index set to 0.6. a i and b i are mechanical parameters, as listed in Table 3.
Based on the aircraft-pavement model, the frequencies at the top of the subgrade were studied for taxiing speeds from 5 to 65 m/s, as shown in Fig. 4. This figure indicates that the frequency at the top of the subgrade is between 1.08 and 4.92 Hz. Therefore, the triaxial tests considered loading frequencies from 0.5 to 5 Hz.

Stress parameters
The stress of the subgrade was calculated as follows [25]: where P 0 is the vertical stress affected by the weight of overlying structures, and h i , γ i are the thickness and the bulk density of the i-th layer, respectively.  where σ 1 is the total vertical stress of the calculated point and σ z is the vertical stress affected by the weight of the aircraft load.
where σ is the confining pressure of the calculated point, σ x,y is the mean lateral stress affected by the weight of the aircraft load, and k 0 is the lateral pressure coefficient, which is equal to 0.85. Based on the above theory and aircraft-pavement model, the influence of the pavement structure on the (16) = x,y + k 0 P 0 stress of the subgrade was studied. The results are presented in Tables 4 and 5.
Tables 4 and 5 show that σ 1 and σ do not vary significantly at 5.0 m away from the top of the subgrade, but do vary significantly at the top of the subgrade. The reason is that P 0 increases and σ z decreases as the depth of the subgrade increases. Moreover, P 0 is only affected by the thickness and the bulk density of the structure layers. Tables 4 and 5 also indicate that σ 1 varies from 0.121 to 0.365 MPa, and σ varies from 0.100 to 0.269 MPa. Therefore, the triaxial tests considered confining pressures from 0.1 to 0.4 MPa.

Materials and methods
The triaxial apparatus used in this study is shown in Fig. 5. It mainly includes a fuel supply system, triaxial test system, sample incubators, loading system, liquid injection cooling circulating system, and numerical control servo system. The maximum values of the axial load and confining   Table 6. The silty clay was divided into cylindrical samples measuring 61.8 mm in diameter and 125 mm in length, as shown in Fig. 6. The sample preparation process included four steps. First, the silty clay was crushed and placed in a drying oven for 12 h at a temperature of 105 °C. Some dry silty clay and distilled water were then mixed to produce wet soil with a predefined moisture content. The wet soil was stored in an airtight container for 24 h to ensure the uniform distribution of water in the silty clay. Second, the quantity of every sample was calculated from the dry density and the moisture content. The wet soil sample was then layered into the mold to ensure it was homogeneous. To check the uniformity of samples, the completed sample was cut into four pieces and the dry density and moisture content were measured. Third, the completed samples were installed in the triaxial pressure chamber and frozen at the desired temperature for 24 h by a refrigerated circulator machine. Fourth, the frozen samples were subjected to dynamic loading. Every test result presented in this paper is the average of three samples with the same test parameters.
Previous research has indicated that the dynamic stress form of the subgrade soil is similar to half of a sine wave. Thus, the sine wave of the unidirectional cyclic load was used to simulate the unidirectional impact stress of the dynamic aircraft load. The detailed parameters of the test conditions are listed in Table 7. The termination index of the dynamic triaxial test was a cycle number of 10,000 or a cumulative strain threshold of 10%.
The dynamic resilient modulus is defined as: where E d is the dynamic resilient modulus, σ d is the applied deviatoric stress, and ε r is the recoverable strain during each loading cycle. The cumulative plastic strain is the irrecoverable strain after the action of dynamic stress. In this study, the cumulative plastic strain is defined as the ratio between the height variation of the sample and the original height of the sample after the action of dynamic stress. Figure 7 shows that the resilient modulus decreases with an increase in temperature. The rate of change in the resilient modulus at temperatures from − 5 to − 9 °C is smaller than that from − 1 to − 5 °C. The reason is that the ice is very sensitive to temperature changes in frozen silty clay. Temperatures from − 5 to 0 °C are the main revulsion area of the ice-water phase change. The weakly binding water that exists in the soil particle surface affects the cementation of the soil and the ice. Moreover, the molecular activity of the ice decreases with any reduction in temperature. The molecular activity of the ice increases the compactness of the soil structure and the value of the resilient modulus. Based on the above test data, the relationship between the temperature and the resilient modulus can be written as follows:

Influence of temperature on resilient modulus
where E d is the dynamic resilient modulus.

Influence of confining pressure on resilient modulus
Because the confining pressure varies with the depth of the subgrade, this parameter has a remarkable effect on the resilient modulus. Figure 8 indicates that the resilient modulus increases with the growth of the confining pressure. The reason is that an increase in the confining pressure enhances the lateral restraint and strength of the frozen silty clay. Moreover, when the confining pressure reaches a certain level, the ice in the voids melts. This weakens the cohesion between the ice and the soil in the voids, which induces a lubrication action that decreases the strength of the frozen silty clay. Based on the above test data, the relationship between the confining pressure and the resilient modulus can be written as follows:    Figure 9 illustrates that the resilient modulus increases as the frequency becomes higher. When the frequency is less than 2 Hz, the rate of change in the resilient modulus is large. The viscoplasticity of the frozen silty clay represents the lower frequency range. However, the viscoplasticity cannot fully represent the behavior of the soil as the frequency increases. Based on the above test data, the relationship between the frequency and the resilient modulus is as follows:

Influence of frequency on resilient modulus
where f is the frequency.

Influence of compaction degree on resilient modulus
The compaction degree is an important index reflecting the performance of the subgrade in airport engineering, and is an indispensable influence factor for the resilient modulus of the subgrade. Figure 10 suggests that the resilient modulus increases with the growth of the compaction degree. Because the growth of the compaction degree reduces the porosity of the frozen silty clay, the structure of the frozen silty clay tends to become more stable, and so the resilient modulus increases. Based on the above test data, the relationship between the compaction degree and the resilient modulus can be written as: where K is the compaction degree and K 0 is the reference compaction degree (taken as 0.80 in this study).  Figure 11 shows that the cumulative plastic strain increases with an increase in temperature. When the temperature is from − 5 to − 9 °C, the cumulative plastic strain increases and tends to remain stable with any increase in the cycle number. However, at temperatures from − 1 to − 3 °C, the cumulative plastic strain increases as the cycle number rises. The reason for this is that successive cycles result in the structural failure of the frozen silty clay at temperatures from − 1 to − 3 °C. When the temperature is below − 5 °C, more weakly bound water turns into ice, increasing the cohesive force between soil particles and ice, and so the elastic modulus of the soil increases to reduce the cumulative plastic strain. Moreover, the structure of the soil sample becomes more compact with successive cycles. In short, a decrease in temperature has a significant effect on decreasing the cumulative plastic strain and increasing the strength of frozen silty clay. Fig. 11 Influence of the temperature on the cumulative plastic strain

Influence of confining pressure on cumulative plastic strain
It can be seen from Fig. 12 that the cumulative plastic strain decreases as the confining pressure increases from 0.1 to 0.4 MPa. For a given confining pressure, the rate of increase in cumulative plastic strain decreases with successive cycles. The cumulative plastic strain increases rapidly when the confining pressure is small. The confining pressure produces an extrusion effect and local stress concentration in the soil particles, causing the ice to melt and the weakly bound water to migrate. The weakly bound water migration reduces the friction of the soil particles and the cohesive force between soil particles and ice crystals. This aids the rearrangement of soil particles and ice crystals to enhance the structural stability and improve the ability to resist deformation. Moreover, an increase in the confining pressure increases the lateral deformation limitation, thus enhancing the ability to resist deformation. Figure 13 indicates that the attenuation of the cumulative plastic strain at frequencies of 0.5-2 Hz is smaller than that at frequencies of 2-5 Hz. When the frequency is small, the voids in the frozen silty clay have enough time to develop and destroy the soil structure. Hence, the strength of the frozen silty clay decreases and the cumulative plastic strain increases. However, frozen silty clay has remarkable instantaneous strength, and the voids have little time to grow and damage the structure of the frozen silty clay. Thus, the cumulative plastic strain is generally small.

Influence of compaction degree on cumulative plastic strain
As shown in Fig. 14, the cumulative plastic strain increases with successive cycles when the compaction degree is from 0.80 to 0.95. Moreover, the cumulative plastic strain decreases as the compaction degree increases. As the compaction degree increases, the void fraction decreases, which improves the stabilization of the soil sample. Thereby, the cumulative plastic strain is reduced.

Discussion
Tretovich divided frozen soil into three main ice-water phase transition zones. (1) Intense phase transition zone: As the temperature decreases by 1 °C, the changes in the unfrozen water and ice content are greater than or equal to 1%. At this point, all the free water in the macropores and capillaries and some of the weakly bound water in the soil have frozen.
(2) Transition zone: As the temperature decreases by 1 °C, the changes in the unfrozen water and ice content are between 0.1 and 1%. At this point, the weakly bound water in the soil has frozen. (3) Frozen zone: As the temperature decreases by 1 °C, the changes in the unfrozen water and ice content do not exceed 0.1%. At this point, only strongly bound water remains in the soil, and the strongly bound water content is close to the maximum moisture absorption of the soil. When the temperature is from − 5 to 0 °C, the dynamic resilient modulus and the cumulative plastic strain have large variation ranges. When the temperature is below − 5 °C, the dynamic resilient modulus and the cumulative plastic strain have small variation ranges. This phenomenon is similar to the findings of Zhu [28], Shi [29], and Wang [30]. Because ice is very sensitive to temperature changes in the permafrost, the temperature range from − 5 to 0 °C is the main intense phase transition zone. When the temperature is between − 5 and 0 °C, part of the weakly binding water exists in the soil particle surface, and affects the cementation of ice and soil. The molecular activity of the ice decreases with the reduction of the temperature, and the structure of the soil becomes more compact. This causes the dynamic resilient modulus of the soil to increase and the cumulative plastic strain to decrease. An increase in the confining pressure improves the lateral deformation limitation, which enhances the ability to resist deformation. Moreover, a higher confining pressure causes the ice to melt and the weakly binding water to migrate, which enables the rearrangement of soil particles and ice crystals so as to enhance the structural stability and the ability to resist deformation. Hence, the dynamic resilient modulus of the soil increases and the cumulative plastic strain decreases as the confining pressure increases, which confirms the findings of Shi [29] and Fan [31]. Compared with the research of Shi [29], there is a critical frequency. When the frequency is less than this critical frequency, the frozen soil exhibits viscoplasticity. The viscoplasticity of the frozen soil can be fully exerted. The dynamic elastic modulus of the frozen soil increases greatly and the cumulative plastic strain decreases greatly as the frequency increases. When the frequency is above the critical frequency, the frozen soil exhibits elasticity. The dynamic elastic modulus of the frozen soil increases slightly and the cumulative plastic strain decreases slightly as the frequency increases.
Because the porosity of the frozen soil decreases with the growth of the compaction degree, the cohesive force and friction force between soil particles increases to enhance the dynamic resilient modulus and reduce the cumulative plastic strain. In addition, Monismith used silty clay as an experimental material to carry out indoor triaxial compression experiments, and proposed a prediction model for the cumulative plastic strain [28,29]. The Monismith model shows that the cumulative plastic strain increases exponentially with the growth of the cycle number [28,29], which is similar to the findings presented in this paper.

Conclusions
Based on the heat balance control equation and the dynamic aircraft load equation, an aircraft-pavement model has been established to analyze the effects of temperature, frequency, and stress on the silty clay soil in the permafrost region of Northeast China. This model was used to determine the parameters for triaxial tests of the subgrade soil. Numerous dynamic triaxial tests were conducted to analyze the influence of temperature, confining pressure, frequency, and compaction degree on the resilient modulus and cumulative plastic strain. The results can be explained as follows: (1) Temperature is an important factor for the resilient modulus and the cumulative plastic strain of frozen silty clay. When the temperature varies from − 1 to -5 °C, the resilient modulus increases significantly and the cumulative plastic strain decreases significantly with successive compression cycles. When the temperature is from − 5 to − 9 °C, the resilient modulus increases slowly and the cumulative plastic strain decreases slowly with successive cycles. (2) As the confining pressure increases, the resilient modulus increases while the cumulative plastic strain decreases. The reason is that an increase in the con-fining pressure enables the rearrangement of soil particles and ice crystals, which enhances the soil's resistance to deformation. (3) The frequency is one of the main factors affecting the resilient modulus and the cumulative plastic strain. The resilient modulus increases faster at frequencies from 0.5 to 2 Hz than at frequencies from 2 to 5 Hz. The cumulative plastic strain decreases significantly when the frequency is lower than 2 Hz. (4) The compaction degree is another important factor affecting the resilient modulus and the cumulative plastic strain. An increase in the compaction degree reduces the porosity and enhances the strength of the sample structure, thus increasing the resilient modulus and decreasing the cumulative plastic strain. (5) This study has only considered silty clay. Future research will examine the resilient modulus and cumulative plastic strain of clay and silty soil, and will attempt to determine empirical equations for the resilient modulus and the cumulative plastic strain.