Constitutive Model of Lateral Unloading Creep of Soft Soil under Excess Pore Water Pressure

School of Civil Engineering and Architecture, Chongqing University of Science and Technology, Chongqing 401331, China Chongqing Key Laboratory of Energy Engineering Mechanics & Disaster Prevention and Mitigation, Chongqing 401331, China Department of Civil and Environmental Engineering, Jackson State University, Jackson, MS 39217, USA CREEC (Chongqing) Survey, Design and Research Co., Ltd., Chongqing 401121, China Shenzhen Water Planning & Design Institute Co., Ltd., Shenzhen 518000, China Department of Civil and Environmental Engineering, Jackson State University, Jackson, MS 39217, USA Department of Civil and Architectural Engineering, Tennessee State University, Nashville, TN 37209, USA


Introduction
Soft soil, being a special kind of regional soil with high water content, low bearing capacity, and salient creep characteristics, is widely distributed in Northern Europe, North America, and coastal areas of Pearl River Delta and Yangtze River Delta of China [1][2][3]. Currently, with the rapid development of urban construction, many soft soil underground spaces such as foundation pits and subways in the above areas are in the development phases [4]. e engineering practices of underground space excavation in soft soil areas showed that, during the excavation process, soft soil surrounding the rock mainly experiences unloading effect. e deformation and failure of soft soil surrounding the rock will take a period, and it generally means that the soft soil unloading creep failure occurs [5]. In addition, due to the relatively high groundwater level in soft soil areas, the construction vibration of large equipment (such as shield, Figure 1) will bring 21.2 to 56.12 kPa of excess pore water pressure to the geological environment of soft soil with a buried depth of 12 to 15 m [6]. As shown in Figure 1, the excess water pressure weakens the mechanical properties of soft soil, increases the probability of unloading creep failure occurrence of soft soil, and even causes rheological disasters. erefore, the study of the constitutive model of unloading creep of soft soil under excess pore water pressure is of great significance to the numerical analysis and excavation optimization of underground space excavation in water-rich soft soil areas.
According to its stress path, the creep constitutive model of soft soil can be divided into the loading creep constitutive model and the unloading creep constitutive model [7]. e creep constitutive models under loading were studied previously, and research results were applied to the numerical analysis of unloading excavation of the underground space in soft soil, resulting in a large gap between the calculated results and the actual situation [8,9]. erefore, it had become a common practice among scholars to study unloading excavation of the underground space of soft soil by constructing the unloading creep constitutive model based on the unloading stress path. Later, some scholars carried out one-dimensional unloading-reloading tests on soft soil and proposed a one-dimensional elastic-viscoplastic constitutive model (1D EVPS) which takes creep/expansion effect into account and a three-element rheological model that can be used to describe the nonlinear creep of soft soil [10,11]. Feng et al. [12] carried out one-dimensional loading and unloading creep tests on marine sedimentary soft soil and established a modified EVPS model to describe the creep deformation. Study results show that, as one-dimensional consolidation creep cannot be fully developed, the creep rate decreases with the increase in soft soil creep degree. With a goal of revealing the triaxial creep law of soft soil more clearly, Zhou and Chen [13], based on the triaxial creep test results of Guangzhou soft soil, established the unloading creep constitutive model of soft soil using the merchant model. Zhou and Chen [14] established an integral constitutive model of soft soil under unloading creep. Although these unloading creep constitutive models take the anisotropy and structure of soft soil into account [11,15], more unloading creep lies in the one-dimensional model, and further study of the triaxial unloading creep constitutive model of soft soil is needed.
In the construction of the underground space in soft soil, excess pore water pressure is generated in large amount, and creep rate of soft soil surrounding the rock is significantly affected [16][17][18]. Currently, research studies on the excess pore water pressure mainly focus on the variation rule and dissipation time of excess pore water pressure on the construction site, e.g., the excess pore water pressure caused by the shield construction vibration is about 21.2 to 56.12 kPa, and the dissipation time of excess pore water pressure is 20 to 40 h [6,19]. Previous studies indicated that the excess pore water pressure weakens the mechanical properties of the soil and causes the peak shear strength and strength parameters to decrease with the increase in the excess pore water pressure [20][21][22]. However, the functional relationship between the excess pore water pressure and the unloading creep constitutive model of soft soil has not been established, which limits the engineering application of the unloading creep constitutive model in underground space excavation analysis of water-rich soft soil. e above analysis indicates that the existing research either focused only on the establishment of the one-dimensional unloading constitutive model or on the variation law of excess pore water pressure. However, there is insufficient research on the soft soil unloading constitutive model under the combined effect of excess pore water pressure and lateral unloading path, which is essentially one of the most important problems in the numerical analysis of the water-rich soft soil underground space. erefore, in this study, the muddy soft soil in the Shenzhen area of China was selected as the research object, the unloading creep characteristic test research studies of soft soil under different excess pore water pressures were carried out, and the unloading creep constitutive model considering excess pore water pressure was established. e results in this study will provide references for the design and numerical analysis of the underground space of water-rich soft soil.

Preparation of the Sample.
e TSZ-2S fully automatic stress-strain control triaxial instrument was used in this study. e sample is taken from a large foundation pit site in Shenzhen of China. e intact sample has a density of 1.82 g/ cm 3 , a saturated water content of 39.6%, a specific gravity of 2.73, a void ratio of 1.094, and a plasticity index of 16.3. e remolded soil samples were prepared with a height of 80.0 mm and a diameter of 39.1 mm according to Geotechnical Test Regulation (SL237-1999). Based on the triaxial tests of isotropic consolidation, the measured cohesion and internal friction angle of the remolded soil sample were 19.95 kPa and 28°, respectively.

Lateral Unloading Stress Path and Excess Pore Water
Pressure. Excavation of the underground space in soft soil involves a variety of unloading stress paths, wherein the lateral unloading stress path plays an important role in the study of lateral deformation and stability of the underground space in soft soil [18]. erefore, this paper mainly focuses Excess pore water pressure = 21.2 -56.12kPa on the creep constitutive model under the lateral unloading stress path. In the test, the backpressure in the field is increased to simulate the excess pore water pressure. Four different excess pore water pressures u, 0, 20, 40, and 60 kPa, are set to study the effect of excess pore water pressure on lateral unloading creep.

Saturation of Soil Samples.
e backpressure saturation system of the triaxial instrument is used to saturate the sample. During the test, under the conditions that the confining pressure σ 3 is 100 kPa and the backpressure u is 90 kPa, Skempton's B-value of the soil reaches 0.95 after 24 hours, which means that the soil is saturated. Refer to Huang et al. [18] for the specific saturation process.

K0 Consolidation.
e sampling depth is 10.5 m, and the corresponding consolidation confining pressure is about 100 kPa. In order to simulate the gravity stress of the soil, K0 consolidation is carried out under the condition that the confining pressure σ 3 is 100 kPa and the axial pressure σ 1 is 189 kPa.

Application of Excess Pore Water Pressure. After K0
consolidation is completed, the following excess pore water pressures u, 0, 20, 40, and 60 kPa, are applied to the inside of the soil sample through the backpressure system of the triaxial instrument, respectively.

Unloading Creep Test.
e axial load is kept unchanged, and the lateral unloading creep test is carried out in 4 to 7 grades with no drainage. e lateral unloading rate of each stage is 1 kPa/min. After unloading of each stage, the confining pressure and axial pressure remain unchanged for 3 to 4 h. At the time when the axial deformation is less than 0.002 mm within 1 h, the next unloading shall be carried out until the soil sample is destroyed. e detailed unloading process is shown in Table 1. e confining pressure, axial pressure, axial deformation, and excess pore water pressure were recorded during the test. e above phenomenon is inconsistent with testing results of Shanghai mucky soft soil conducted by Fu et al. [23]. It is found from Figure 2(c) A5 load that the soil sample in this study directly transits from constant velocity creep to accelerated creep. is could be due to that the excess pore water pressure shortens the time of stable creep. Meanwhile, excess pore water pressure can obviously accelerate unloading creep of soft soil, and high excess pore water pressure is more likely to lead to unloading creep failure. As seen in an example at A4 load in Figures 2(a)-2(d), the creep deformation increases from 0.33% to 9.74% as excess pore water pressure increases from 0 to 60 kPa. In actual engineering, therefore, the excess pore water pressure should be reduced to prevent the unloading rheological failure of soft soil.

Viscoelastic Model of Unloading Creep.
e unloading creep stress-strain isochronous curves of soft soil under different excess pore water pressures are shown in Figure 3. It can be seen that, with the increase of excess pore water pressure from 0 to 60 kPa, the yield deviatoric stress (σ 1 -σ 3 ) s decreases from 34 kPa to 20 kPa. When the deviatoric stress σ 1 -σ 3 is lower than the yield deviatoric stress (σ 1 -σ 3 ) s , the creep of soft soil is mainly viscoelastic deformation. As an approximate straight-line segment shown in Figure 3, the slope of the straight-line segment is defined as viscoelastic modulus E ve (t), which shows that viscoelastic modulus E ve (t) gradually decreases with the increase of time. e relationship between viscoelastic modulus E ve (t) and time is shown in Figure 4.
According to the data distribution of viscoelastic modulus and time in Figure 4, the following function expressions of viscoelastic modulus and time can be obtained by adopting nonlinear fitting: where E ve (t) is the viscoelastic modulus of unloading creep, E 01 + E 02 is the initial viscoelastic modulus, B is the fitting coefficient, and t is the creep time. Fitted values of E 01 , E 02 , and B are shown in Table 2.
From Table 2, it can be seen that when the excess pore water pressure is 0 to 60 kPa, the above parameters all change to different degrees. Taking the fitting parameter of 0 kPa of excess pore water pressure as the standard, the variation range of initial viscoelastic modulus E 01 is 5% to 82%, the variation range of E 02 is 12% to 81%, and the variation range of fitting parameter B is 210% to 414%, which shows that excess pore water pressure had a significant influence on fitting parameter B. For simplicity, we mainly considered the influence of excess pore water pressure on parameter B, and the variation diagram of parameter B with excess pore water pressure is plotted, as shown in Figure 5. We found that the parameter B is almost linear with the excess pore water pressure. By fitting it with equation (2), the corresponding fitting parameters a and b are 0.30951 and 0.01574, respectively, and the square R 2 of the parameter fitting correlation coefficient is 0.94416; then, where a and b are fitting parameters, and u is the excess pore water pressure. By substituting equation (2) into equation (1), the following relationship between unloading creep viscoelastic modulus of soft soil and excess pore water pressure u and time t can be obtained: Strain (%)   Currently, the creep model of soft soil is mainly established by the element model wherein the elastic element, viscous element, and plastic element are used to simulate elastic deformation, viscous deformation, and plastic deformation, respectively. e combination of different elements can be used to simulate more complex deformation. e merchant model consists of the Kelvin model and series elastic elements. As shown in Figure 6(a), it presents better simulations in the viscoelastic deformation [24,25], but because all parameters in the merchant model are constant, the influence of excess pore water pressure was beyond the consideration of the model parameters.
erefore, assuming that the elastic modulus E of the spring in the merchant model conforms to the functional relation of equation (3), a modified merchant model can be obtained, as shown in Figure 6(b). By solving the differential equation, the following constitutive model of the modified merchant model can be obtained: where ε 1 is the viscoelastic strain, E 01 + E 02 is the initial viscoelastic modulus, E H is the instantaneous elastic modulus, η is the viscosity coefficient, σ is the stress level, a and b   Mathematical Problems in Engineering are the fitting parameters, u is the excess pore water pressure, and t is the creep time.

Viscoplastic Model of Unloading Creep.
According to the stress-strain isochronous curve of Figure 3, the viscoplastic deformation value can be obtained by subtracting the viscoelastic deformation value, as shown in Figure 7, which shows that when the stress level is higher than the yield deviatoric stress (σ 1 -σ 3 ) s , the unloading creep is mainly viscoplastic deformation. To describe the viscoplastic deformation as accurately as possible, a horizontal line is made along the stress [(σ 1 -σ 3 ) max -(σ 1 -σ 3 ) s ]/2, which intersects with the stress-strain curve at different times. e intersection point and the vertex at different times are connected into a straight line, and the slope of the straight line is defined as the viscoplastic modulus, as shown in Figure 7. e relationship between viscoplastic modulus and time at different times is plotted as shown in Figure 8. Similarly, the relationship between viscoplastic modulus and time can be expressed as a function of equation (1): where E vp (t) is the viscoplastic modulus of unloading creep, E P01 + E P02 is the initial viscoplastic modulus, B′ is the fitting coefficient that reflects the change rule of elastic modulus with time, and t is the creep time. e fitted values of E P01 , E P02 , and B′ are shown in Table 3. According to Table 3, various fitting parameters have changed to different degrees under the influence of excess pore water pressure. Taking the fitting parameters under 0 kPa as the standard, the change range of initial viscoplastic modulus E P01 is 33% to 61%, the change range of E P02 is 15% to 43%, and the change range of fitting parameter B′ is 30% to 135%. For simplicity, in this study, only the influence of excess pore water pressure on fitting parameter B′ was considered. By fitting it with equation (6), the corresponding fitting parameters a′ and b′ are 0.34235 and 0.00812, respectively, and the square R 2 of the parameter fitting correlation coefficient is 0.89955 (Figure 9); then, where a′ and b′ are the fitting parameters and u is the excess pore water pressure.

Mathematical Problems in Engineering
By substituting equation (6) into equation (5), the following relationship between unloading creep viscoplastic modulus and excess pore water pressure can be obtained: us, according to equation (7), the relationship between viscoplastic modulus and excess pore water pressure is established. As viscoplastic deformation occurs only when the stress exceeds the yield stress, in this study, the combined model of the merchant model and parallel plastic elements is used to describe viscoplastic deformation in this paper, as shown in Figure 10(a).
Since the parameter E vp in Figure 10(a) is constant, the influence of excess pore water pressure cannot be taken into consideration. erefore, by assuming that the parameter E vp in Figure 10(a) conforms to the functional relation of equation (7), the following constitutive model of the viscoplastic model can be obtained:  Deviatoric stress (σ Deviatoric stress (σ modulus, η P is the viscosity coefficient, and a′ and b′ are the fitting parameters; then, the following can be obtained by fitting test data: u is the excess pore water pressure; σ s is the yield stress; σ is the stress level; and t is the creep time. By connecting Figure 6(b) with Figure 10(b) in series, a viscoplastic model that describes unloading creep of soft soil can be obtained, as shown in Figure 11. e constitutive model of the viscoplastic model shown in Figure 11 is added by equations (4) and (8) Since equation (9) is the unloading creep constitutive model obtained under one-dimensional condition, while in the unloading creep test of soft soil, the soil is in a threedimensional stress state, equation (9) is deduced to obtain the unloading creep constitutive model under the threedimensional stress state.
Under a three-dimensional stress, the stress tensor σ ij on soft soil can be expressed by the following equation: where σ ij is the stress tensor; σ m is the spherical stress tensor; and S ij is the deviatoric stress tensor. It is generally recognized that the spherical stress tensor σ m changes only the volume, while the deviatoric stress tensor S ij changes only the shape. Accordingly, the strain tensor can be expressed by the following equation: where ε ij is the strain tensor; ε m is the spherical strain tensor; and e ij is the deviatoric strain tensor. For simplicity, on the assumption that the creep of soft soil is only generated by the deviatoric stress tensor [17], the following constitutive model of unloading creep under the triaxial stress state can be obtained: (1) When the deviatoric stress (σ 1 − σ 3 ) is less than the yield deviatoric stress (σ 1 − σ 3 ) s , (2) When the deviatoric stress (σ 1 − σ 3 ) is greater than the yield deviatoric stress (σ 1 − σ 3 ) s ,

Identification and Verification of Model Parameters.
By adopting the quasi-Newton algorithm (BFGS algorithm) of 1stOpt mathematical optimization analysis software, combined in the triaxial unloading creep test results of soft soil in Figure 2, the parameters of equations (12) and (13) Tables 4 and 5.
As the fitting curve shown in Figure 12, since the squares of the fitting correlation coefficients of the model parameters are all above 0.89, the fitting curve of the model is in good agreement with testing data, and the overall fitting effect is good. Zhou and Chen [14] established the merchant model that can be used to describe the unloading of soft soil under the one-dimensional state, but the model was not developed to describe the three-dimensional stress state, and the effect of excess pore water pressure was not considered. Further analysis shows that the fitting accuracy of the model is higher at low stress level than at high stress level. Especially, when soft soil creep changes from attenuation creep to constant velocity creep at high stress level, the curvature of the test data curve changes greatly. Limited by the equation, the simple multivariate model cannot accurately simulate this process; thus, the fitting value of the model is slightly lower than the test data value. In order to simulate this process more accurately, more Kelvin models in series are needed to be connected in the merchant model. e larger the curvature of the model equation, the more accurately it can   reflect the real characteristics of soft soil unloading creep. However, the series connection of more Kelvin models will undoubtedly increase the number of parameters to be determined in the models, which will inhibit their application in actual projects. erefore, considering the application of the models in projects, the merchant model (including a Kelvin model) is adopted in this paper.

Conclusion
In this study, the triaxial unloading creep tests of marine sedimentary soft soil in Shenzhen were conducted under different excess pore water pressures in undrained conditions, and an unloading creep constitutive model of soft soil is deduced. e following conclusions can be drawn: (1) Excess pore water pressure can significantly increase unloading creep of soft soil. Under the same partial stress, the unloading creep deformation of soft soil increases with the increase of excess pore water pressure. erefore, in the actual project, excess pore water pressure caused by construction shall be reduced as much as possible in order to prevent unloading rheological failure.
(2) Viscoelastic and viscoplastic moduli of soft soil exponentially reduced with creep time. By fitting the relationship between functional parameters and excess pore water pressure, the functional relationships between viscoelastic modulus, viscoplastic modulus, creep time, and excess pore water pressure are established, respectively. e modified merchant model can be used not only to describe the unloading creep under the condition of lower than yield deviatoric stress but also to represent the unloading creep under the condition of higher than yield deviatoric stress after parallel connection with plastic elements. (3) e unloading creep constitutive model of soft soil in the three-dimensional state is obtained through deduction. e accuracy of the model is related to the number of Kelvin models connected in series. e more Kelvin models are used, the higher the accuracy of the model is, but the parameters of the model also increase accordingly. erefore, the application of the unloading creep constitutive model of soft soil in engineering shall be considered. By comparing the test curve with the model curve, it shows that the two are in good agreement, which is consistent with the accuracy and applicability of the unloading creep constitutive model.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that there are no conflicts of interest.