# A Study of Coupled Creep Damaged Constitutive Model of Artificial Frozen Soil

Authors: Dongwei Li; Junhao Chen; Yan Zhou
Journal: Advances in Materials Science and Engineering
 (2018)

---

## Abstract

Artificial frozen soil is a kind of typical creep material, and the frozen clay under the unloading stress paths of high-confining pressure conforms to the improved the Zienkiewicz–Pande parabola-type yield criterion, and the Mohr–Coulomb yield function can describe the shear yield surface of artificial frozen clay under low-confining pressure. Based on the results of triaxial creep and shear tests for artificial frozen soil, the viscoplastic damage variable and evolution rule of artificial frozen clay were obtained by using the theory of viscoelastic-plastic mechanics and damage mechanics. An improved Zienkiewicz–Pande parabola-type yield criterion was used instead of a linear Newton body to obtain a coupled constitutive model of viscoelastic-plastic damage in the frozen soil under the unloading stress paths and to derive the coupling flexibility matrix for viscoelastic and viscoplastic damage. A finite element program of artificial frozen soil considering creep damage was written in the Visual Fortran 6.6A environment and embedded into the nonlinear finite element software ADINA as a user subroutine. The results of numerical simulation and laboratory testing were identical, with a maximum error of no more than 4.8%. This work shows that it is reasonable to describe the creep constitutive model of frozen soil with the viscoelastic-plastic-coupled constitutive model.


---

## Body

## 1. Introduction

Andersland and Douglas [1] successively proposed the creep theory of frozen soil. A series of constitutive models describing the creep yield process of frozen soil were presented, and a nonlinear yield criterion for frozen soil was proposed by Ladanyi [2, 3]. In the design of an artificial frozen wall during this period, attention was paid to the creep properties of frozen soil. It is found that the relationship between the uniaxial compressive strength of frozen soil and its temperature conforms to the power function [4], and the power function is applicable to both conditions of compression and tension [5–7].Frozen soil has special static-dynamic characteristics, so it is very important to study its basic mechanical properties [8]. Many scholars have conducted experiments, numerical simulations, and constitutive parameter inversions to study the creep damaged characteristics of frozen soil. However, there have been relatively few studies of the creep damaged characteristics of deep artificial frozen soil under unloading stress paths, especially the numerical calculation of the creep damaged model of artificial frozen clay under these paths [9, 10]. The damage variable of artificial frozen soil and its law of damage evolution can be derived based on the triaxial creep test, and its creep constitutive behavior can be studied using the finite element analysis.


## 2. Damage Constitutive Model of Artificial Frozen Soil Based on Viscoelastic-Plastic Theory

### 2.1. Creep Damaged Constitutive Model of Artificial Frozen Soil

Creep tests performed under unloading stress paths of artificial frozen soil show the following [11]:(1)
The artificial frozen soil has transient elastic strain under the unloading stress paths, and the instantaneous elastic strain is about 1.0%, less than 10% of the total creep deformation generally.(2)
The creep test of deep artificial frozen soil under the lower deviatoric stress shows two stages of unsteady creep and steady creep.(3)
The accelerating creep stage occurs when the stress level is above a certain critical value and will end after about 3 h.The Nishihara model can fully reflect the viscoelastic-plastic constitutive model in geotechnical creep material. Therefore, this model was selected according to the experimental results of artificial frozen soil under unloading stress paths and as the basis for improving the instantaneous elastic deformation and viscoelastic-plastic items, to apply the complex creep deformation law of deep artificial frozen soil. According to the rule of shear and creep deformation of artificial frozen soil under unloading stress path, the modified Zienkiewicz–Pande parabola-type yield surface was analyzed and used instead of the Saint-Venant body model. The damage variable was added for a new viscoelastic-plastic constitutive model suitable for deep artificial frozen soil under unloading stress paths of high-confining pressure, as shown in Figure1.Figure 1
Coupled creep damaged constitutive model of artificial frozen soil.For calculation, according to classical rheological finite element theory, the deep artificial frozen soil was regarded as a nonlinear continuum medium. The total strainεij includes a nonlinear viscoelastic strain component εijve and a nonlinear viscoplastic damage strain component εijvp. The nonlinear viscoelastic strain includes transient elastic strain and viscoelastic creep, and the nonlinear viscoplastic damage strain includes viscoplastic creep and damage creep. Therefore, the total strain of the coupled constitutive model under unloading stress paths can be expressed as follows:(1)εij=εijve+εijvp.The nonlinear viscoelastic strain componentεijve and the nonlinear viscoplastic damage strain component εijvp are solved below. For numerical calculation, the nonlinear ADINA finite element program was used to simulate the mechanical process of deep shaft excavation and unloading construction. In this program, the secondary development of the constitutive model requires the element flexibility matrix of each component. Therefore, the element flexibility matrix of the nonlinear viscoelastic components and the element flexibility matrix of nonlinear viscoplastic damage were derived below.


### 2.2. Element Flexibility Matrix of Nonlinear Viscoelastic Body

This section presents the viscoelastic-plastic constitutive model of artificial frozen soil under unloading stress paths of high-confining pressure, based on the test results of artificial frozen soil, in which the model of rheological elements includes the viscoelastic deformation property. The expression of the nonlinear viscoelastic strain componentεijve can be obtained directly according to the principle of rheological mechanics:(2)εijve=1E0+1E11−exp−E1tη1σij,where E1 is the viscoelastic modulus of deep artificial frozen soil (MPa); E0 is the transient modulus (MPa); η1 is the viscosity coefficient (MPa·h); t is the load time (h).The viscoelastic properties of deep artificial frozen soil can be treated as exhibiting instantaneous elasticity in the numerical calculation. By comparing with the generalized Hooke’s law, the viscoelastic flexibility matrixCve can be given directly by(3)Cve=1ET1−μT−μT000−μT1−μT000−μT−μT100021+μT00Symmetry21+μT021+μT,where ET is the equivalent nonlinear viscoelastic modulus and is the sum of the instantaneous elastic modulus and viscoelastic modulus. Based on elastic-plastic mechanics, the nonlinear viscoelastic modulus can be directly given by the following formula:(4)ET=11/E0T+1/E1T1−e−E1T1/η1Tt.


### 2.3. Element Flexibility Matrix of the Nonlinear Viscoplastic Body of Artificial Frozen Soil

The last item in the viscoelastic-plastic constitutive model of artificial frozen soil under unloading stress paths of high-confining pressure is the term for nonlinear viscoelastic-plastic damage, as shown in Figure1. The expression of the nonlinear viscoplastic strain component εijvp can be obtained directly from the rheomechanics principle:(5)εijvp=FF0tη2∂F∂σij,where F is the parabolic yield function of artificial frozen soil, as determined by (6). For frozen soil, F0=1; η2 is the viscosity coefficient (MPa·h); and t is load time (h):(6)F=sinφI1−3ccosφ+3sinφJ22,where φ is the internal friction angle of artificial frozen soil under unloading stress paths (°), determined by the triaxial test, and c is the cohesive force (MPa), determined by triaxial test.In order to improve the calculation efficiency and reliability, the implicit integration method was used to calculate the nonlinear viscoplastic flexibility matrix of a large geotechnical underground project. The basic idea was to determine the strain increment generated in a time increment intervalΔt=tn+1−tn by assuming the displacements δn and nodal stresses σn in the frozen wall at the moment of time tn:(7)Δεijn=Δtnε˙ijn+CvpΔσijn,where(8)Cvp=∂ε˙ijvp∂σij.According to the theory of geotechnical rheological mechanics, the viscoplastic strain rate can be expressed as follows:(9)∂ε˙ijvp∂σij=Δtη2∂F∂σij∂F∂σijT+F∂2F∂σij2,where(10)∂F∂σij=∂F∂I1∂I1∂σij+∂F∂J2∂J2∂σij+∂F∂J3∂J3∂σij,∂I1∂σij=P3σxσyσzτzyτxzτxyT,∂J2∂σij=12J2SxSySz2τzy2τxz2τxyT=Q2J2σxσyσzτzyτxzτxyT,∂J3∂σij=2J23PσxσyσzτzyτxzτxyT,C1=∂F∂I1=∂sinφ/3I1−ccosφ+sinφJ22∂I1=sinφ3,C2=∂F∂J2=∂sinφ/3I1−ccosφ+sinφJ22∂J22=36J23−sinφ2,C3=∂F∂J3=∂sinφ/3I1−ccosφ+sinφJ22∂J3=1833+sinφ,where P, R, and S are the coefficient matrixes of the following constitutive model:(11)P=100000100001000Symmery000000,Q=132−1−10002−10002000Symmetry600606,R=132−1−10002−10002000Symmetry120012012,S=132−1−10002−10002000Symmetry000000.Substitution of the above formulas into (9) yields(12)∂F∂σij=C1P3+C2Q2J2+C32J2P3.SinceP, R, and S coefficient matrixes are symmetrical matrixes, the transposed matrix of ∂F/∂σijT and coefficients in ∂F/∂σij are the same, as shown in the following formula:(13)∂F∂σijT=C1P3+C2Q2J2+C32J2P3.Therefore, the first item in (9) can be expressed as follows:(14)∂F∂σij∂F∂σijT=C12P9+C22R4J2+C324J22P9+C1C2S3J2+C1C34J2P9+C2C32J2S3=C12+4J22C32+4J2C1C39P+C224J2R+C1C2+2J2C2C33J2S.Similarly, sinceP, R, and S coefficient matrixes are symmetrical matrixes, the second item in (9) can be expressed as follows:(15)∂F∂σij2=C12+4J22C32+4J2C1C39P+C224J2R+C1C2+2J2C2C33J2S,Cvp=Δt1+Fη2C12+4J22C32+4J2C1C39P+C224J2R+C1C2+2J2C2C33J2S.The specific expressions of the nonlinear viscoplastic element flexibility matrix of deep artificial frozen soil under unloading stress paths are given below. It is assumed that(16)d1=Δt1+Fη2C12+4J22C32+4J2C1C39,(17)d2=Δt1+Fη2C224J2,(18)d3=Δt1+Fη2C1C2+2J2C2C33J2.Therefore, the specific nonlinear viscoplastic element flexibility matrix of the artificial frozen soil under unloading stress paths can be expressed in 3D space as follows:(19)Cvp=d1+2d2+d33−d2+d33−d2+d33000d1+2d2+d33−d2+d33000d1+2d2+d33000Symmetry4d2004d204d2.From (18), the nonlinear viscoplastic element flexibility matrix contains three parameters for frozen soil, that is, c,φ,andη2, which can be obtained through the triaxial creep test.


### 2.4. Element Flexibility Matrix of Artificial Frozen Soil Considering Creep Damage

For deep shaft freezing engineering, the strength and creep characteristics of artificial frozen soil are closely related to the internal damage of frozen wall. Due to the rheology of artificial frozen soil in the deep shaft freezing method, microcracks are produced inside the frozen wall and can develop continuously with time and gradually reduce the elastic modulus, Poisson’s ratio, and the strength of the frozen wall. The artificial viscoplastic creep exhibits obvious nonlinearity, which is reflected by the decrease of the viscoplastic viscosity coefficient. The decrease of the plastic viscosity coefficient with time is mainly caused by microelement damage to the artificial frozen soil. Therefore, the frozen soil microelement strength criterion (F∗) is assumed to conform to the Weibull-type probability distribution, and the probability density function can be expressed as follows:(20)PF∗=mF0F∗F0m−1×exp−F∗F0m,where F∗ is the Weibull-type probability distribution variable for microelement damage in artificial frozen soil and m and F0 are the Weibull probability distribution parameters, as determined by the creep test for artificial frozen soil.According to the theory of probability theory, the damage probability of artificial frozen soil, represented by the damage variableD can be expressed as follows:(21)D=∫0F∗Pydy=1−exp−F∗F0m.The first stress invariable and second deviatoric stress invariable that are expressed by effective stress and overstress definitions are as follows:(22)σm=σ^m1−D,J2=J^21−D.Therefore, the element flexibility matrix of viscoplastic damage coupling of artificial frozen soil under unloading stress paths is as follows:(23)Cvp=f1+2f2+f33−f2+f33−f2+f33000f1+2f2+f33−f2+f33000f1+2f2+f33000Symmetry4f2004f204f2,where(24)f1=Δt1+Fη2k12+4J22/1−D2k32+41−Dk1k39,f2=Δt1+Fη2k221−D4J2,f3=Δt1+Fη2k1k2+2J2/1−Dk2k33J2/1−D,k1=sinφ3,k2=36J2/1−D3−sinφ2,k3=1833+sinφ.When the artificial frozen soil conforms to the improved Zienkiewicz–Pande yield function within the main stress space, the viscoplastic creep begins. As the time of creep increases, microcracks in the frozen soil begin to expand and damage to the frozen soil appears. The damage gradually increases at a growing rate with the increase of deviatoric stress. The frozen soil structure is destroyed when many microcracks are formed, and the artificial frozen soil is completely damaged.


## 2.1. Creep Damaged Constitutive Model of Artificial Frozen Soil

Creep tests performed under unloading stress paths of artificial frozen soil show the following [11]:(1)
The artificial frozen soil has transient elastic strain under the unloading stress paths, and the instantaneous elastic strain is about 1.0%, less than 10% of the total creep deformation generally.(2)
The creep test of deep artificial frozen soil under the lower deviatoric stress shows two stages of unsteady creep and steady creep.(3)
The accelerating creep stage occurs when the stress level is above a certain critical value and will end after about 3 h.The Nishihara model can fully reflect the viscoelastic-plastic constitutive model in geotechnical creep material. Therefore, this model was selected according to the experimental results of artificial frozen soil under unloading stress paths and as the basis for improving the instantaneous elastic deformation and viscoelastic-plastic items, to apply the complex creep deformation law of deep artificial frozen soil. According to the rule of shear and creep deformation of artificial frozen soil under unloading stress path, the modified Zienkiewicz–Pande parabola-type yield surface was analyzed and used instead of the Saint-Venant body model. The damage variable was added for a new viscoelastic-plastic constitutive model suitable for deep artificial frozen soil under unloading stress paths of high-confining pressure, as shown in Figure1.Figure 1
Coupled creep damaged constitutive model of artificial frozen soil.For calculation, according to classical rheological finite element theory, the deep artificial frozen soil was regarded as a nonlinear continuum medium. The total strainεij includes a nonlinear viscoelastic strain component εijve and a nonlinear viscoplastic damage strain component εijvp. The nonlinear viscoelastic strain includes transient elastic strain and viscoelastic creep, and the nonlinear viscoplastic damage strain includes viscoplastic creep and damage creep. Therefore, the total strain of the coupled constitutive model under unloading stress paths can be expressed as follows:(1)εij=εijve+εijvp.The nonlinear viscoelastic strain componentεijve and the nonlinear viscoplastic damage strain component εijvp are solved below. For numerical calculation, the nonlinear ADINA finite element program was used to simulate the mechanical process of deep shaft excavation and unloading construction. In this program, the secondary development of the constitutive model requires the element flexibility matrix of each component. Therefore, the element flexibility matrix of the nonlinear viscoelastic components and the element flexibility matrix of nonlinear viscoplastic damage were derived below.


## 2.2. Element Flexibility Matrix of Nonlinear Viscoelastic Body

This section presents the viscoelastic-plastic constitutive model of artificial frozen soil under unloading stress paths of high-confining pressure, based on the test results of artificial frozen soil, in which the model of rheological elements includes the viscoelastic deformation property. The expression of the nonlinear viscoelastic strain componentεijve can be obtained directly according to the principle of rheological mechanics:(2)εijve=1E0+1E11−exp−E1tη1σij,where E1 is the viscoelastic modulus of deep artificial frozen soil (MPa); E0 is the transient modulus (MPa); η1 is the viscosity coefficient (MPa·h); t is the load time (h).The viscoelastic properties of deep artificial frozen soil can be treated as exhibiting instantaneous elasticity in the numerical calculation. By comparing with the generalized Hooke’s law, the viscoelastic flexibility matrixCve can be given directly by(3)Cve=1ET1−μT−μT000−μT1−μT000−μT−μT100021+μT00Symmetry21+μT021+μT,where ET is the equivalent nonlinear viscoelastic modulus and is the sum of the instantaneous elastic modulus and viscoelastic modulus. Based on elastic-plastic mechanics, the nonlinear viscoelastic modulus can be directly given by the following formula:(4)ET=11/E0T+1/E1T1−e−E1T1/η1Tt.


## 2.3. Element Flexibility Matrix of the Nonlinear Viscoplastic Body of Artificial Frozen Soil

The last item in the viscoelastic-plastic constitutive model of artificial frozen soil under unloading stress paths of high-confining pressure is the term for nonlinear viscoelastic-plastic damage, as shown in Figure1. The expression of the nonlinear viscoplastic strain component εijvp can be obtained directly from the rheomechanics principle:(5)εijvp=FF0tη2∂F∂σij,where F is the parabolic yield function of artificial frozen soil, as determined by (6). For frozen soil, F0=1; η2 is the viscosity coefficient (MPa·h); and t is load time (h):(6)F=sinφI1−3ccosφ+3sinφJ22,where φ is the internal friction angle of artificial frozen soil under unloading stress paths (°), determined by the triaxial test, and c is the cohesive force (MPa), determined by triaxial test.In order to improve the calculation efficiency and reliability, the implicit integration method was used to calculate the nonlinear viscoplastic flexibility matrix of a large geotechnical underground project. The basic idea was to determine the strain increment generated in a time increment intervalΔt=tn+1−tn by assuming the displacements δn and nodal stresses σn in the frozen wall at the moment of time tn:(7)Δεijn=Δtnε˙ijn+CvpΔσijn,where(8)Cvp=∂ε˙ijvp∂σij.According to the theory of geotechnical rheological mechanics, the viscoplastic strain rate can be expressed as follows:(9)∂ε˙ijvp∂σij=Δtη2∂F∂σij∂F∂σijT+F∂2F∂σij2,where(10)∂F∂σij=∂F∂I1∂I1∂σij+∂F∂J2∂J2∂σij+∂F∂J3∂J3∂σij,∂I1∂σij=P3σxσyσzτzyτxzτxyT,∂J2∂σij=12J2SxSySz2τzy2τxz2τxyT=Q2J2σxσyσzτzyτxzτxyT,∂J3∂σij=2J23PσxσyσzτzyτxzτxyT,C1=∂F∂I1=∂sinφ/3I1−ccosφ+sinφJ22∂I1=sinφ3,C2=∂F∂J2=∂sinφ/3I1−ccosφ+sinφJ22∂J22=36J23−sinφ2,C3=∂F∂J3=∂sinφ/3I1−ccosφ+sinφJ22∂J3=1833+sinφ,where P, R, and S are the coefficient matrixes of the following constitutive model:(11)P=100000100001000Symmery000000,Q=132−1−10002−10002000Symmetry600606,R=132−1−10002−10002000Symmetry120012012,S=132−1−10002−10002000Symmetry000000.Substitution of the above formulas into (9) yields(12)∂F∂σij=C1P3+C2Q2J2+C32J2P3.SinceP, R, and S coefficient matrixes are symmetrical matrixes, the transposed matrix of ∂F/∂σijT and coefficients in ∂F/∂σij are the same, as shown in the following formula:(13)∂F∂σijT=C1P3+C2Q2J2+C32J2P3.Therefore, the first item in (9) can be expressed as follows:(14)∂F∂σij∂F∂σijT=C12P9+C22R4J2+C324J22P9+C1C2S3J2+C1C34J2P9+C2C32J2S3=C12+4J22C32+4J2C1C39P+C224J2R+C1C2+2J2C2C33J2S.Similarly, sinceP, R, and S coefficient matrixes are symmetrical matrixes, the second item in (9) can be expressed as follows:(15)∂F∂σij2=C12+4J22C32+4J2C1C39P+C224J2R+C1C2+2J2C2C33J2S,Cvp=Δt1+Fη2C12+4J22C32+4J2C1C39P+C224J2R+C1C2+2J2C2C33J2S.The specific expressions of the nonlinear viscoplastic element flexibility matrix of deep artificial frozen soil under unloading stress paths are given below. It is assumed that(16)d1=Δt1+Fη2C12+4J22C32+4J2C1C39,(17)d2=Δt1+Fη2C224J2,(18)d3=Δt1+Fη2C1C2+2J2C2C33J2.Therefore, the specific nonlinear viscoplastic element flexibility matrix of the artificial frozen soil under unloading stress paths can be expressed in 3D space as follows:(19)Cvp=d1+2d2+d33−d2+d33−d2+d33000d1+2d2+d33−d2+d33000d1+2d2+d33000Symmetry4d2004d204d2.From (18), the nonlinear viscoplastic element flexibility matrix contains three parameters for frozen soil, that is, c,φ,andη2, which can be obtained through the triaxial creep test.


## 2.4. Element Flexibility Matrix of Artificial Frozen Soil Considering Creep Damage

For deep shaft freezing engineering, the strength and creep characteristics of artificial frozen soil are closely related to the internal damage of frozen wall. Due to the rheology of artificial frozen soil in the deep shaft freezing method, microcracks are produced inside the frozen wall and can develop continuously with time and gradually reduce the elastic modulus, Poisson’s ratio, and the strength of the frozen wall. The artificial viscoplastic creep exhibits obvious nonlinearity, which is reflected by the decrease of the viscoplastic viscosity coefficient. The decrease of the plastic viscosity coefficient with time is mainly caused by microelement damage to the artificial frozen soil. Therefore, the frozen soil microelement strength criterion (F∗) is assumed to conform to the Weibull-type probability distribution, and the probability density function can be expressed as follows:(20)PF∗=mF0F∗F0m−1×exp−F∗F0m,where F∗ is the Weibull-type probability distribution variable for microelement damage in artificial frozen soil and m and F0 are the Weibull probability distribution parameters, as determined by the creep test for artificial frozen soil.According to the theory of probability theory, the damage probability of artificial frozen soil, represented by the damage variableD can be expressed as follows:(21)D=∫0F∗Pydy=1−exp−F∗F0m.The first stress invariable and second deviatoric stress invariable that are expressed by effective stress and overstress definitions are as follows:(22)σm=σ^m1−D,J2=J^21−D.Therefore, the element flexibility matrix of viscoplastic damage coupling of artificial frozen soil under unloading stress paths is as follows:(23)Cvp=f1+2f2+f33−f2+f33−f2+f33000f1+2f2+f33−f2+f33000f1+2f2+f33000Symmetry4f2004f204f2,where(24)f1=Δt1+Fη2k12+4J22/1−D2k32+41−Dk1k39,f2=Δt1+Fη2k221−D4J2,f3=Δt1+Fη2k1k2+2J2/1−Dk2k33J2/1−D,k1=sinφ3,k2=36J2/1−D3−sinφ2,k3=1833+sinφ.When the artificial frozen soil conforms to the improved Zienkiewicz–Pande yield function within the main stress space, the viscoplastic creep begins. As the time of creep increases, microcracks in the frozen soil begin to expand and damage to the frozen soil appears. The damage gradually increases at a growing rate with the increase of deviatoric stress. The frozen soil structure is destroyed when many microcracks are formed, and the artificial frozen soil is completely damaged.


## 3. Finite Element Method to Solve the Problem of Artificial Frozen Soil for Viscoelastic Damage

The basic equation of the finite element method is the general equilibrium differential equation. At timet, this equation is expressed as follows:(25)∫ΩBtTσtdΩ+ft=0,where Bt is the creep damaged matrix of artificial frozen soil at time t. The artificial frozen soil is a kind of viscoelastic-plastic damage material, so the creep matrix is related to the stress state; σt is the total stress in artificial frozen soil at time t; and ft is the equivalent nodal load of artificial frozen soil at the time of t.Formula (25) is expressed in the incremental form as(26)∫ΩBtTΔσtdΩ+Δft=0,where σt is the total stress increment in artificial frozen soil during the period of Δt and ft is the equivalent nodal load change during the period of Δt.Calculation of theΔσt stress during the period of Δt with creep increment Δεt can be done using(27)Δσt=CsvpΔεt.The Euler integral method can then be applied:(28)Δε˙vpn=Δt1−Θεvpn+Θε˙vpn+1.The Taylor series expansion can then be done for viscoplastic creep damage of artificial frozen soil, to obtain(29)ε˙vpn+1=ε˙vpn+HΔσt.Substitution of (28) into (29) gives the arithmetic expression of the strain increment for viscoplastic creep damage as follows:(30)εvpn=Δtnε˙vpn+ΘΔtnHΔσt.Substitution of (30) into (27) yields:(31)Δσt=D−1+ΘΔtnH−1BΔd−ε˙vpΔt,H=∂ε˙vp∂σt.Therefore, in the incremental form, the equilibrium equation of the artificial frozen soil can be expressed as(32)∫ΩBtTD−1+ΘΔtnH−1BΔd−ε˙vpΔtdΩ+Δft=0.Therefore, the expression for stress increment of artificial frozen soil is as follows:(33)Δσt=D−1+ΘΔtnH−1BΔd−ε˙vpΔt.The expression for strain increment of viscoplastic damage in artificial frozen soil is as follows:(34)εvpn=BΔd−D−1Δσt.Formula (23) shows the flexibility matrix of finite element calculation that uses the viscoelastic-plastic damage coupled constitutive model of artificial frozen soil under unloading stress paths. Numerical calculation can be performed using the nonlinear finite element program ADINA, operating in the Visual Fortran 6.6A secondary development environment to generate the file of the dynamic link library [14]. The viscoelastic-plastic damage constitutive model of artificial frozen soil presented herein is written into an identifiable constitutive program in Fortran language by using the incremental finite element format in (34).


## 4. Determination and Verification of Parameters of the Constitutive Model of Artificial Frozen Soil

### 4.1. Determination of Parameters

The coupled constitutive model for creep damage presented herein contains eight parameters:c,φ,μ,η1,η2,E0,E1,andm.


#### 4.1.1. Determination of the Elastic ParameterscT,φT,Rf,andE0T

Based on the shear tests performed at different negative temperatures and axial pressures under the triaxial unloading stress paths, the curveJ2−σm was fitted on the coordinate plane of the first stress invariable and the second stress variable of deviatoric stress, to obtain the strength parameters cTandφT of artificial frozen soil at different temperatures.μ
T
and
E
0
T represent the instantaneous Poisson’s ratio and the modulus of elasticity of the artificial frozen soil, respectively. The method to determine the parameters is determined consistent with the determined parameters method of the Duncan–Chang model, a hyperbolic constitutive model considering the nonlinearity of instantaneous deformation [12, 13]. Another important elastic parameter of artificial frozen soil is Poisson’s ratio, which is determined by the ratio of axial strain and lateral strain in triaxial shear tests.


#### 4.1.2. Determination of the Viscoelastic Parametersη1andE1

Based on the uniaxial or triaxial shear creep test under unloading stress paths, the relationship between axial viscoelastic creep and time was obtained. The nonlinear fitting of the curve with the Gauss–Newton method according to the principle of least squares was used to obtain three viscoelastic parameters of artificial frozen soil:η1,η2,andE1.


#### 4.1.3. Determination of the Viscoplastic Damage Parametersη2andm

Based on the uniaxial or triaxial shear creep test under unloading stress path, the relationship between axial viscoplastic creep and time was obtained. The nonlinear fitting of the curve with the Gauss–Newton method according to the principle of least squares was used to obtain two viscoplastic parameters of artificial frozen soil:η2andm.The creep test for frozen soil samples in the air shaft at Huainan coal mine was performed under triaxial unloading conditions, and the following parameters are obtained for the constitutive model of artificial frozen soil herein:c,φ,μ,η1,η2,E0,andE1.


### 4.2. Verification of the Triaxial Test for the Constitutive Model of Artificial Frozen Soil

The creep deformation law of artificial frozen clay under test conditions was obtained through numerical simulation, with the axial stress-time relation curves, shown in Figure2.Figure 2
Comparison between measured and theoretical values from triaxial creep tests of artificial frozen soil at different freezing temperature and consolidation confining pressure: (a) −10°C, 4.0 MPa; (b) −10°C, 5.0 MPa; (c) −10°C, 6.0 MPa; (d) −10°C, 7.0 MPa.
(a)
(b)
(c)
(d)


## 4.1. Determination of Parameters

The coupled constitutive model for creep damage presented herein contains eight parameters:c,φ,μ,η1,η2,E0,E1,andm.


### 4.1.1. Determination of the Elastic ParameterscT,φT,Rf,andE0T

Based on the shear tests performed at different negative temperatures and axial pressures under the triaxial unloading stress paths, the curveJ2−σm was fitted on the coordinate plane of the first stress invariable and the second stress variable of deviatoric stress, to obtain the strength parameters cTandφT of artificial frozen soil at different temperatures.μ
T
and
E
0
T represent the instantaneous Poisson’s ratio and the modulus of elasticity of the artificial frozen soil, respectively. The method to determine the parameters is determined consistent with the determined parameters method of the Duncan–Chang model, a hyperbolic constitutive model considering the nonlinearity of instantaneous deformation [12, 13]. Another important elastic parameter of artificial frozen soil is Poisson’s ratio, which is determined by the ratio of axial strain and lateral strain in triaxial shear tests.


### 4.1.2. Determination of the Viscoelastic Parametersη1andE1

Based on the uniaxial or triaxial shear creep test under unloading stress paths, the relationship between axial viscoelastic creep and time was obtained. The nonlinear fitting of the curve with the Gauss–Newton method according to the principle of least squares was used to obtain three viscoelastic parameters of artificial frozen soil:η1,η2,andE1.


### 4.1.3. Determination of the Viscoplastic Damage Parametersη2andm

Based on the uniaxial or triaxial shear creep test under unloading stress path, the relationship between axial viscoplastic creep and time was obtained. The nonlinear fitting of the curve with the Gauss–Newton method according to the principle of least squares was used to obtain two viscoplastic parameters of artificial frozen soil:η2andm.The creep test for frozen soil samples in the air shaft at Huainan coal mine was performed under triaxial unloading conditions, and the following parameters are obtained for the constitutive model of artificial frozen soil herein:c,φ,μ,η1,η2,E0,andE1.


## 4.1.1. Determination of the Elastic ParameterscT,φT,Rf,andE0T

Based on the shear tests performed at different negative temperatures and axial pressures under the triaxial unloading stress paths, the curveJ2−σm was fitted on the coordinate plane of the first stress invariable and the second stress variable of deviatoric stress, to obtain the strength parameters cTandφT of artificial frozen soil at different temperatures.μ
T
and
E
0
T represent the instantaneous Poisson’s ratio and the modulus of elasticity of the artificial frozen soil, respectively. The method to determine the parameters is determined consistent with the determined parameters method of the Duncan–Chang model, a hyperbolic constitutive model considering the nonlinearity of instantaneous deformation [12, 13]. Another important elastic parameter of artificial frozen soil is Poisson’s ratio, which is determined by the ratio of axial strain and lateral strain in triaxial shear tests.


## 4.1.2. Determination of the Viscoelastic Parametersη1andE1

Based on the uniaxial or triaxial shear creep test under unloading stress paths, the relationship between axial viscoelastic creep and time was obtained. The nonlinear fitting of the curve with the Gauss–Newton method according to the principle of least squares was used to obtain three viscoelastic parameters of artificial frozen soil:η1,η2,andE1.


## 4.1.3. Determination of the Viscoplastic Damage Parametersη2andm

Based on the uniaxial or triaxial shear creep test under unloading stress path, the relationship between axial viscoplastic creep and time was obtained. The nonlinear fitting of the curve with the Gauss–Newton method according to the principle of least squares was used to obtain two viscoplastic parameters of artificial frozen soil:η2andm.The creep test for frozen soil samples in the air shaft at Huainan coal mine was performed under triaxial unloading conditions, and the following parameters are obtained for the constitutive model of artificial frozen soil herein:c,φ,μ,η1,η2,E0,andE1.


## 4.2. Verification of the Triaxial Test for the Constitutive Model of Artificial Frozen Soil

The creep deformation law of artificial frozen clay under test conditions was obtained through numerical simulation, with the axial stress-time relation curves, shown in Figure2.Figure 2
Comparison between measured and theoretical values from triaxial creep tests of artificial frozen soil at different freezing temperature and consolidation confining pressure: (a) −10°C, 4.0 MPa; (b) −10°C, 5.0 MPa; (c) −10°C, 6.0 MPa; (d) −10°C, 7.0 MPa.
(a)
(b)
(c)
(d)


## 5. Conclusions

The creep damaged constitutive program for viscoelastic-plastic damage of artificial frozen soil under unloading stress paths was developed into a large nonlinear finite element program using ADINA. The following conclusions can be obtained:(1)
The triaxial shear and creep deformation laws of deep artificial frozen soil were obtained by simulating stress paths after freezing and unloading to establish the constitutive equation of artificial frozen soil under complex stress paths.(2)
The value from the triaxial creep test of deep artificial frozen soil under unloading stress paths was consistent with the simulated value, with a maximum error of 4.8%. A constitutive model can reflect the primary creep, steady creep, or the accelerating creep of artificial frozen soil was proposed. These results demonstrate that it is reasonable to describe the nonlinear creep of artificial frozen soil with the viscoelastic-plastic damage coupled model using the improved Zienkiewicz–Pande parabola yield criterion.(3)
A combination of basic rheological elements and yield surface models was used to describe the creep damaged deformation law of artificial frozen soil. This approach allowed development of a viscoelastic-plastic damage coupled constitutive model of deep artificial frozen soil under unloading stress paths. The calculation parameters can be obtained through a simple test, and the model parameters are few and easy to be determined. Overall, this model can be easily and widely applied during artificial freezing engineering.(4)
The constitutive model provides an effective research approach to study the mechanical behavior of field frozen soil structure and presents theoretical guiding significance for the long-term stability analysis and the prediction of artificial frozen soil structure.



---

*Source: 7458696-2018-06-25.xml*