# The Constitutive Model of Rockfill Based on Property-Dependent Plastic Potential Theory for Geomaterials

Authors: Xuefeng Li; Ruijie Li; Junhui Zhang; Qi Wang
Journal: Advances in Civil Engineering
 (2020)

---

## Abstract

To better control the strength and deformation of the roadbed, a constitutive model of rockfill was established based on property-dependent plastic potential theory for geomaterials. The effect of the particle gradation on the anisotropy was described in the model. According to the effect of the particle grading and crushing on the fractal dimension, the fractal theory and fabric tensor were introduced to establish the yield and failure criteria of the rockfill. By combining the property-dependent concepts of the materials and the results of the rockfill strength test, a critical state line considering the microstructure, fractal dimension, particle breakage, and stress state of the rockfill was established. The dilatancy equation was derived based on the novel potential theory and the hardening criterion affected by the critical state was established. A constitutive model of the rockfill in the general stress space was established under the framework of the novel potential theory. The 3D strength and its intensity change in theπ plane were simulated through the drainage strength test results, which verified the description of the critical state under various stress paths. By simulating the stress-strain relationship, the validity and rationality of the model were verified.


---

## Body

## 1. Introduction

Rockfill is a common engineering material that is largely used in railway, transportation, and highway roadbeds construction projects due to its high filling density, strong permeability, good compacting performance, and high strength. Studying more about the stress, strain, and deformation characteristics of rockfill materials under various loads in detail and establishing a reasonable and effective constitutive model were of great significance to the theoretical analysis of rockfill engineering application and the numerical calculation and analysis of transportation geotechnics [1–3].The constitutive relationship is not only the core of modern soil mechanics, but also an important entry point for studying the mechanical properties of rockfill materials. Compared with sand, the particle size and the void of rockfill material are larger. Under complex stress conditions, the spatial fabric and particle size distribution of rockfill are usually more prone to slipping and breaking. The complex macroscopic mechanical properties often show a tendency for nonlinear forms of strength and volume change, which affect the stability of a rockfill fabric [4]. This has brought some difficulties in studying the constitutive relationship of rockfill. Some researchers have done a large amount of work and achieved positive research results such as the strain-softening and dilatancy of rockfill [5], material constants and particle characteristics [6], load and creep [7, 8], particle breakage and the relationship between intermediate principal stress coefficient b value and particle breakage [9], loading path [10], wetting deformation [11], the coupling effect of mean effective stress p and deviatoric stress q on deformation [12], and so on.In the context of generalized plasticity theory, Wang et al. [13] and Liu et al. [14] established a constitutive model of rockfill to describe the particle breakage considering critical state and dilatancy. Xiao and Liu [15] proposed a critical state line for particle crushing called the breakage critical state line (BCSL) and established an elastoplastic model for rockfill materials considering the state dependence and particle crushing based on critical state soil mechanics. Liu et al. [16] and Fang et al. [17] proposed an elastoplastic constitutive model and a state-dependent 3D multimechanism boundary surface model, respectively, by introducing state-related parameters, but the models still had certain limitations in the adaptation and expansion of complex stress paths. Liu and Chen [18] established an exponential-parabolic nonlinear elastic constitutive model reflecting the porosity and density of rockfill material based on the improved Nanshui double yield surface model. Besides, Brito et al. [19] established a new model specifically describing soil-rockfill mixtures (SRM).In summary, many researchers have established a constitutive model that could describe the hardness, dilatancy, loading, and unloading of rockfill materials suitably under complex stress paths based on elastic-plastic theory, subloading surface theory, generalized plasticity theory, boundary surface theory, and others. These constitutive models have achieved positive and effective results in describing the basic strength and deformation characteristics of materials, such as friction and dilatancy. However, the existing models still have the shortcomings of narrow application range, many model parameters, and complicated forms. Meanwhile, these models cannot accurately describe the isotropy and anisotropy of rockfill materials. Furthermore, some models also cannot reflect the plastic deformation that is caused by the rotation of the principal stress. Therefore, determining how to correctly describe the deformation characteristics of rockfill materials, quoting the small amount of necessary mechanical parameters, and unifying the anisotropy and isotropic characteristics of the material characteristics have become the research difficulties. This study tries to adopt the potential theory to solve the above problems. Compared with the traditional constitutive models, the mathematical principle is clear and is not based on a plastic postulate. It connected physical assumptions and mathematical foundations established by the constitutive model. Thus, it forms a complete set for a theoretical system. The main principle of potential theory involves taking the principal stress and principal strain of the main space or its increment as a mathematical vector and using the idea of vector fitting to fit the known vector [20, 21].Under the framework of the theory of continuum mechanics, Li et al. [22] linked the strain distribution rule of geotechnical materials with its fabric properties and proposed the property-dependent plastic potential theory. The proposed theory was verified and applied to sand constitutive models with a good description result. Compared to the noncoaxial plastic theory proposed by Gutierrez et al., the theory of Li et al. has a better description of the noncoaxial characteristics. However, it is necessary to consider the differences in the fabrics and the mechanical properties between sand and rockfill, so it is not easy to directly apply the potential theory of sand to rockfill. The particle fractal and fabric of rockfill material can be described as the microscopic characteristics of the material, and the size effect builds a bridge between the sand research and the rockfill material. The key idea of the plasticity theory related to material properties in this research was that, according to the material’s meso-structural properties, it was assumed that the material properties could be described by a fabric tensor. Then, the strain distribution of the material could be affected by the material properties. Based on this, in this research, the property-dependent plastic potential theory will be extended to describe anisotropy characteristics of rockfill, and a constitutive model is proposed for rockfill with general stress space. Based on the fractal dimension, the novel state variables will be introduced. Additionally, in this research, the yield failure criterion of rockfill will be proposed, the dilatancy equation of rockfill will be derived by proposed potential theory, hardening criterion considering the influence of critical state will be established, and the rationality of its 3D strength description will be verified through large-scale triaxial test results and the determination of effectiveness.


## 2. Fractal Theory

Rockfill has the particular properties of coarse particles with large sizes and its easily broken characteristics. There is a significant relationship among the changing gradation along with the stress, dilatancy, and strength-deformation characteristics in rockfill materials [23, 24]. Therefore, it is important to study the mechanical properties of rockfill materials by quantitatively describing its gradation of rockfill materials which approaches the true properties and status of raw materials. Zhang and Zhang [25] found that the particle yield and compression state were closely related to the particle-size distribution and fractal curve characteristics of particle breakage through a one-dimensional compression test. Besides, the study found that the fractal dimension D could be used to describe the rockfill particle-size distribution. This was especially true for the two indicators Cu and Cc used to characterize the particle-size distribution of the large particle rockfill material of the earth-rock dam and highway roadbed. The indicators Cu and Cc could accurately determine the size distribution and the optimal content of fine particles. This means that the fractal theory had a strong theoretical and practical significance for the study of constitutive relationships [26].The fractal relationship for the particle size of the rockfill and the particle mass and particle volume was described by the fractal model. It could be expressed as(1)lgMr<dmaxMT=3−DlgdmaxλV,where M is the mass of the certain particles, MT is the total mass of the particles, and dmax and λV stand for the maximum particle size and the characteristic particle size, respectively. D is the fractal dimension of the rockfill particle size, which was defined as follows:(2)D=3−lnMr<di/MTlndi/dmax.This shows that as the nonuniformity coefficient of the rockfill gradation becomes larger, the number of particlesNr that are less than r from equation (2) increases. Thus, the dimension D was larger. Formula (2) better reflected the grading of rockfill, so it could be used to represent the grading of rockfill.Studies have shown that the crushing of rockfill particles becomes increasingly severe with an increase in confining pressure [27]. The reason for this is that the particles inside the rockfill are easily broken with the multiple effects occurring under higher confining pressures. This usually caused the quality of particles with larger particle sizes to decrease. However, the relative density of the rockfill material had little effect on the particle crushing. The specific behavior was as follows. The finer the initial rockfill gradation was (i.e., the larger the initial fractal dimension was), the smaller the crushing degree was. The experiment confirmed that the fractal dimension D had a good linear relationship with σ3/Pa, which could be obtained by binary linear regression analysis:(3)D=l+κσ3Pa+βD0,where l, κ, and β are model parameters that could be calculated from the experimental data.


## 3. Constitutive Model Construction

### 3.1. Elastic Description

The elastic shear modulus of rockfill that describes its pressure hardness characteristics based on the empirical equation (4) can be found in the work by Richart et al.:(4)G=G02.97−e21+ep⋅pa0.5,where e is the current void ratio, Pa is the atmospheric pressure (Pa = 101 kPa), and G0 is the material constant. Therefore, the bulk modulus of elasticity can be expressed as(5)K=G021+ν31−2ν⋅pa2.97−e21+eppa0.5,where v is Poisson’s ratio. Equations (4) and (5) indicate that the strength of a material increases with the enhancement of the stress level. This conforms to the pressure hardness characteristic for which the strength of the rockfill increases with the confining pressure.


### 3.2. Yield Surface

Based on the work of the failure criterion for anisotropic sand and rockfill by Li et al. [28, 29], the novel anisotropic state variable was introduced, and the yield criterion in this research could be written as(6)f=q−1+ςAMgθσp=0,where M is the internal variable of hardening and ς is the weight coefficient, which was usually equal to 0.2. p and q are the mean stress and the generalized deviatoric stress, respectively. p=σkk/3, q=3J2, J2=sijsij/2, and sij=σij−δijp. δij is the Kronecker tensor (when i≠j, δij = 0 and when i=j, δij = 1). θσ=sin−133J3/2J23/3, J3=sijsjkski/3, and θσ is the stress Lode angle. A is the anisotropic state variable. gθσ is a function of θσ that reflects the influence of the intermediate principal stress. The elliptic function of gθσ in equation (6) was proposed by William and Warnke [30]:(7)gθσ=21−c2cosπ/6+θσ+2c−141−c2cos2π/6+θσ+c5c−441−c2cos2π/6+θσ+2c−12,where c is the ratio of the triaxial tension strength to the triaxial compression strength. To ensure the convexity of the yield surface, the law of compression and elongation was used. The yield surface had to satisfy the following conditions: when θσ was equal to π/6 or −π/6, dgθσ/dθσ was equal to 0. When g−π/6 was equal to c, gπ/6 was equal to 1. Mf is the peak failure stress ratio, which stands for equation (8) with the peak internal friction angle φf:(8)Mf=6sinφf3−sinφf.It was found that the strength envelope of the rockfill materials exhibited a nonlinear trend. With the increase of the confining pressure, the slope of the envelope gradually decreased and the longitudinal intercept increased. Therefore, the Mohr–Coulomb yield surface had some limitations in describing the strength characteristics of the rockfill. The test showed that the relationship between the ratio of the deviatoric stress and the atmospheric pressureτ/Pa and the ratio of normal stress and the atmospheric pressure σ/Pa followed the power function in(9)τpa=AcσPaBc.The peak internal friction angle decreased with the increase of the confining pressure, which could be expressed as(10)φf=arctanAcσPaBc−1,where Ac and Bc are the strength parameters. The amplitude parameters of the x1−x2 and x2−x3 planes could be obtained similarly as follows:(11)Ac=fD0,e,Bc=fe.Ac was proportional to the fractal dimension D0 and inversely proportional to the void ratio e, while the fractal dimension D0 had little effect on Bc.


### 3.3. Orthotropic Fabric Tensor

The anisotropic state variable based on the research of the quantitative description of mesoscopic fabric by Li et al. [31] in equation (12) was expressed as(12)A=σ¯ijσmF¯ij−σ¯ijσmF¯ij0=σ¯ijσmF¯ij−A0,where σ¯ij is stress tensor. In particular, σ¯ij is principal stress tensor and σm is the partial tensor of the stress tensor when the fabric tensor is expressed in terms of the main fabric space. In equation (12), the equation σ¯ij/σmF¯ij=A¯ represents the secondary anisotropic state variable induced by an external load. Thus, A0 is the native anisotropic state variable. In conclusion, according to the work by Li et al. [31], the mesoscopic fabric tensor was given by(13)F¯ij=1+a1+a2+a1a23+a1+a2−a1a21+a1−a2−a1a23+a1+a2−a1a21−a1+a2−a1a23+a1+a2−a1a2,where a1 and a2 are orthotropic amplitude parameters that could be expressed as(14)a1=12N∑k=12Ncos2θ1k−sin2θ1kcos2αk2+∑k=12Nsin2θ1kcosαk2,a2=12N∑k=12Ncos2θ1k−sin2θ1ksin2αk2+∑k=12Nsin2θ1ksinαk2,a3=12N∑k=12Nsin2θ1kcos2αk2+∑k=12Nsin2θ1ksin2αk2,where a1, a2, and a3 are the Fx1−Fx3, Fx1−Fx2, and Fx2−Fx3 orthotropic amplitude parameters on three orthogonal surfaces. According to the microscopic test, all three parameters could be determined. The range of the values allowed the amplitude parameters to be in the range of 0 to 1. It reflected the anisotropy degree of the materials on the three deposition surfaces. θk is the angle between the long axis orientation of the kth particle and the corresponding coordinate axis. N is the total number of measured rockfill particles and αk is the kth angle between the projection of the long axis of the particle on the horizontal plane and the xk axis.The definition of the 3D fabric equation (13) required that Fij be equal to 1. There were only two independent variations for the determined fabric amplitude parameters on the a1, a2, and a3 surfaces. The expressions of the three orthogonal fabric tensors (F1, F2, and F3) obtained by (14) were functions of a1 and a2, a1 and a3, and a2 and a3, respectively. Therefore, only the information for the microscopic arrangement of the rockfill on any two of the three orthogonal surfaces of the particles was obtained through experiments, so the orthotropic fabric tensors of the rockfill material could be determined completely.


## 3.1. Elastic Description

The elastic shear modulus of rockfill that describes its pressure hardness characteristics based on the empirical equation (4) can be found in the work by Richart et al.:(4)G=G02.97−e21+ep⋅pa0.5,where e is the current void ratio, Pa is the atmospheric pressure (Pa = 101 kPa), and G0 is the material constant. Therefore, the bulk modulus of elasticity can be expressed as(5)K=G021+ν31−2ν⋅pa2.97−e21+eppa0.5,where v is Poisson’s ratio. Equations (4) and (5) indicate that the strength of a material increases with the enhancement of the stress level. This conforms to the pressure hardness characteristic for which the strength of the rockfill increases with the confining pressure.


## 3.2. Yield Surface

Based on the work of the failure criterion for anisotropic sand and rockfill by Li et al. [28, 29], the novel anisotropic state variable was introduced, and the yield criterion in this research could be written as(6)f=q−1+ςAMgθσp=0,where M is the internal variable of hardening and ς is the weight coefficient, which was usually equal to 0.2. p and q are the mean stress and the generalized deviatoric stress, respectively. p=σkk/3, q=3J2, J2=sijsij/2, and sij=σij−δijp. δij is the Kronecker tensor (when i≠j, δij = 0 and when i=j, δij = 1). θσ=sin−133J3/2J23/3, J3=sijsjkski/3, and θσ is the stress Lode angle. A is the anisotropic state variable. gθσ is a function of θσ that reflects the influence of the intermediate principal stress. The elliptic function of gθσ in equation (6) was proposed by William and Warnke [30]:(7)gθσ=21−c2cosπ/6+θσ+2c−141−c2cos2π/6+θσ+c5c−441−c2cos2π/6+θσ+2c−12,where c is the ratio of the triaxial tension strength to the triaxial compression strength. To ensure the convexity of the yield surface, the law of compression and elongation was used. The yield surface had to satisfy the following conditions: when θσ was equal to π/6 or −π/6, dgθσ/dθσ was equal to 0. When g−π/6 was equal to c, gπ/6 was equal to 1. Mf is the peak failure stress ratio, which stands for equation (8) with the peak internal friction angle φf:(8)Mf=6sinφf3−sinφf.It was found that the strength envelope of the rockfill materials exhibited a nonlinear trend. With the increase of the confining pressure, the slope of the envelope gradually decreased and the longitudinal intercept increased. Therefore, the Mohr–Coulomb yield surface had some limitations in describing the strength characteristics of the rockfill. The test showed that the relationship between the ratio of the deviatoric stress and the atmospheric pressureτ/Pa and the ratio of normal stress and the atmospheric pressure σ/Pa followed the power function in(9)τpa=AcσPaBc.The peak internal friction angle decreased with the increase of the confining pressure, which could be expressed as(10)φf=arctanAcσPaBc−1,where Ac and Bc are the strength parameters. The amplitude parameters of the x1−x2 and x2−x3 planes could be obtained similarly as follows:(11)Ac=fD0,e,Bc=fe.Ac was proportional to the fractal dimension D0 and inversely proportional to the void ratio e, while the fractal dimension D0 had little effect on Bc.


## 3.3. Orthotropic Fabric Tensor

The anisotropic state variable based on the research of the quantitative description of mesoscopic fabric by Li et al. [31] in equation (12) was expressed as(12)A=σ¯ijσmF¯ij−σ¯ijσmF¯ij0=σ¯ijσmF¯ij−A0,where σ¯ij is stress tensor. In particular, σ¯ij is principal stress tensor and σm is the partial tensor of the stress tensor when the fabric tensor is expressed in terms of the main fabric space. In equation (12), the equation σ¯ij/σmF¯ij=A¯ represents the secondary anisotropic state variable induced by an external load. Thus, A0 is the native anisotropic state variable. In conclusion, according to the work by Li et al. [31], the mesoscopic fabric tensor was given by(13)F¯ij=1+a1+a2+a1a23+a1+a2−a1a21+a1−a2−a1a23+a1+a2−a1a21−a1+a2−a1a23+a1+a2−a1a2,where a1 and a2 are orthotropic amplitude parameters that could be expressed as(14)a1=12N∑k=12Ncos2θ1k−sin2θ1kcos2αk2+∑k=12Nsin2θ1kcosαk2,a2=12N∑k=12Ncos2θ1k−sin2θ1ksin2αk2+∑k=12Nsin2θ1ksinαk2,a3=12N∑k=12Nsin2θ1kcos2αk2+∑k=12Nsin2θ1ksin2αk2,where a1, a2, and a3 are the Fx1−Fx3, Fx1−Fx2, and Fx2−Fx3 orthotropic amplitude parameters on three orthogonal surfaces. According to the microscopic test, all three parameters could be determined. The range of the values allowed the amplitude parameters to be in the range of 0 to 1. It reflected the anisotropy degree of the materials on the three deposition surfaces. θk is the angle between the long axis orientation of the kth particle and the corresponding coordinate axis. N is the total number of measured rockfill particles and αk is the kth angle between the projection of the long axis of the particle on the horizontal plane and the xk axis.The definition of the 3D fabric equation (13) required that Fij be equal to 1. There were only two independent variations for the determined fabric amplitude parameters on the a1, a2, and a3 surfaces. The expressions of the three orthogonal fabric tensors (F1, F2, and F3) obtained by (14) were functions of a1 and a2, a1 and a3, and a2 and a3, respectively. Therefore, only the information for the microscopic arrangement of the rockfill on any two of the three orthogonal surfaces of the particles was obtained through experiments, so the orthotropic fabric tensors of the rockfill material could be determined completely.


## 4. Model Framework

### 4.1. Property-Dependent Plastic Potential Theory

To reflect the anisotropy of the rockfill fabric and the uniqueness of the rockfill critical state, the consideration of rockfill mesoscopic fabric material properties related to the plastic potential theory was put forward based on material status-related dilatancy theory. It was assumed that the material properties could be represented by the fabric tensor and that the rockfill material properties would influence the strain distribution. For an incremental property of the geotechnical materials, the general expression of the property-dependent plastic potential related to rockfill materials can be written as(15)dεijp=dλ∂g∂σilF¯lj,where g is the plastic potential function, dλ is the increment of the plastic factor, and F¯ij is the tensor of the fabric tensor Fij after the transformation of equation (15). According to the transformation characteristics of equation tensors, equation (15) was further simplified to(16)Fij=pFδij+s˜ij=pFF¯ij,where pF and s˜ij are the spherical tensor and the partial tensor of the fabric tensor, respectively.Equations (15) and (16) could be rearranged as another expression of the material property potential theory:(17)dεijp=dλ∂g∂σij+dλ∂g∂σils˜lj,where dλ∂g/∂σij is traditional plastic potential theory. When the material was isotropic, dλ∂g/∂σijs˜lj was equal to 0, and the plastic factor dλ was a nonnegative proportional coefficient in equation (17).


### 4.2. Hardening Law

The hardening law of internal variables was adopted in the incremental hyperbolic form proposed by Li and Dafalias [32]. It was given by (18)dεsp=MhsGMp−MdM,where G is the shear modulus. Mp is the peak stress ratio in a conventional triaxial compression test, and hs is the material constant.It must be noted that the peak stress ratioMp had to not be a constant, but rather a quantity related to the state parameters ψ. In this research, the equation for rockfill was adopted, Mp=Mcsexp−kpψ, as suggested by Li and Dafalias, and kp was the constant of the model. The state parameter ψ was expressed as(19)ψ=e−ec,where ec is the critical void ratio which could be determined from the following equation:(20)ec=eΓ−λcppaξ,The results showed that there was a good linear relationship between the critical void ratio intercepteГ and the initial void ratio e0, which satisfied the relationship(21)eΓ=c′+b′e0.Moreover, there was a linear strip relationship between the critical void ratio intercepteГ and the fractal dimension D, so the critical void ratio intercept eГ was expressed as(22)eΓ=c′−a′D+b′e0,where a′, b′, and c′ are the model parameters that could be calculated from the experimental data.For the value ofhs, Li and Dafalias [32] considered that it was related to the void ratio, taken as hs = h1−h2e, where h1 and h2 were the constants of the model. In addition, the model parameters were also affected by anisotropy, and its expression could be rearranged as(23)hs=h1−h2e1+khA,where kh is the model parameter and A is the defined anisotropic state variable in equation (12), when the material is isotropic, A≡0.To consider the effect of the anisotropy, the expression of the peak stress ratioMp could be changed into the function of an anisotropic state variable:(24)Mp=Mcsgθσe−kpψA,where ψA is the state parameter related to the initial anisotropy described in equation (24) and kp is the model parameter. The expression of the peak stress ratio showed that the peak stress ratio was not a constant but rather an obvious function of the anisotropic state variable.


### 4.3. Dilatancy Equation

The dilatancy equation could be deduced from property-dependent plastic potential theory. Based on the basic theory of rockfill mechanics in a critical state, when the rockfill reached the critical state, the volume change was equal to zero and the stress ratio reached the critical state pressure ratioηc. Lü et al. found that stress dilatation is a function related to the state of geomaterials [33]. At this point, the deformation of the rockfill only had shear strain, so its energy relationship satisfied(25)dWp=pdεvp+qdγpcosθσ−θdε=ηcpdγp,where ηc is the stress ratio under the critical state. According to the energy relationship of the critical state given in equation (25), the dilatancy equation based on the potential theory related to material properties could be obtained as follows:(26)dεvpdγp=ηc−cosθσ−θdεqp.The dilatancy equation is usually determined by conventional triaxial tests, and it can be expressed as follows:(27)d=d0ηc−cosθσ−θdεqp,where d0 is the dilatation coefficient, which was determined by the experiment, θσ is the stress Lode angle, and θdε is the strain increment Lode angle in principal stress space. When the stress ratio was greater than the critical state (d < 0) for the dilatancy equation, rockfill dilatancy occurred. When the stress ratio was less than the critical stress (d > 0), the rockfill contracted. When the stress ratio was equal to the critical state (d = 0), the rockfill was in the state of phase, as shown in Figure 1.Figure 1
Relationship between the critical state stress ratio and dilatation.In the principal stress space, the strain incrementθdε of the Lode angle was defined as(28)θdε=arctan132dε2−dε1−dε3dε1−dε3,where the main strain was determined by property-dependent plastic potential theory in equation (25). Assuming q/p=η, equation (29) could be rearranged as(29)dq=ηdp+pdη,By substituting equation (29) into equation (27), we obtained the following:(30)dpp=−dηd+η,The plastic potential function was obtained by integrating equation (30):(31)g=q−pd0ηcd0cψ−11−pp0d0cψ−1=0,where cψ is the state parameter, the value could be obtained from the equation:cψ=cosθσ−θdε under the critical state of rockfill.Considering the influence of anisotropy on the critical state stress ratioηc, the expression of ηc in 3D space was defined as follows:(32)ηc=1+ςAMcsgθσekdψA,where kd is the model parameter.


### 4.4. Incremental Stress-Strain Relationship

Based on the expression method of the strain increment, the total strain could be expressed as(33)dεij=dεije+dεijp,where dεij is the total strain increment and dεije is the elastic strain increment, the expression for which is shown in equation (34). According to the potential theory related to the material properties, dεijp is the plastic strain increment, which is expressed in equation (17), as follows:(34)dεije=Cijkledσkl.According to equations (17), (33), and (34), the total strain increment was derived as dσij=Dijkledεkl−dλ∂g/∂σisFsj:(35)dεij=Cijkledσkl+dλ∂g∂σikFkj.Equation (35) was simplified as follows:(36)dσij=Dijkledεkl−dλ∂g∂σisFsj.Equation (40) could be rearranged as(41)dλ=∂f/∂σijDijkledεklAp+∂f/∂σmnDmnpqe∂g/∂σsuFut.where Dijkle is the elastic stiffness matrix. dλ is also related to the hardening parameter Ap, so it could be rearranged as(37)dλ=1Ap∂f∂σijdσij.According to the consistency equation df equal to 0 on the yield surface, the equation became(38)df=∂f∂σijdσij+∂f∂MdM=0,where f is the yield function and M is the hardening parameter.Based on Equations (37) and (38), this led to the expression below:(39)dλAp=−∂f∂M∂M∂εijpdεijp.After combining equations (36)–(39), the expression was given by(40)∂f∂σijDijkledεkl−∂f∂σmnDmnpqedλ∂g∂σsupFF¯ut−Apdλ=0.Equations (37) and (38) were combined as follows:(42)Ap=−∂f∂MdMdλ=−∂f∂M∂M∂εijp∂g∂σikFkj.Because the hardening law took the plastic shear strain as the hardening parameter in the form of an incremental hyperbola, it could be obtained according to the flow law:(43)Ap=−∂f∂MdMdλ=−∂f∂M∂M∂εqp∂g∂qFkj2+1q∂g∂θσFkj2.The hardening law equation (18) was substituted into equation (43) to obtain the hardening function as(44)Ap=gθσphsGMp−M2pMp2∂g∂qFkj2+1q∂g∂θσFkj2.


### 4.5. Model Parameter Determination

There were four groups and seventeen model parameters in this research. The elastic parameters and the fractal parameters could be measured using the basic test data of the rockfill materials. The critical state parameters were tested with a triaxial tensile test. The state-related parameters were determined with the drained triaxial test. All model parameters are shown in Table1. The constitutive model parameters were divided into four groups: elastic parameter, fractal parameter, critical parameter, and state parameter. The specific methods for determining these parameters were as follows.Table 1
Parameters of the state-dependent constitutive model for rockfill.
Elastic parameterFractal parameterCritical parameterState parameterG0=190l = 0.744MCS = 1.727d0 = 2.267v = 0.3k = 0.008a′  = 0.032cψ = 0.98β = 12.8b′ = 0.527h1 = 0.46c′ = 0.313h2 = 0.78kh = 0.12kp = 4.95λc = 0.14kd = 0.458


#### 4.5.1. Elastic Parameter

G0 was the material constant. It can be obtained from the three-dimensional stress-strain curve of the rockfill material. v was Poisson’s ratio reflecting the transverse deformation of the material and can be obtained by the triaxial compression test.


#### 4.5.2. Fractal Parameter

l, κ, and β were the fractal parameters. The initial fractal dimension D0 and fractal dimension D could be obtained from equation (2). Then according to equation (3), the fractal parameters could be obtained from fitting the fractal dimension D and σ3/Pa through the binary linear regression analysis.


#### 4.5.3. Critical Parameter

According to the measured initial void ratioe0, the values of a′, b′, and c′ were obtained from equations (21) and (22). For rockfill materials, the value of the parameter λc was obtained by normalizing mean effective stress after obtaining the critical state line through triaxial test data. Regarding the other critical parameters kh and Mcs, for details, please refer to relevant research for determination [32, 34].


#### 4.5.4. State Parameter

d0 was the expansion coefficient which can be obtained from the dilatancy equation and conventional triaxial test of rockfill material. kp and kd were model parameters that can be obtained from equations (24) and (30) in the test, respectively. cψ was defined as cψ=cosθσ−θdε. It can be obtained from the stress Lode angle in the critical state and the principal stress space strain increment Lode angle. h1 and h2 were model parameters, please refer to the relevant research for details [34].


## 4.1. Property-Dependent Plastic Potential Theory

To reflect the anisotropy of the rockfill fabric and the uniqueness of the rockfill critical state, the consideration of rockfill mesoscopic fabric material properties related to the plastic potential theory was put forward based on material status-related dilatancy theory. It was assumed that the material properties could be represented by the fabric tensor and that the rockfill material properties would influence the strain distribution. For an incremental property of the geotechnical materials, the general expression of the property-dependent plastic potential related to rockfill materials can be written as(15)dεijp=dλ∂g∂σilF¯lj,where g is the plastic potential function, dλ is the increment of the plastic factor, and F¯ij is the tensor of the fabric tensor Fij after the transformation of equation (15). According to the transformation characteristics of equation tensors, equation (15) was further simplified to(16)Fij=pFδij+s˜ij=pFF¯ij,where pF and s˜ij are the spherical tensor and the partial tensor of the fabric tensor, respectively.Equations (15) and (16) could be rearranged as another expression of the material property potential theory:(17)dεijp=dλ∂g∂σij+dλ∂g∂σils˜lj,where dλ∂g/∂σij is traditional plastic potential theory. When the material was isotropic, dλ∂g/∂σijs˜lj was equal to 0, and the plastic factor dλ was a nonnegative proportional coefficient in equation (17).


## 4.2. Hardening Law

The hardening law of internal variables was adopted in the incremental hyperbolic form proposed by Li and Dafalias [32]. It was given by (18)dεsp=MhsGMp−MdM,where G is the shear modulus. Mp is the peak stress ratio in a conventional triaxial compression test, and hs is the material constant.It must be noted that the peak stress ratioMp had to not be a constant, but rather a quantity related to the state parameters ψ. In this research, the equation for rockfill was adopted, Mp=Mcsexp−kpψ, as suggested by Li and Dafalias, and kp was the constant of the model. The state parameter ψ was expressed as(19)ψ=e−ec,where ec is the critical void ratio which could be determined from the following equation:(20)ec=eΓ−λcppaξ,The results showed that there was a good linear relationship between the critical void ratio intercepteГ and the initial void ratio e0, which satisfied the relationship(21)eΓ=c′+b′e0.Moreover, there was a linear strip relationship between the critical void ratio intercepteГ and the fractal dimension D, so the critical void ratio intercept eГ was expressed as(22)eΓ=c′−a′D+b′e0,where a′, b′, and c′ are the model parameters that could be calculated from the experimental data.For the value ofhs, Li and Dafalias [32] considered that it was related to the void ratio, taken as hs = h1−h2e, where h1 and h2 were the constants of the model. In addition, the model parameters were also affected by anisotropy, and its expression could be rearranged as(23)hs=h1−h2e1+khA,where kh is the model parameter and A is the defined anisotropic state variable in equation (12), when the material is isotropic, A≡0.To consider the effect of the anisotropy, the expression of the peak stress ratioMp could be changed into the function of an anisotropic state variable:(24)Mp=Mcsgθσe−kpψA,where ψA is the state parameter related to the initial anisotropy described in equation (24) and kp is the model parameter. The expression of the peak stress ratio showed that the peak stress ratio was not a constant but rather an obvious function of the anisotropic state variable.


## 4.3. Dilatancy Equation

The dilatancy equation could be deduced from property-dependent plastic potential theory. Based on the basic theory of rockfill mechanics in a critical state, when the rockfill reached the critical state, the volume change was equal to zero and the stress ratio reached the critical state pressure ratioηc. Lü et al. found that stress dilatation is a function related to the state of geomaterials [33]. At this point, the deformation of the rockfill only had shear strain, so its energy relationship satisfied(25)dWp=pdεvp+qdγpcosθσ−θdε=ηcpdγp,where ηc is the stress ratio under the critical state. According to the energy relationship of the critical state given in equation (25), the dilatancy equation based on the potential theory related to material properties could be obtained as follows:(26)dεvpdγp=ηc−cosθσ−θdεqp.The dilatancy equation is usually determined by conventional triaxial tests, and it can be expressed as follows:(27)d=d0ηc−cosθσ−θdεqp,where d0 is the dilatation coefficient, which was determined by the experiment, θσ is the stress Lode angle, and θdε is the strain increment Lode angle in principal stress space. When the stress ratio was greater than the critical state (d < 0) for the dilatancy equation, rockfill dilatancy occurred. When the stress ratio was less than the critical stress (d > 0), the rockfill contracted. When the stress ratio was equal to the critical state (d = 0), the rockfill was in the state of phase, as shown in Figure 1.Figure 1
Relationship between the critical state stress ratio and dilatation.In the principal stress space, the strain incrementθdε of the Lode angle was defined as(28)θdε=arctan132dε2−dε1−dε3dε1−dε3,where the main strain was determined by property-dependent plastic potential theory in equation (25). Assuming q/p=η, equation (29) could be rearranged as(29)dq=ηdp+pdη,By substituting equation (29) into equation (27), we obtained the following:(30)dpp=−dηd+η,The plastic potential function was obtained by integrating equation (30):(31)g=q−pd0ηcd0cψ−11−pp0d0cψ−1=0,where cψ is the state parameter, the value could be obtained from the equation:cψ=cosθσ−θdε under the critical state of rockfill.Considering the influence of anisotropy on the critical state stress ratioηc, the expression of ηc in 3D space was defined as follows:(32)ηc=1+ςAMcsgθσekdψA,where kd is the model parameter.


## 4.4. Incremental Stress-Strain Relationship

Based on the expression method of the strain increment, the total strain could be expressed as(33)dεij=dεije+dεijp,where dεij is the total strain increment and dεije is the elastic strain increment, the expression for which is shown in equation (34). According to the potential theory related to the material properties, dεijp is the plastic strain increment, which is expressed in equation (17), as follows:(34)dεije=Cijkledσkl.According to equations (17), (33), and (34), the total strain increment was derived as dσij=Dijkledεkl−dλ∂g/∂σisFsj:(35)dεij=Cijkledσkl+dλ∂g∂σikFkj.Equation (35) was simplified as follows:(36)dσij=Dijkledεkl−dλ∂g∂σisFsj.Equation (40) could be rearranged as(41)dλ=∂f/∂σijDijkledεklAp+∂f/∂σmnDmnpqe∂g/∂σsuFut.where Dijkle is the elastic stiffness matrix. dλ is also related to the hardening parameter Ap, so it could be rearranged as(37)dλ=1Ap∂f∂σijdσij.According to the consistency equation df equal to 0 on the yield surface, the equation became(38)df=∂f∂σijdσij+∂f∂MdM=0,where f is the yield function and M is the hardening parameter.Based on Equations (37) and (38), this led to the expression below:(39)dλAp=−∂f∂M∂M∂εijpdεijp.After combining equations (36)–(39), the expression was given by(40)∂f∂σijDijkledεkl−∂f∂σmnDmnpqedλ∂g∂σsupFF¯ut−Apdλ=0.Equations (37) and (38) were combined as follows:(42)Ap=−∂f∂MdMdλ=−∂f∂M∂M∂εijp∂g∂σikFkj.Because the hardening law took the plastic shear strain as the hardening parameter in the form of an incremental hyperbola, it could be obtained according to the flow law:(43)Ap=−∂f∂MdMdλ=−∂f∂M∂M∂εqp∂g∂qFkj2+1q∂g∂θσFkj2.The hardening law equation (18) was substituted into equation (43) to obtain the hardening function as(44)Ap=gθσphsGMp−M2pMp2∂g∂qFkj2+1q∂g∂θσFkj2.


## 4.5. Model Parameter Determination

There were four groups and seventeen model parameters in this research. The elastic parameters and the fractal parameters could be measured using the basic test data of the rockfill materials. The critical state parameters were tested with a triaxial tensile test. The state-related parameters were determined with the drained triaxial test. All model parameters are shown in Table1. The constitutive model parameters were divided into four groups: elastic parameter, fractal parameter, critical parameter, and state parameter. The specific methods for determining these parameters were as follows.Table 1
Parameters of the state-dependent constitutive model for rockfill.
Elastic parameterFractal parameterCritical parameterState parameterG0=190l = 0.744MCS = 1.727d0 = 2.267v = 0.3k = 0.008a′  = 0.032cψ = 0.98β = 12.8b′ = 0.527h1 = 0.46c′ = 0.313h2 = 0.78kh = 0.12kp = 4.95λc = 0.14kd = 0.458


### 4.5.1. Elastic Parameter

G0 was the material constant. It can be obtained from the three-dimensional stress-strain curve of the rockfill material. v was Poisson’s ratio reflecting the transverse deformation of the material and can be obtained by the triaxial compression test.


### 4.5.2. Fractal Parameter

l, κ, and β were the fractal parameters. The initial fractal dimension D0 and fractal dimension D could be obtained from equation (2). Then according to equation (3), the fractal parameters could be obtained from fitting the fractal dimension D and σ3/Pa through the binary linear regression analysis.


### 4.5.3. Critical Parameter

According to the measured initial void ratioe0, the values of a′, b′, and c′ were obtained from equations (21) and (22). For rockfill materials, the value of the parameter λc was obtained by normalizing mean effective stress after obtaining the critical state line through triaxial test data. Regarding the other critical parameters kh and Mcs, for details, please refer to relevant research for determination [32, 34].


### 4.5.4. State Parameter

d0 was the expansion coefficient which can be obtained from the dilatancy equation and conventional triaxial test of rockfill material. kp and kd were model parameters that can be obtained from equations (24) and (30) in the test, respectively. cψ was defined as cψ=cosθσ−θdε. It can be obtained from the stress Lode angle in the critical state and the principal stress space strain increment Lode angle. h1 and h2 were model parameters, please refer to the relevant research for details [34].


## 4.5.1. Elastic Parameter

G0 was the material constant. It can be obtained from the three-dimensional stress-strain curve of the rockfill material. v was Poisson’s ratio reflecting the transverse deformation of the material and can be obtained by the triaxial compression test.


## 4.5.2. Fractal Parameter

l, κ, and β were the fractal parameters. The initial fractal dimension D0 and fractal dimension D could be obtained from equation (2). Then according to equation (3), the fractal parameters could be obtained from fitting the fractal dimension D and σ3/Pa through the binary linear regression analysis.


## 4.5.3. Critical Parameter

According to the measured initial void ratioe0, the values of a′, b′, and c′ were obtained from equations (21) and (22). For rockfill materials, the value of the parameter λc was obtained by normalizing mean effective stress after obtaining the critical state line through triaxial test data. Regarding the other critical parameters kh and Mcs, for details, please refer to relevant research for determination [32, 34].


## 4.5.4. State Parameter

d0 was the expansion coefficient which can be obtained from the dilatancy equation and conventional triaxial test of rockfill material. kp and kd were model parameters that can be obtained from equations (24) and (30) in the test, respectively. cψ was defined as cψ=cosθσ−θdε. It can be obtained from the stress Lode angle in the critical state and the principal stress space strain increment Lode angle. h1 and h2 were model parameters, please refer to the relevant research for details [34].


## 5. Test Verification

### 5.1. Introduction to the Experiment

The experiment adopted the consolidated drained triaxial test results of Cai [23] by using a large triaxial instrument. The test specimen was dolomitic limestone. The density of the particle was 2.77 g/cm3 and the particle size was less than 800 mm. Furthermore, the grading was good. The nonuniformity coefficients Cu and Cc were 35.48 and 1.35, respectively. There were four types of gradation for the rockfill samples and the relative densities were 0.65, 0.75, 0.90, and 1.0. To characterize the strength characteristics of the rockfill under different relative compactness and confining pressures, the triaxial shear test was conducted under four confining pressures: 300 kPa, 600 kPa, 1000 kPa, and 1500 kPa.


### 5.2. Experimental Characteristics and 3D Description of the Fractal Dimension

Tapias et al. found that grain crushing largely resulted in the deformation of rockfill particles under different pressures [35]. The experiment showed that when the contact force was small, contact-related crushing dominated, although equivalent volume division dominated beyond given yield stress. The experimental characteristics of the fractal dimension reproduced these observations and trends. Figure 2 shows the relationship curves between the grain gradation and the fractal dimension under the confining pressure of 1000 kPa for Dr equal to 0.90. The initially graded 1 to 4 samples showed obvious particle fragmentation after the test. The particle fragmentation was mainly concentrated in the range of 80 mm to 40 mm. The results of particle crushing were about 5%–7% higher than the average particle mass percentage before crushing. The fractal dimension of crushing was proportional to the confining pressure. The dependence of stress path was also one of the important characteristics of rockfill material properties [36, 37]. The large-scale triaxial compression test showed that the stress path had a hierarchical evolution effect on the critical state and the particle fragmentation of the rockfill materials.Figure 2
Relationship between the gradation and the fractal dimension. (a) The ratio of the particle size and the mass percentage content. (b) The curve of the fractal dimension withσ3/Pa.
(a)(b)Besides, Guo et al. found that the particle breakage and particle size distribution had an important effect on the critical state line of the rockfill under stress [38–40]. The critical state line of rockfill was very similar to sand. The effective mean normal stress p′, deviatoric stress q, and volume changes of the rockfill were no longer changing after reaching the critical state. Then, the continuous increment of the shear deformation would cause different degrees of damage. However, because of the different initial material properties (particle size, initial gradation, etc.) of the rockfill, the critical state was difficult to determine accurately. Therefore, we would approximate achieve the critical state of the rockfill when the deviatoric stress and volume deformation were stable in the shear test. So, it is generally believed that the stable state was when the shear strain was greater than 20%. Figure 3 shows that the critical state point of the rockfill material is linear and unique in the p′∼q space. This is consistent with the research results on critical state line of rockfill material [40–42].Figure 3
Critical state point and trend line inp′∼q space.Figure4 shows the e∼p′/paξ curve of the rockfill. The slope of the critical void ratio of the four-gradation rockfill was essentially the same, but the y-axis intercepts were not the same, which reflected the compaction of the rockfill with different levels and different relative densities. The experimental results showed that the critical void ratio intercept was closely related to the initial void ratio and the fractal dimension. Figure 5 shows that the critical void ratio intercept had a good linear relationship with the initial void ratio and the fractal dimension, which further reflected the necessity of the fractal dimension when describing the void of the rockfill particles.Figure 4
The critical state point and trend line ine∼p′/paξ space. (a) Gradation 1. (b) Gradation 2. (c) Gradation 3. (d) Gradation 4.
(a)(b)(c)(d)Figure 5
The properties of the granular materials. (a) The relation between the intercepteГ and the initial void ratio e0. (b) The relation between the intercept eГ and the fractal dimension D.
(a)(b)Figure6(a) reflects the failure surfaces of the rockfill materials at different levels. The figure indicates that the gradation, fractal dimension, and stone-piling strength were correlated. The greater the inhomogeneity coefficient of the rockfill gradation was, the greater the fractal dimension was and the greater the yield strength was. Figures 6(b) and 6(c) show the variation rules of the strength parameters and the strength failure surface on the plane, respectively. The figures show that the strength failure lines on the π plane were similar shapes in equilateral curved triangles and they were symmetrical for the condition in which the orthotropic amplitude parameters a1 and a2 were equal. With an increase of the Ac or Bc, the strength failure line on the π plane increased along the coordinate axis. The increase of the strength parameter Ac essentially reflected the compactness of the arrangement of particles with the decrease of the void ratio. The aggregation of the particle action showed the characteristic of the stronger yield strength macroscopically. However, the increase of the strength parameter Bc made the yield strength decrease because the parameter Bc was directly proportional to the void ratio. When the amplitude a1 and the amplitude a2 changed, the yield strength in the vicinity of σ3 and σ2 changed in response.Figure 6
Schematic of various failure surfaces. (a) Failure surface with different gradations. (b) Failure surface with variations of the strength parameterAc. (c) Failure surface with variations of the strength parameter Bc.
(a)(b)(c)Figure7 shows the relationship between the void ratio and the intermediate principal stress coefficient. It is shown in Figure 7(a) that the degree of anisotropy of the rockfill materials decreased with the increase of the intermediate principal stress coefficient from the microscopic characteristics of the rockfill particles. However, on the macro level, the critical void ratio decreased and the yield strength increased. Figure 7(b) shows the denser arrangement of rockfill particles on the microscopic level when the gradation nonuniformity coefficient and the initial fractal dimension were larger under the same pressure. The test results showed that the critical void ratio increased with the yield strength.Figure 7
(a) Void ratio at various coefficient of intermediate principal stress. (b) Void ratio at various gradations.
(a)(b)


### 5.3. Deformation Simulation of the Triaxial Drainage Test

The test results and the model simulation results of three kinds of rockfill gradation were compared atDr = 0.75, as shown in Figure 8. It can be seen from the stress-strain curve that the constitutive equation established in this research could better predict the dilatancy characteristics of the rockfill. The gradation unevenness coefficient and the initial fractal dimension had a positive correlation to the impact of the rockfill swelling. Moreover, as deviatoric stress σ1−σ3 increased at the peak state, the softening of the rockfill became more obvious. The test results showed that the confining pressure σ3 had an important influence on the process of the softening to the hardening of the rockfill. Specifically, the volume from contractive to dilative was increasingly obvious with the increase of the peak value of the deviatoric stress. The finer the gradation was (the initial fractal dimension D0 was larger), the greater the peak value of the deviatoric stress σ1−σ3 was. At that time, it changed from the hardening type to the softening type, and the volume changed from contractive to dilative.Figure 8
The comparison of test results and the simulation results (Dr = 0.75). (a) Gradation 1. (b) Gradation 2. (c) Gradation 3.
(a)(b)(c)


## 5.1. Introduction to the Experiment

The experiment adopted the consolidated drained triaxial test results of Cai [23] by using a large triaxial instrument. The test specimen was dolomitic limestone. The density of the particle was 2.77 g/cm3 and the particle size was less than 800 mm. Furthermore, the grading was good. The nonuniformity coefficients Cu and Cc were 35.48 and 1.35, respectively. There were four types of gradation for the rockfill samples and the relative densities were 0.65, 0.75, 0.90, and 1.0. To characterize the strength characteristics of the rockfill under different relative compactness and confining pressures, the triaxial shear test was conducted under four confining pressures: 300 kPa, 600 kPa, 1000 kPa, and 1500 kPa.


## 5.2. Experimental Characteristics and 3D Description of the Fractal Dimension

Tapias et al. found that grain crushing largely resulted in the deformation of rockfill particles under different pressures [35]. The experiment showed that when the contact force was small, contact-related crushing dominated, although equivalent volume division dominated beyond given yield stress. The experimental characteristics of the fractal dimension reproduced these observations and trends. Figure 2 shows the relationship curves between the grain gradation and the fractal dimension under the confining pressure of 1000 kPa for Dr equal to 0.90. The initially graded 1 to 4 samples showed obvious particle fragmentation after the test. The particle fragmentation was mainly concentrated in the range of 80 mm to 40 mm. The results of particle crushing were about 5%–7% higher than the average particle mass percentage before crushing. The fractal dimension of crushing was proportional to the confining pressure. The dependence of stress path was also one of the important characteristics of rockfill material properties [36, 37]. The large-scale triaxial compression test showed that the stress path had a hierarchical evolution effect on the critical state and the particle fragmentation of the rockfill materials.Figure 2
Relationship between the gradation and the fractal dimension. (a) The ratio of the particle size and the mass percentage content. (b) The curve of the fractal dimension withσ3/Pa.
(a)(b)Besides, Guo et al. found that the particle breakage and particle size distribution had an important effect on the critical state line of the rockfill under stress [38–40]. The critical state line of rockfill was very similar to sand. The effective mean normal stress p′, deviatoric stress q, and volume changes of the rockfill were no longer changing after reaching the critical state. Then, the continuous increment of the shear deformation would cause different degrees of damage. However, because of the different initial material properties (particle size, initial gradation, etc.) of the rockfill, the critical state was difficult to determine accurately. Therefore, we would approximate achieve the critical state of the rockfill when the deviatoric stress and volume deformation were stable in the shear test. So, it is generally believed that the stable state was when the shear strain was greater than 20%. Figure 3 shows that the critical state point of the rockfill material is linear and unique in the p′∼q space. This is consistent with the research results on critical state line of rockfill material [40–42].Figure 3
Critical state point and trend line inp′∼q space.Figure4 shows the e∼p′/paξ curve of the rockfill. The slope of the critical void ratio of the four-gradation rockfill was essentially the same, but the y-axis intercepts were not the same, which reflected the compaction of the rockfill with different levels and different relative densities. The experimental results showed that the critical void ratio intercept was closely related to the initial void ratio and the fractal dimension. Figure 5 shows that the critical void ratio intercept had a good linear relationship with the initial void ratio and the fractal dimension, which further reflected the necessity of the fractal dimension when describing the void of the rockfill particles.Figure 4
The critical state point and trend line ine∼p′/paξ space. (a) Gradation 1. (b) Gradation 2. (c) Gradation 3. (d) Gradation 4.
(a)(b)(c)(d)Figure 5
The properties of the granular materials. (a) The relation between the intercepteГ and the initial void ratio e0. (b) The relation between the intercept eГ and the fractal dimension D.
(a)(b)Figure6(a) reflects the failure surfaces of the rockfill materials at different levels. The figure indicates that the gradation, fractal dimension, and stone-piling strength were correlated. The greater the inhomogeneity coefficient of the rockfill gradation was, the greater the fractal dimension was and the greater the yield strength was. Figures 6(b) and 6(c) show the variation rules of the strength parameters and the strength failure surface on the plane, respectively. The figures show that the strength failure lines on the π plane were similar shapes in equilateral curved triangles and they were symmetrical for the condition in which the orthotropic amplitude parameters a1 and a2 were equal. With an increase of the Ac or Bc, the strength failure line on the π plane increased along the coordinate axis. The increase of the strength parameter Ac essentially reflected the compactness of the arrangement of particles with the decrease of the void ratio. The aggregation of the particle action showed the characteristic of the stronger yield strength macroscopically. However, the increase of the strength parameter Bc made the yield strength decrease because the parameter Bc was directly proportional to the void ratio. When the amplitude a1 and the amplitude a2 changed, the yield strength in the vicinity of σ3 and σ2 changed in response.Figure 6
Schematic of various failure surfaces. (a) Failure surface with different gradations. (b) Failure surface with variations of the strength parameterAc. (c) Failure surface with variations of the strength parameter Bc.
(a)(b)(c)Figure7 shows the relationship between the void ratio and the intermediate principal stress coefficient. It is shown in Figure 7(a) that the degree of anisotropy of the rockfill materials decreased with the increase of the intermediate principal stress coefficient from the microscopic characteristics of the rockfill particles. However, on the macro level, the critical void ratio decreased and the yield strength increased. Figure 7(b) shows the denser arrangement of rockfill particles on the microscopic level when the gradation nonuniformity coefficient and the initial fractal dimension were larger under the same pressure. The test results showed that the critical void ratio increased with the yield strength.Figure 7
(a) Void ratio at various coefficient of intermediate principal stress. (b) Void ratio at various gradations.
(a)(b)


## 5.3. Deformation Simulation of the Triaxial Drainage Test

The test results and the model simulation results of three kinds of rockfill gradation were compared atDr = 0.75, as shown in Figure 8. It can be seen from the stress-strain curve that the constitutive equation established in this research could better predict the dilatancy characteristics of the rockfill. The gradation unevenness coefficient and the initial fractal dimension had a positive correlation to the impact of the rockfill swelling. Moreover, as deviatoric stress σ1−σ3 increased at the peak state, the softening of the rockfill became more obvious. The test results showed that the confining pressure σ3 had an important influence on the process of the softening to the hardening of the rockfill. Specifically, the volume from contractive to dilative was increasingly obvious with the increase of the peak value of the deviatoric stress. The finer the gradation was (the initial fractal dimension D0 was larger), the greater the peak value of the deviatoric stress σ1−σ3 was. At that time, it changed from the hardening type to the softening type, and the volume changed from contractive to dilative.Figure 8
The comparison of test results and the simulation results (Dr = 0.75). (a) Gradation 1. (b) Gradation 2. (c) Gradation 3.
(a)(b)(c)


## 6. Conclusions

Because of the research status of large particle sizes, compression crushing characteristics, and lack of relatively large-scale triaxial test results and the difficult description on the constitutive equation, this research established a constitutive model based on the test results of rockfill. Additionally, the research method for combining the macro and the micro effects of particle gradation with its anisotropy was described, along with the plastic state theory of the material state. The main conclusions were as follows:(1)
A constitutive model of rockfill was established based on property-dependent plastic potential theory for geomaterials. The analysis of the model could be used to characterize the particle gradation of the rockfill and the law of particle breakage under loading conditions. At the same time, the microscopic arrangement and macroscopic yield characteristics of the rockfill were connected by the fabric tensor. The hardening law could better describe the anisotropic mechanical response of the mesoscopic fabric, fractal dimension, particle breakage, and stress state to the rockfill. Modeling based on property-dependent plastic potential theory for geomaterials made the model give a better description of the mechanical characteristics of the principal stress axis rotation caused by complex stress conditions.(2)
A conventional rockfill triaxial drainage test was used to verify the model. The new model established the relationship between the mesoscopic particles and the stress state of the rockfill and used anisotropic state variables to reflect the anisotropy caused by the particle arrangement, such as initial anisotropy and stress-induced anisotropy. The model also considered the influences of the fractal dimension, particle fragmentation, gradation, and other factors on the critical state line of the rockfill. The performance of the model was verified by the test results in describing the critical state under various stress paths. In this research, the intensity change in three dimensions and its intensity change rule on theπ plane were simulated, and the validity and the rationality of the model were verified with the simulation.(3)
The triaxial tests and simulations confirmed that the fractal dimension could reflect the crushing of the rockfill particles and the scale of the rockfill during the loading of the stress space internal force, as well as the gradation change of the rockfill and its critical void ratio and yield strength deformation. In this research, based on the fractal dimension theory, the two material properties of sand and rockfill were organically linked, which could not only quantitatively describe the particle grading, but also clearly express the void ratio strength parametersAc and Bc.(4)
The experimental results showed that the gradation and the initial fractal dimension determined the material properties of the rockfill, the peak value of the deviatoric stress had an important effect on the softening of the rockfill, and the influence of the confining pressure was mainly in the process from contractive to dilative. The constitutive model accurately simulated this result. However, based on the limited experimental data and the need for practical application, the constitutive equation needs to be further improved.



---

*Source: 8868600-2020-12-03.xml*