Study on Creep Characteristics and Constitutive Relation of Fractured Rock Mass

In order to effectively describe the whole creep process of fractured rock mass, triaxial unloading creep tests were carried out on prefractured coal samples using constant axial pressure and graded unloading confining pressure, and the axial and lateral creep laws of fractured coal samples with different dip angles were studied. Combined with the characteristics of creep curve and based on Kachanov’s creep damage theory, the damage variable is introduced into the constitutive relation and creep equation, and the evolution equation of damage variable with time in the whole creep process is derived. At the same time, a newmethod to calculate the initial damage is proposed. *e elastoplastic body with damage variable is connected with the Burgers model in series. Meanwhile, lade criterion and switch element are introduced into the creepmodel to establish a new fracture damage creepmodel. *e one-dimensional and three-dimensional damage creep equations are derived. *e damage creep equation is obtained according to the superposition principle. A simple and feasible method for parameter identification of the model is given based on the characteristics of creep curve. *e applicability of the model is verified by comparing the creep test curve of fractured coal sample with the theoretical curve.*e results show that the two models are in good agreement.*e model can not only accurately reflect the nonlinear characteristics of creep curves in the attenuation and isokinetic stages but also describe the accelerated creep characteristics of fractured rocks.


Introduction
Rock burst refers to the local high stress concentration in coal (rock) caused by the mining of working face and the tunneling of roadway, which makes the local elastic energy in coal (rock) accumulate stably in time and transfer to weak surface through easy propagation path, resulting in continuous generation of weak structure in coal and rock, causing local fracture, damage, instability, and sliding of coal (rock). e phenomenon that the elastic potential energy stored in coal (rock) is released instantaneously after some external disturbance [1][2][3][4][5]. is phenomenon can only last for a few seconds to tens of seconds, which will produce severe vibration, accompanied by a large number of coal (rock) thrown out, seriously threatening the safety of coal production [6][7][8]. With the increase of coal mining depth, when the rock stratum is relatively soft and the coal seam with strong risk is impacted, the strength, failure range, and destructive force of shock wave of rock burst have obvious change characteristics. Rock burst will also occur under the condition of no disturbance. e sudden release of this type of energy is because of creep-type rock burst [9][10][11]. erefore, it is of great significance for the prevention and control of coal body impact disaster that the creep-type impact risk can be predicted by quickly judging the creep deformation of coal body on-site. e existing prevention and control methods of rock burst mainly include the following: (1) energy releasing drilling method; (2) high-pressure water injection in coal seam; (3) deep hole controlled blasting; and (4) hydraulic slotting pressure relief method [12][13][14][15]. Li Chao et al. used the method of cutting coal body with high-pressure water jet to prevent rock burst and monitored the test results with KBD5 electromagnetic radiation monitor. It shows that the effect of cutting on stress unloading of shallow coal is the best, and the stress is redistributed around the seam, so that the high stress is transferred from shallow to deep. e existing prevention methods of rock burst are to release the accumulated energy in coal and rock mass by releasing confining pressure, so as to reduce the stress concentration in coal and rock mass and prevent rock burst [16][17][18][19][20]. e rheological property of rock is an important factor affecting the long-term stability of rock mass engineering.
ere are many factors affecting the creep process of rock, and one of the important factors is the influence of rock damage on creep. In the process of rock creep, the internal damage cracks of rock are produced continuously, which makes the bearing capacity of rock decrease continuously, making the rock to enter the accelerated creep stage. In addition, under the coupling effect of engineering disturbance and rock damage fracture penetration and fusion, roof falling, rock burst, and other disasters are likely to occur. erefore, it is necessary and significant to introduce damage body to simulate the whole creep process curve of rock. At present, the existing rock creep model simulation has made great achievements, but there are still many problems in the simulation of accelerated creep process [21][22][23][24][25]. erefore, this paper establishes a creep model to simulate accelerated creep process by introducing a damage body.
In order to simulate the accelerated creep process of rock, scholars at home and abroad have done a lot of research on nonlinear rheological theory and achieved good results: one is that three different material element models (such as Hooke body, Saint Venant body, and Newton body) can be combined in series and parallel to establish a new creep constitutive model, and a lot of rock creep can be established by this method of constitutive models, such as Maxwell model, Kelvin model, and Xiyuan model; Zhao Yanlin constructed a new rheological model, which was composed of components in series, and could well describe the creep law under the state of graded stress; Xia Caichu et al. proposed a unified rheological mechanical model including four basic rheological mechanical forms: viscoelastic, viscoplastic, viscous, and viscoelastic-plastic, and gave the identification method of the model.
Although there are many methods to construct traditional components, they still have some limitations. ey can only reflect the characteristics of viscoelastic-plasticity and cannot simulate the process of accelerated creep well. erefore, some nonlinear models are proposed. e others are based on the damage theory because the basic parameters of creep can be changed, so it can be realized by using the damage factor. e new component is combined with the traditional element in series, and the creep constitutive equation which can simulate the whole process of creep is established. Murakammi and Ohno proposed a damage model for jointed rock mass with multiple sets of planar fractures based on the description of the geometric size of the fracture system; Kawamoto et al. further hypothesized and improved the abovementioned model for jointed rock mass; Zhu Weishen and Qi Yinping introduced a secondorder tensor to describe the damage of rock mass based on the geometric damage principle and the definition of the abovementioned damage variables. Based on Kachanov's damage theory, Wu Zhulin et al. simulated the initial damage factor in series with Burgers element model and achieved good results. However, most of the existing creep models are based on the combination of original or improved components. Either there are many components, or the formula is complex, or cannot be solved. It can be seen that the simulation method of accelerated creep stage still has some shortcomings, which needs to be improved. erefore, the damage factor is introduced into the creep model to simulate the whole creep process. ere are two methods of damage factor: the first is to define the damage factor from the angle of energy damage; the second is to define the damage factor based on the effective bearing area of the structure in terms of geometric damage. However, the nonlinear characteristics of rock creep are not only in the accelerated creep stage but also in the whole creep process. However, most of the creep damage models focus on the accelerated creep process after the stress level exceeds the long-term strength of the rock, and there is little research on the nonlinear characteristics of the early stable creep. In this paper, the initial damage is defined on the basis of existing damage mechanics, and the evolution equation of damage variable is derived according to the Kachanov damage rate. At the same time, based on the lade failure criterion, the elastoplastic damage element is introduced to simulate the plastic damage body and then is combined with the Burgers model in series, a constitutive model with few parameters, easy to determine, uncomplicated components, and able to simulate accelerated creep is constructed. e research in this paper can be used for reference for the stability of underground engineering in the future, especially for the stability of support in deep mining and the prevention of rock burst [26][27][28][29].

Specimen Preparation.
e coal sample was taken from Wang JiaLing coal mine in Shanxi Province. e coal is hard and anthracite. In order to avoid the influence of anisotropy on the test results, rock samples were processed along the same direction when preparing coal samples. At the same time, in order to ensure the authenticity and comparability of the test results, the appearance of the processed coal samples was carefully observed before the test to determine that there were no obvious weak planes such as joints and cracks so as to ensure that there is no obvious difference between the test rock samples in macroscopic. According to the requirements of the International Society for Rock Mechanics (ISRM), the rock samples are processed into standard cylindrical rock samples with a diameter of 50 mm × 100 mm.
e precast fractured rock sample is shown in Figure 1, in which the fracture length is 2a and the fracture dip angle is α. e high-speed electric cutting machine is used to process three-dimensional cracks. e cutting wheel is 0.3 mm thick ultrathin diamond saw blade. e crack thickness is 0.3-0.5 mm, and the crack is filled with gypsum. Change the fracture dip angle α (30°, 45°, and 60°), fracture length 2a � 20 mm, and fracture spacing 2b � 30 mm.

2
Shock and Vibration

Test System and Loading Process.
e triaxial unloading creep test was performed using the ROCK 600-50 triaxial test system manufactured by Topouan, France ( Figure 2). e function and application scope of the equipment are as follows: it can be used for mechanical experiment of rock, rock-soil, and concrete under normal condition and temperature bar with a sample size of φ 50 mm × 100 mm, rheological test of rock mechanics, permeability test of rock mechanics, and mechanical coupling test of rock under conventional condition and temperature condition. e main technical indexes of the equipment are as follows: the maximum axis is 380 MPa, the maximum confining pressure is 60 MPa, the maximum pore pressure is 40 MPa; the axial range of deformation measurement is 0-10 mm, and the circumferential range is 0-4 mm; and the measurement resolution is 0.00001 mm, which fully meets the working requirements of long-time high load.
In the test, constant axial pressure is used to release confining pressure by stages, and load control mode is adopted. Firstly, according to the hydrostatic pressure condition σ 1 � σ 2 � σ 3 � 10 MPa, confining pressure is applied to 10 MPa at the loading rate of 0.5 MPa/s, and then the axial load is continuously applied to 75% of the predicted ultimate bearing capacity under the confining pressure, which exceeds the uniaxial compressive strength of the coal sample and then maintains the principal stress difference σ 1 -σ 3 . Under the same condition, the confining pressure is gradually reduced by 0.5 MPa/s, and the unloading is carried out step by step from 10 MPa to 8 MPa, 6 MPa, 4 MPa, and 2 MPa; each stage is maintained for 24 h until the coal sample is completely destroyed. e specific control parameters are shown in Table 1.

Test Results.
In this section, the unloading creep test results of prefractured coal samples obtained by using the triaxial rheometer are analyzed, and the rheological deformation and failure law of coal samples under different fracture dip angles are studied. In order to ensure the reliability of the test data, the creep test is divided into three groups, four specimens are tested in each group, and the creep time of each stage is 24 h. e test results are shown in Table 2.
e creep law of the prefractured coal sample is obtained through the creep test of unloading confining pressure step by step under constant axial pressure, as shown in Figure 3. e creep curve obtained under the condition of step unloading creep is ladder shape, which cannot be directly used and needs to be transformed. According to Chen's loading treatment method, the whole curve is transformed into a creep curve under various loads, as shown in  It can be seen from     Complete coal sample  45  10  8  6  4  2  A30  30°angle  40  10  8  6  4  2  A45  45°angle  35  10  8  6  4  2  A60  60°angle  33  10  8  6  4 2 Table 2: Triaxial unloading creep test data of prefabricated cracked coal samples.

Shock and Vibration
Number  1 and σ 3 are the axial pressure and confining pressure; ε i1 , ε i3 , and ε iv are the axial, lateral, and volumetric instantaneous strain increments, respectively; ε c1 , ε c3 , and ε cv are the axial, lateral, and volumetric creep strain increments, respectively.  Shock and Vibration

Analysis on Creep Deformation Law of Complete Coal
Sample. Taking A00 sample as an example, the instantaneous deformation and creep deformation of A00 sample increased with the increase of deviator stress (Figures 3 and 4). e larger the deviator stress is, the greater the increase range of creep deformation is. When the confining pressure is 10 MPa, obvious instantaneous strain is produced. With the decrease of confining pressure, the axial and lateral instantaneous strain values increase gradually. When the confining pressure is reduced from 10 MPa to 2 MPa, the increment of axial instantaneous strain is about 1.77 × 10 − 2 , and the lateral increment is about 3.5 × 10 − 4 . Since the second stage unloading, the axial and lateral creep deformation of each stress level is significantly increased, and the axial strain increase is about 0.5 × 10 − 2 . When the confining pressure is 2 MPa, the axial and lateral creep deformation increases significantly, and the accelerated creep phenomenon appears. Within 1 hour, the axial strain increased by 4.05 × 10 − 4 in the direction almost parallel to the axial strain axis, and the lateral deformation increased by 0.2 × 10 − 4 . e results show that the internal bearing capacity of the sample is weakening with time, and the through fracture surface under the level of fracture confining pressure leads to the failure of the sample.

Damage Analysis.
In the process of rock formation, there are a lot of microdefects (joints, holes, etc.) in the process of rock formation, so the initial damage is not zero before the creep test. Damage variable is used to measure the damage state of materials, and it is an important link to establish the damage mechanics model of materials. It can be established from microscopic (number of voids, length, area, volume, etc.) and macroscopic (elastic constant, yield stress, tensile strength, etc.). Considering that the angle and length of some joints and fissures in rock mass can be observed in practical engineering, this paper describes the initial damage of rock mass from the geometric damage principle. e fracture geometry distribution of the rock sample is shown in Figure 8. e damage calculation method of the fractured rock sample is as follows: firstly, each fracture is projected on the plane perpendicular to x 1 , and the projection area S 30 is calculated; secondly, the distance L 30 between the center normal vectors of the two projection planes is calculated; thirdly, considering the influence of the fracture spacing on the mechanical properties of the rock sample, the spacing L 0 with the minimum influence of the fracture on the mechanical properties of the coal sample is calculated; and fourthly, the damage coefficient λ is introduced. e specific calculation formula is where D 30 is the damage value of 30°fractured coal sample; V 30 is the projected volume of crack in the x 1 direction; λ is the damage coefficient, λ � L 30 /L 0 ; and V is the volume of the rock sample. According to the nonlinear characteristics of rock creep stage, many researches on the rock creep model have been carried out in recent years. ere are four methods to establish the nonlinear creep model: ① According to the indoor or field creep test and measured data, the empirical model can be obtained by using the creep formula or empirical formula fitting regression, but its universality is poor ② e parameters of the linear element are treated unsteady to make it a nonlinear element, and then the modified nonlinear creep model is obtained ③ Based on the damage theory of creep and damage quantity, the damage model is established ④ e creep model is obtained considering the external environment, such as the nonlinear creep model coupled with temperature and chemical conditions Under the action of stress, the accumulation and evolution of internal crack damage is the essential reason for the final deformation and failure of rock. However, the fractured rock mass itself has some damage before the test, and then the creep test is carried out to further expand the damage. erefore, the creep damage of fractured rock mass does not only exist in the accelerated creep stage but also exist throughout the whole creep process. e research idea of this section is based on the time-dependent evolution mechanism of unloading creep damage of fractured rock mass, introducing the damage factor and combining with the existing creep element model, to establish a creep model which can describe the nonlinear deformation characteristics of each stage of unloading creep process of fractured rock mass, including attenuation creep, steady-state creep, and accelerated creep.
According to the Kachanov damage rate [30], the evolution equation of damage variable is as follows: where C and n are the material parameters. By integrating the above formula, we can get e fractured rock mass itself contains different degrees of initial cracks, and the rock mass contains initial damage during loading. erefore, the initial damage is not 0 at the moment of loading.
By introducing equation (4) into (3), we can get Shock and Vibration 7 By introducing equation (5) into (3), we can get e evolution law of damage variable of fractured rock mass in the whole creep process is as follows: e Burgers model is composed of the Maxwell model and the Kelvin model in series. e Kelvin body describes the viscoelastic properties of materials and reflects the attenuation creep of rocks; the Maxwell volume describes the elastic and viscous characteristics of the material, reflecting the instantaneous deformation of rock in the loading stage and the constant rate creep deformation in the steady-state creep stage. e Burgers model has the characteristics of elastic viscoelastic viscosity, which can better describe the attenuation and steady-state creep process of general rock mass under low stress level. In this paper, the Lade criterion considering plastic deformation is introduced into the Burgers model, and the damage variable D and switching element are introduced to construct the damage creep model which can describe the instantaneous elastic strain, viscoelastic strain, viscous strain, and nonlinear elastic-plastic strain of rock. e mechanical model is shown in Figure 9.

Damage Creep Model of Fractured Rock Mass.
According to the series parallel relationship of the model, we can know that where σ e , σ ve , σ v , and σ vp are the stress of the elastic body, Kelvin body, viscous body, and damaged elastoplastic body, respectively, and ε e , ε ve , ε v , and ε vp are the strain of the elastic body, Kelvin body, viscous body, and damaged elasticplastic body, respectively. In order to compare with the test results, the following reasonable assumptions are made firstly: ① e rock is isotropic, and the damage is assumed to be consistent in all directions ② e damage time and creep time are consistent in the accelerated creep stage e evolution law of the one-dimensional Burgers model based on damage variable in the whole creep process is as follows:

ree-Dimensional Creep Equation of Damage Creep Model of Fractured Rock Mass.
In the three-dimensional stress state, the stress tensor inside the material can be decomposed into the spherical stress tensor σ m and the deviatoric stress tensor S ij : where δ ij is the Kronecker function. It is generally believed that the spherical stress tensor can only change the volume of the material, but not its shape, while the partial stress tensor can only cause the shape change. Similarly, the strain tensor can be decomposed into spherical strain 8 Shock and Vibration tensor and deviatoric strain tensor, as shown in the following formula: For the Hooke body under three-dimensional stress state, where K is the bulk modulus and G is the shear modulus. Under the condition of isotropic material, it is assumed that the elastic strain is caused by the stress sphere tensor and the creep is caused by the stress deviator. e creep equation of the material under three-dimensional stress state is where K is the bulk modulus and G 1 and G 2 are the shear modulus. K, G 1 , G 2 , η 1 , and η 2 are all the model parameters, which are determined by fitting the test results: , where μ is Poisson's ratio. erefore, the parameters K, G 1 , and G 2 that need to be fitted can be transformed into E 1 and E 2 : e three-dimensional creep equation of the damaged plastic body is as follows: where function, F is the yield function of rock; F 0 is the initial value of the rock yield function, which is generally taken as 1 for the convenience of calculation; Φ(F/F 0 ) is the power function, Φ(F/F 0 ) � (F/F 0 ) m ; M is the experimental constant, generally take 1; and Q is the plastic potential function.
For the convenience of calculation, the relevant flow rule is adopted, F � Q. e Mohr-Coulomb yield criterion and the von Mises yield criterion are widely used. However, the Mohr-Coulomb criterion cannot consider the effect of intermediate principal stress on rock yield, and von Mises yield criterion ignores the influence of spherical stress on rock creep characteristics. In this paper, the yield criterion is adopted, which has the following advantages: (1) e results show that the three curves are smooth on the meridian plane, which is conducive to the numerical calculation and the determination of the direction of plastic strain increment (2) e S-D effect of different tensile and compressive strengths of geotechnical materials is considered (3) e material parameters are few and easy to be tested (4) It can be widely used in sandy soil, rock-soil, concrete, and over consolidated clay (5) It reflects the nonlinear relationship between the shear expansion yield curve and the hydrostatic pressure (6) e effect of hydrostatic pressure P on the yield of geotechnical materials and the nonlinear characteristics of yield and failure are considered e yield criterion is expressed as where I 1 and I 3 are the first and third invariants of the stress tensor; p a is the atmospheric pressure; and m and k are the stress level parameters. For conventional creep tests, Shock and Vibration e creep equation of the nonlinear creep model of fractured rock mass under three-dimensional stress state is obtained by using simultaneous equations (15)- (17):

Parameter Identification.
In the process of model parameter identification, the selection of the mathematical model identification method is very important. At present, the main methods to determine rock creep model parameters are as follows: least squares method, regression inversion method, and displacement back analysis method. e least squares fitting curve is widely used to determine the relevant parameters, but the disadvantage is that the unreasonable selection of the initial value will lead to divergence of results or slow convergence speed. In this paper, the parameters of the rock creep model are identified by the mathematical optimization analysis software 1stOpt. e software does not need the user to provide appropriate initial parameters, which overcomes the shortcomings of ordinary least squares method. Levenberg-Marquardt + general global optimization algorithm is adopted to find the optimal model parameters. Before fitting, equation (19) is complex, which is not conducive to the preparation of the user-defined function, so it is necessary to simplify the coefficient. e results are as follows: , σ ≤ σ s ,  where According to the fitting results, the coefficients A∼G are obtained, and then the model parameters are converted.

Model Validation.
In order to verify the effect of the creep damage model established in this paper on the whole creep stage, the traditional Burgers model and the F-D model considering the initial damage are used to fit the test curve at the same time.
e fitting results are shown in Table 3, and the comparison of fitting curves is shown in Figure 10.
Combined with Table 3 and Figure 10, it can be seen that both Burgers model and F-D model can accurately describe viscoelastic creep of fractured rock mass. In terms of fitting effect, adding damage elements to the Burgers model in series and introducing the Lade criterion can make the simulation effect of the improved model (F-D) more accurate and better reflect the creep deformation characteristics of fractured rock mass. Moreover, the correlation coefficient R 2 of model identification is more than 0.99, and the fitting accuracy of the model meets the requirements. e experimental curve and the theoretical curve are in good agreement; the theoretical curve of the model is basically consistent with the test curve.
e F-D model can not only accurately describe the creep characteristics of attenuation and constant velocity stage but also describe the accelerated creep characteristics of fractured rock mass. e deformation characteristics of three typical stages of rock rheology can be well described, which verifies the applicability of the model.

Conclusions
Rock burst is one of the major disasters in coal mines because of its great destructiveness. Because of the mining depth, creep-type rock burst occurs frequently, and the number of coal mines and coal seams is increasing. is kind of rock burst disaster has seriously restricted the safety production of coal mines in China. e prediction and prevention of this kind of rock burst is different from the conventional type. If the mining boundary is complex, the self-weight stress field is large, and the residual tectonic stress field exists, and the superposition of these multiple factors will make the coal (rock) body produce "strong" creep failure. If disturbed by the dynamic stress field of the working face, it is easy to produce large-scale and serious rock burst. e existing prevention methods of rock burst are to release the energy accumulated in coal and rock mass by releasing confining pressure, so as to reduce the stress concentration in coal and rock mass and prevent rock burst phenomenon. Based on this, triaxial unloading confining pressure creep tests are carried out on fractured coal samples with different angles, and the fracture damage creep model is established. e conclusions are as follows: (1) Under the same stress state, the different fracture dip angle significantly affects the creep deformation of coal samples. With the decrease of confining pressure, the axial instantaneous strain and creep strain increase with the increase of fracture dip angle, each curve is basically straight line, and the deformation increase amount of the same type coal sample under each stress level is basically the same. e order of axial instantaneous strain and creep strain from large to small are 60°fracture dip angle, 45°fracture dip angle, 30°fracture dip angle, and complete coal sample. Under the same stress state, with the decrease of confining pressure, the order of lateral instantaneous strain from large to small is 60°f racture angle, 45°fracture angle, 30°fracture angle, and complete coal sample; compared with complete coal sample, due to the existence of prefabricated cracks, the bearing capacity of fractured coal sample in creep process decreases, and the lateral expansion deformation increases. e larger the fracture angle, the greater the lateral strain. (2) Based on the geometric damage principle, a secondorder tensor is introduced to describe the initial damage of rock mass. According to Kachanov's creep damage theory, the damage variable with initial damage is introduced into the constitutive relation and creep equation, and the evolution equation of damage variable with time in the whole creep process is derived.

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

Conflicts of Interest
e authors declare that they have no conflicts of interest.