# Improved Nonlinear Nishihara Shear Creep Model with Variable Parameters for Rock-Like Materials

Authors: Hang Lin; Xing Zhang; Yixian Wang; Rui Yong; Xiang Fan; Shigui Du; Yanlin Zhao
Journal: Advances in Civil Engineering
 (2020)

---

## Abstract

Creep property is an important mechanical property of rocks. Given the complexity of rock masses, mechanical parameters change with time in the creep process. In this work, a nonlinear function for describing the time-dependent change of parameters was introduced and an improved variable-parameter nonlinear Nishihara shear creep model of rocks was established. By creating rock-like materials, the mechanical properties of rocks under the shear creep test condition were studied, and the deformation characteristics and long-term shear strength of rocks during creep were analyzed. The material parameters of the model were identified using the creep test results. Comparison of the model’s calculated values and experimental data indicated that the model can describe the creep characteristics of rocks well, thus proving the correctness and rationality of the improved model. During shear creep, the mechanical properties of rocks have an aging effect and show hardening characteristics under low shear stress. Furthermore, according to the fact thatGk of the nonlinear model can characterize the creep deformation resistance, a method to determine the long-term shear strength is proposed.


---

## Body

## 1. Introduction

Creep property is an important mechanical property of rocks and the main factor affecting long-term stability of geotechnical engineering [1–7]. Obtaining the stress and deformation characteristics of rock masses from the perspective and method of their rheology is consistent with reality [8–10]. Thus, studying the creep characteristics of rocks is crucial. The creep behavior of rocks can be divided into three stages, namely, attenuated, stable, and accelerated creep. In the initial creep stage, the strain increases with time. In the steady-state creep stage, the strain increases uniformly with time, and the duration of this stage is relatively long. In the accelerated creep stage, the strain rate increases with time, and the rock is destroyed after reaching the critical fracture point. Numerous studies have shown that rocks exhibit typical nonlinear characteristics during creep, in which the accelerated creep stage is the most obvious [4, 11, 12]. To reveal the mechanical properties of rock rheological behavior, scholars have conducted targeted in-depth research and obtained abundant research results [13–15]. By summarizing existing rheological models, Xia et al. [16] proposed a unified rheological mechanical model. Relevant studies have shown that model parameters vary during rock creep and are related to stress, time, temperature, and other factors [17, 18]. Referring to damage mechanics theory [19, 20], Jia et al. [21] proposed a nonlinear elastic-viscoplastic damage model based on the modified Mohr–Coulomb criterion. Moghadamaab [22] used -constitutive model to describe expansion and short- and long-term damages in the creep process of rock salt and established a numerical model of the rock salt excavation underground chamber as it changes with time. Yang and Li [11] replaced the viscous element of the Nishihara model with a new element by using fractional calculus and established a nonlinear Nishihara model that describes rock rheology. Shao et al. [23] established a full stress-strain creep model from the perspective of microstructure evolution by considering various factors. These models effectively describe the nonlinear characteristics of model parameters during rock creep but focus on the creep mechanical properties of rocks under uniaxial and triaxial compression. Research on rock characteristics under direct shear creep conditions is relatively scarce. The primary failure mode of rocks under natural conditions is shear failure [24–30]. Therefore, the shear creep failure behavior of rocks should be studied. The classical Nishihara model can describe the mechanical properties of rocks in decay and stable creep stages [11, 31, 32], but its parameters are constant and cannot reflect the nonlinear characteristics of creep. To reveal the mechanism of rock mechanical behavior under shear creep conditions and establish a corresponding shear creep model, this study investigated the shear rheological mechanical characteristics of rocks with plastic-elastic deformation by performing a shear test on a similar sandstone material. A time-varying factor was introduced to the Nishihara model to describe the functional relationship between the viscosity coefficient and time. On the basis of this factor, an improved nonlinear Nishihara shear creep model that can reveal the unsteady and nonlinear characteristics of parameters was established. The model’s rationality was verified through a comparison with experimental values.


## 2. Specimen Preparation and Test Methods

Natural rock masses contain many micro- and macrocracks, which lead to different physical and mechanical characteristics [33–41]. Miller divided the stress-strain curves of rocks with different physical and mechanical properties into six categories [42, 43]. The third category is called the plastic-elastic deformation feature, in which the stress-strain curve is concave in the initial stage and then becomes a straight-line segment until destruction [44–46]. Rocks with many pores and microfractures, such as granite and sandstone, reflect this deformation characteristic [47, 48]. In this study, the mechanical behavior and deformation characteristics under shear creep conditions were studied by fabricating similar materials with the plastic-elastic deformation feature of rocks [49].Experiments were performed on a RYL-600 computer-controlled rock shear rheometer, which is usually used for rheological testing of materials, such as rock and concrete, and composed of three components: the host, the test and control system, and the computer system. The rheometer can simultaneously apply vertical-axis compression test force and horizontal-axis shear test force to a specimen and control the loading by force or deformation. In the shear creep test, the rock-like specimens were under constant temperature and humidity conditions. The loading rate at each stage was 300 N/s. The specimens adopted the same type of material, mix ratio, curing time, and mold specification (10  cm × 10  cm × 10 cm). Through the uniaxial and triaxial compression tests, the mechanical parameters of rock materials with different cement, river sand, and water ratio were obtained. Finally, the sandstone materials with cement, river sand, and water ratio of 1 : 1 : 0.406 were selected as research objects. The mechanical parameters of the rock-like specimens are shown in Table1, which is similar to the mechanical properties of sandstones with elastoplastic deformation characteristics. During pouring, the material was mixed according to the ratio, and the surface was smoothed after tamping in the steel mold to minimize adverse effects in the manufacturing process. The mold was taken apart after modeling for 24 hours and maintained for 28 days. All specimens are intact rock-like materials.Table 1
The mechanical parameters of the rock-like specimens.
Parameter
σ
c (MPa)
μ
φ (°)
c (MPa)
Specimen
15.6
0.21
60
3.12As the loading model shown in Figure1, the specimens were subjected to a fast shear test to compare the strength characteristics of the rock’s fast shearing and shear creep conditions. The test results are shown in Figure 2. At the beginning of loading, the slope of the stress-strain curve increased with the loading stress, and this process lasted for a long time. After the microcracks were compacted completely, the slope of the curve entered the linear elastic stage with a constant size and then maintained this feature until failure. Significant plastic elastic deformation characteristics were exhibited. The fast shear strengths of specimens 1–4 were 4.1 MPa, 5.6 MPa, 6.45 MPa, and 7.33 MPa, respectively. In the grading loading creep test, normal stress was taken as 5%, 10%, 15%, and 20% of the uniaxial compressive strength of a specimen. The shear stress levels of specimen 1 were 30%, 40%, 50%, 60%, and 70% of the fast shear strength. The shear stress levels of specimen 2 were 30%, 40%, 50%, and 60% of the fast shear strength until destruction. The shear stress levels of specimen 3 were 30%, 40%, 60%, 70%, and 80% of the fast shear strength. The shear stress levels of specimen 4 were 30%, 50%, 60%, 70%, and 80% of the fast shear strength. The loading method is presented in Table 2, and the process is shown in Figure 3.Figure 1
(a) Specimen loading model and (b) test equipment diagram.
(a)
(b)Figure 2
Shear stress-shear displacement relationship in the fast shear test under 3.12 MPa normal stress.Table 2
Test grading loading method.
Specimen number
Normal stress (MPa)
Shear strength (MPa)
Shear stress (MPa)
1
2
3
4
5
1
0.78
4.1
1.23
1.64
2.05
2.46
2.87
2
1.56
5.6
1.68
2.24
2.8
3.36
Broken
3
2.34
6.45
1.935
2.58
3.87
4.515
5.16
4
3.12
7.33
2.199
3.665
4.398
5.131
5.864Figure 3
Schematic diagram of grading loading process.


## 3. Test Results and Analysis

### 3.1. Mechanical Properties with Time

The creep curves of graded shear loading under the same normal stress level were obtained through the creep shear test (Figure4). In the creep test curve obtained, a straight line parallel to the longitudinal axis was made to intersect with the creep curve under the same normal stress and different shear forces, and the isochronal curve (Figure 5) was drawn by the stress and displacement values at each intersection point obtained. As can be seen from Figure 5, the slope of the isochronous shear stress-shear displacement curve varies with stress. Stress is the factor that affects the mechanical properties of rocks with plastic-elastic deformation characteristics, which are different from the characteristics of the microfracture compaction stage completed immediately after loading in most rocks. To study the time deterioration effect on the rock creep process, the deformation modulus under different creep time conditions was calculated through the slope of the curve section of the adjacent shear stress levels, and the relationship between deformation modulus and time was plotted and is shown in Figure 6.Figure 4
Shear displacement-time curves under different normal stresses: (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)Figure 5
Isochronal graphs under different normal stresses: (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)Figure 6
Deformation modulus-time diagrams under different normal stresses: (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)The formula for calculating the deformation modulus can be written as(1)E=τ1−τ2γτ1t−γτ2t,where τ1 and τ2 characterize the adjacent two-stage shear stress and γτ1tandγτ2t indicate the shear displacement at different moments under the respective shear stress conditions.Figure6 indicates that for the rocks with plastic-elastic deformation characteristics, under the condition of low shear stress, the deformation modulus of rocks in creep behavior increased with time, which is contrary to the deterioration of the mechanical properties of most rocks with time [50, 51]. This result is due to the existence of many microcracks in plastic-elastic rocks [52, 53]. Low shear stress cannot easily compact the microcracks inside the rock completely in a short time, and the microcracks that are not compacted continue to be compacted in the following creep behavior, thus showing a hardening effect. When the shear stress level was increased to a certain value, the deformation modulus weakened with time and showed aging deterioration. This phenomenon occurred because the shear stress at this time exceeded the critical shear stress when the microcrack was completely compacted, which not only stopped the hardening effect brought by microcrack compaction but also developed the internal crack of the specimen with time and produced significant aging degradation effects. When the shear stress continues to increase, the mechanical properties of the rock continued to deteriorate with time because the unfinished area of the rock when shear stress was applied to the upper stage continued to deteriorate under the new first-order shear stress condition. Meanwhile, the increase in shear stress led to the deterioration of the area that is not easy to deteriorate under the condition of the shear stress of the upper level, resulting in continuous internal damage of the specimen.


### 3.2. Determination of Long-Term Strength under Different Normal Stress Conditions

According to the isochronous curve method which is commonly used in rheological mechanics [54], the stress corresponding to the second turning point of the isochronous curve after the fully compacted linear elastic phase was used as the long-term strength value of the rock. In Figures 5(a)–5(d), the shear stress value corresponding to the second turning point is 2.04 MPa, 2.8 MPa, 4.515 MPa, and 5.131 MPa, respectively. Thus, the long-term shear strength of the specimen under the normal stress condition of 0.78 MPa, 1.56 MPa, 2.34 MPa, and 3.12 MPa is 2.04 MPa, 2.8 MPa, 4.515 MPa, and 5.131 MPa, respectively. Correspondingly, the ratio of long-term shear strength to fast shear strength is 0.5, 0.5, 0.7, and 0.702, respectively. The phenomenon wherein the specimen had different ratios of long-term strength to fast shear strength under different normal stress conditions could be caused by the numerous microcracks in the rock with plastic-elastic deformation characteristics. Under low normal stress, the internal damage of the rock accumulated continuously overtime and with the increasing shear stress. Damage occurred when the local damage reached the limit and lost the bearing capacity. Under high normal stress, the microcracks in the normal direction of the rock were compacted closely, thereby making the hardening effect caused by the compaction of microcracks more obvious than that under low stress. Meanwhile, under the condition of dense microfractures, the mechanical properties and internal continuity of the rock improved, thus strengthening the capability to transfer stress inside the rock, increasing the effective bearing area, and strengthening the capability to resist damage. Therefore, the ratio of long-term shear strength to fast shear strength of the rock under high normal stress is greater than that under low normal stress.


## 3.1. Mechanical Properties with Time

The creep curves of graded shear loading under the same normal stress level were obtained through the creep shear test (Figure4). In the creep test curve obtained, a straight line parallel to the longitudinal axis was made to intersect with the creep curve under the same normal stress and different shear forces, and the isochronal curve (Figure 5) was drawn by the stress and displacement values at each intersection point obtained. As can be seen from Figure 5, the slope of the isochronous shear stress-shear displacement curve varies with stress. Stress is the factor that affects the mechanical properties of rocks with plastic-elastic deformation characteristics, which are different from the characteristics of the microfracture compaction stage completed immediately after loading in most rocks. To study the time deterioration effect on the rock creep process, the deformation modulus under different creep time conditions was calculated through the slope of the curve section of the adjacent shear stress levels, and the relationship between deformation modulus and time was plotted and is shown in Figure 6.Figure 4
Shear displacement-time curves under different normal stresses: (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)Figure 5
Isochronal graphs under different normal stresses: (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)Figure 6
Deformation modulus-time diagrams under different normal stresses: (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)The formula for calculating the deformation modulus can be written as(1)E=τ1−τ2γτ1t−γτ2t,where τ1 and τ2 characterize the adjacent two-stage shear stress and γτ1tandγτ2t indicate the shear displacement at different moments under the respective shear stress conditions.Figure6 indicates that for the rocks with plastic-elastic deformation characteristics, under the condition of low shear stress, the deformation modulus of rocks in creep behavior increased with time, which is contrary to the deterioration of the mechanical properties of most rocks with time [50, 51]. This result is due to the existence of many microcracks in plastic-elastic rocks [52, 53]. Low shear stress cannot easily compact the microcracks inside the rock completely in a short time, and the microcracks that are not compacted continue to be compacted in the following creep behavior, thus showing a hardening effect. When the shear stress level was increased to a certain value, the deformation modulus weakened with time and showed aging deterioration. This phenomenon occurred because the shear stress at this time exceeded the critical shear stress when the microcrack was completely compacted, which not only stopped the hardening effect brought by microcrack compaction but also developed the internal crack of the specimen with time and produced significant aging degradation effects. When the shear stress continues to increase, the mechanical properties of the rock continued to deteriorate with time because the unfinished area of the rock when shear stress was applied to the upper stage continued to deteriorate under the new first-order shear stress condition. Meanwhile, the increase in shear stress led to the deterioration of the area that is not easy to deteriorate under the condition of the shear stress of the upper level, resulting in continuous internal damage of the specimen.


## 3.2. Determination of Long-Term Strength under Different Normal Stress Conditions

According to the isochronous curve method which is commonly used in rheological mechanics [54], the stress corresponding to the second turning point of the isochronous curve after the fully compacted linear elastic phase was used as the long-term strength value of the rock. In Figures 5(a)–5(d), the shear stress value corresponding to the second turning point is 2.04 MPa, 2.8 MPa, 4.515 MPa, and 5.131 MPa, respectively. Thus, the long-term shear strength of the specimen under the normal stress condition of 0.78 MPa, 1.56 MPa, 2.34 MPa, and 3.12 MPa is 2.04 MPa, 2.8 MPa, 4.515 MPa, and 5.131 MPa, respectively. Correspondingly, the ratio of long-term shear strength to fast shear strength is 0.5, 0.5, 0.7, and 0.702, respectively. The phenomenon wherein the specimen had different ratios of long-term strength to fast shear strength under different normal stress conditions could be caused by the numerous microcracks in the rock with plastic-elastic deformation characteristics. Under low normal stress, the internal damage of the rock accumulated continuously overtime and with the increasing shear stress. Damage occurred when the local damage reached the limit and lost the bearing capacity. Under high normal stress, the microcracks in the normal direction of the rock were compacted closely, thereby making the hardening effect caused by the compaction of microcracks more obvious than that under low stress. Meanwhile, under the condition of dense microfractures, the mechanical properties and internal continuity of the rock improved, thus strengthening the capability to transfer stress inside the rock, increasing the effective bearing area, and strengthening the capability to resist damage. Therefore, the ratio of long-term shear strength to fast shear strength of the rock under high normal stress is greater than that under low normal stress.


## 4. Improved Nonlinear Shear Creep Model Based on the Nishihara Model

The test results above illustrate that during the process of creep from the decay creep stage to the stable creep stage, the internal microcracks are continuously compacted, which causes irreversible plastic deformation and a hardening effect. Therefore, for rocks with an obvious hardening effect under low shear stress conditions, the traditional Nishihara model with linear parameters is no longer applicable. The model parameters must be described in a nonlinear manner. To reflect the mechanical properties of the rock in the creep process accurately, a time-varying factor was introduced to improve the traditional Nishihara model, and an improved variable-parameter nonlinear Nishihara shear creep model of rock was established. The rationality of the model was verified by the creep test result of specimens.


### 4.1. Description of the Characteristics of the Time-Dependent Parameter in the Nishihara Model

The Nishihara model is a rheological model composed of Hooke, Kelvin, and viscoplastic bodies in series, as shown in Figure7. The model expression is(2)γ=τGM+τGK1−e−GK/ηKt,τ<τs,τGM+τGK1−e−GK/ηKt+τ−τsη2t,τ>τs,where GM is the Hooke body coefficient, which is the instantaneous elastic modulus; GK is the Kelvin model elastic coefficient; ηK is the Kelvin model viscosity coefficient; η2 is the viscous element coefficient in the viscoplastic body; τs is the long-term shear strength of the rock; and ts is the moment when the accelerated creep stage begins.Figure 7
Nishihara rheological model.To describe the nonlinearity of parameters in the Nishihara model under shear loading, time-varying factorD was introduced to the Nishihara model by referring to the method of establishing fractional-derivative models [55], to describe the variation in element coefficients in the Nishihara model with time, and a nonlinear shear creep model with variable parameters was constructed (Figure 8).Figure 8
Improved nonlinear shear creep model with variable parameters.A large number of studies have shown that the internal damage/hardening variable is closely related to time as an exponential function [56]. Therefore, on the basis of previous studies, the time-varying factor D is defined as(3)D=μe−αt,where D ∈ [0, 1]; μ and α are material parameters representing the variation in element coefficients μ ∈ [−∞, 1]; and α ∈ [0, +∞].Compared with the nonlinear function in the fractional-derivative model, which describes the aging characteristics of the model parameters with only one parameter [55, 57], in this paper, the function D used to describe the nonlinear characteristics of model parameters has two parameters, where μ represents the comparison of model parameters before and after creep occurs, and α reflects the trend of the model parameters over time. It is not as simple and easy to calculate as the function in fractional-derivate models, but it can describe the evolution characteristics of model parameters during the creep process more comprehensively.In the Nishihara model, the Kelvin body represents the attenuation and stable creep stages, and the viscoplastic body represents the accelerated creep stage. In the actual creep process, all parameters change with time. Time-varying factorD was introduced to the Kelvin and viscoplastic bodies to describe the relationship between the viscosity coefficient and time.The Kelvin model is composed of an elastic element and a viscous element in parallel, and its constitutive equation is(4)τ=τ1+τ2=GKγ+ηKγ˙,(5)γ=γ1=γ2=τ1GK=τ2ηKt,(6)γt=τGK1−e−GK/ηKt.After introducing time-varying factorD, the Kelvin viscous element with variable parameters is expressed as(7)ηKt=ηK∞1−D1=ηK∞1−μ1e−α1t,where ηKt is the actual viscosity coefficient that changes over time, when t = ∞, ηKt = ηK∞, ηK∞ is the viscosity coefficient of the final moment of rock stable creep; μ1 is the material parameter characterizing the damage/hardening degree of ηK; when μ1 is positive, ηK hardens with time during rock creep; when μ1 ⟶ 1, ηK hardens into a plastomer at t = ∞; when μ1 = 0, ηK is constant, which is the original Kelvin model; when μ1 is negative, ηK deteriorates with time during rock creep; when μ1 ⟶ ∞, ηK softens into an ideal fluid at t = ∞; α1 is a parameter reflecting the trend of ηK over time in the creep process; the larger α1 is, the longer it takes for the ηK to stabilize.Equation (7) indicates that when t = 0, the initial viscous element coefficient ηK0 is(8)ηK0=ηK∞1−D1=ηK∞1−μ1.In accordance with (3)–(7), the constitutive equation of nonlinear Kelvin model with parameter changes is obtained as follows:(9)τ=τ1+τ2=GK+ηKtγ˙,ηKtγ˙=ηK∞1−μ1e−α1tγ˙=τ−GKγ.Combined with the initial conditionst = 0 and γ = 0, the Kelvin nonlinear shear creep model equation with parameter variation is obtained as follows:(10)γkel=τGK1−e−GK/ηK∞lneα1t−μ1−ln1−μ1/α1.Similarly, for the viscoplastic body introduced with time-varying factorD, the viscous element with variable parameters is expressed as(11)η2t=η2∞1−D2=η2∞1−μ2e−α2t,where η2t is the actual viscosity coefficient that changes over time, when t = ∞, η2t = η2∞, η2∞ is the viscosity coefficient of the final moment of rock accelerated creep; when μ2 > 0, η2 hardens with time, and the larger μ2 represents the higher degree of hardening; when μ2 < 0, η2 softens with time, and the smaller μ2 represents the higher degree of softness.The constitutive equation of viscoplastic body of Nishihara model is(12)τ−τs=η2γ˙vp.Combining (11) and (12), the constitutive equation of the viscous element with variable parameters is obtained as follows:(13)τ−τs=η21−D2γ˙vp=η21−μ2e−α2t−tsγ˙vp,(14)γvp=τ−τsα2η2lneα2t−ts−μ2.Combining (10) and (14), the constitutive model of improved nonlinear Nishihara shear creep model with variable parameters can be obtained as follows:(15)γt=τGM+τGK1−e−GK/ηK∞lneα1t−μ1−ln1−μ1/α1,τ<τs,τGM+τGK1−e−GK/ηK∞lneα1t−μ1−ln1−μ1/α1+τ−τsα2η2∞lneα2t−ts−μ2,τ>τs.


### 4.2. Model Verification

To test the correctness of the improved nonlinear shear creep model established in this study, the model was used to fit the test curve with the experimental data of specimens, which was compared with the fitting curve of the original Nishihara model to determine the advantages and disadvantages of the improved model. The parameters of the improved model are given in Table3. The fitting results are shown in Figure 9, and it indicates that the improved Nishihara model was closer to the creep test curve than the original Nishihara model was, thereby reflecting the mechanical behavior of the specimen creep process and verifying the reasonableness of the model.Table 3
Model parameters under various shear stresses conditions.
Normal stress (MPa)
Shear stress (MPa)
G
M(MPa)
G
K(MPa)
η
K
∞(MPa·h)
η
2
∞(MPa·h)
μ
1
α
1
τ
s(MPa)
t
s (h)
μ
2
α
2
Correlation coefficientR2
0.78
1.23
0.8
11.60
44.866
—
0.985
3.861
—
—
—
—
0.998
1.64
0.93
24.35
174.511
—
0.977
0.601
—
—
—
—
0.999
2.05
1.01
22.36
137.981
—
0.996
1.096
—
—
—
—
0.997
2.46
1.06
21.59
105.654
—
0.993
1.324
—
—
—
—
0.996
2.87
1.08
19.67
196.583
—
0.991
0.703
—
—
—
—
0.999
1.56
1.68
1.06
17.58
260.98
—
0.999
0.435
—
—
—
—
0.999
2.24
1.57
116.67
1212.06
—
0.977
0.628
—
—
—
—
0.999
2.8
1.76
39.27
982.89
—
0.971
0.390
—
—
—
—
0.999
3.36
0.74
7.10
3.9693
0.087
0.970
20.32
2.8
0.2634
−14.70
3.64
0.996
2.34
1.935
1.16
17.06
37.564
—
0.984
3.868
—
—
—
—
0.998
2.58
1.26
32.09
1223.15
—
0.983
0.145
—
—
—
—
0.997
3.87
1.48
42.2
603.86
—
0.990
0.787
—
—
—
—
0.999
4.515
1.55
17.94
291.13
—
0.986
0.461
—
—
—
—
0.996
5.16
1.57
9.52
245.49
—
0.977
0.192
—
—
—
—
0.998
3.12
2.199
0.54
22.15
212.33
—
0.993
1.894
—
—
—
—
0.999
2.99
0.65
22.46
370.49
—
0.979
0.246
—
—
—
—
0.998
4.398
0.86
25.86
287.12
—
0.997
0.769
—
—
—
—
0.997
5.131
0.97
21.73
186.73
—
0.999
1.678
—
—
—
—
0.999
5.864
1.06
21.76
290.39
—
0.999
0.821
—
—
—
—
0.998Figure 9
Comparison of the test and calculated values of the creep model under normal stress (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)


### 4.3. Relationship between Creep Model Parameters and Stress

In rock rheology, constitutive model parameters have important physical and mechanical significance. Studying the relationship between model parameters and loading conditions (stress, loading rate, etc.) can effectively reveal the mechanical mechanism of rock rheology [58–60]. In order to better understand the mechanical properties characterized by D, on the basis of Table 3, the relationship between μ1, α1, and shear stress under the normal stress 1.56 MPa is plotted in Figure 10. It could be seen from Figure 10 that μ1 is greater than 0 under each shear stress condition, and the maximum value occurs when the shear stress is 1.68 MPa and then decreases with the increase of shear stress and finally stabilizes to 0.97, indicating that the nonlinear Kelvin viscosity coefficient has hardening under various shear stress conditions, and the hardening degree is the highest under low shear stress condition, and the hardening degree tends to be stable after shear stress increases. As shown in Figure 10, as the shear stress increases, α1 increases first, then decays, and then continues to increase. Combined with (7), it can be seen that the time required for the viscosity coefficient in the nonlinear Kelvin body to reach stability is relatively long under the condition of low shear stress and then fluctuates with the increase of shear stress, but all remain stable at a certain value. Finally, when the shear stress reaches the value beyond the long-term strength, it decreases rapidly, and this phenomenon is caused by the compaction of rock-like specimen under low shear stress conditions. The relationship between ηKt and time in accordance with (7) and Table 4 is plotted in Figure 11, which can better show the change process of ηKt represented by μ1 and α1.Figure 10
The relationship betweenμ1, α1, and shear stress diagram.Table 4
Initial viscosity coefficientη0 under different shear stress conditions.
Shear stress (MPa)
1.68
2.24
2.8
3.36
η
0 (MPa·h)
23.56
67.65
46.84
0.098Figure 11
The relationship betweenηKt and time diagram.
(a)
(b)The relationship between Kelvin model elastic coefficientGK and shear stress is plotted in Figure 12(b), and it shows that when the range of the shear stress was from 1.68 MPa to 2.24 MPa, Kelvin model elastic coefficient GK increased with stress, and this result indicates the obvious internal hardening effect of the specimen and the characteristics of the plastic-elastic deformation of the rock-like specimen that not only hardens under low-stress fast shear but also under low-stress creep. When the shear stress was loaded to 2.8 MPa, the GK was attenuated, indicating that under this shear stress condition, the internal microcracks of the specimen began to develop and damage began to occur. After the shear stress reached 3.36 MPa, the specimen quickly entered the accelerated creep stage. The viscosity parameter of the viscous element in the viscoplastic body that describes the accelerated creep stage was negative, indicating that the internal damage of the specimen was continuously intensified and the creep rate was continuously accelerated until failure.Figure 12
The relationship between viscosity modulusGK-shear stress diagram under normal stress (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)


### 4.4. Long-Term Strength Determination Method Based on Improved Nonlinear Model

It can be concluded from (10) that the rock-like specimen creep deformation is controlled by GK; when GK is small, the mass of rock creep deformation is large; that is, the larger GK is, the smaller the creep deformation is. GK characterizes the ability to resist creep deformation. Therefore, in the improved nonlinear model, the author proposes a method to determine the long-term strength by using the curve of GK-shear stress. This method uses the nonlinear Nishihara model proposed in this paper to fit shear displacement-time curve to obtain the curve of GK−τ and then determines the shear stress of GK decay in the curve as the long-term strength. According to the parameters in Table 3, the GK−τ curves of the specimen under different normal stress levels were drawn (Figure 12) and were compared with the long-term shear strength determined by the isochronous curve method, as shown in Table 5.Table 5
Comparison of two methods for determining long-term shear strength.
Specimen number
Normal stress (MPa)
Long-term shear stress (MPa)
Isochronous curve method
Modulus inflection point method
1
0.78
2.04
2.04
2
1.56
2.8
2.8
3
2.34
4.515
4.515
4
3.12
5.131
5.131Table5 shows that the long-term shear strength determined by these two methods is basically the same, which proves the effectiveness of the long-term strength determination method proposed in this paper. And compared with the isochronous curve method, the modulus inflection point method only needs to fit the experimental data to determine the GK under each stress condition. Then, the long-term strength of the rock can be obtained according to the attenuation of GK, which has a simpler and faster operation process.


## 4.1. Description of the Characteristics of the Time-Dependent Parameter in the Nishihara Model

The Nishihara model is a rheological model composed of Hooke, Kelvin, and viscoplastic bodies in series, as shown in Figure7. The model expression is(2)γ=τGM+τGK1−e−GK/ηKt,τ<τs,τGM+τGK1−e−GK/ηKt+τ−τsη2t,τ>τs,where GM is the Hooke body coefficient, which is the instantaneous elastic modulus; GK is the Kelvin model elastic coefficient; ηK is the Kelvin model viscosity coefficient; η2 is the viscous element coefficient in the viscoplastic body; τs is the long-term shear strength of the rock; and ts is the moment when the accelerated creep stage begins.Figure 7
Nishihara rheological model.To describe the nonlinearity of parameters in the Nishihara model under shear loading, time-varying factorD was introduced to the Nishihara model by referring to the method of establishing fractional-derivative models [55], to describe the variation in element coefficients in the Nishihara model with time, and a nonlinear shear creep model with variable parameters was constructed (Figure 8).Figure 8
Improved nonlinear shear creep model with variable parameters.A large number of studies have shown that the internal damage/hardening variable is closely related to time as an exponential function [56]. Therefore, on the basis of previous studies, the time-varying factor D is defined as(3)D=μe−αt,where D ∈ [0, 1]; μ and α are material parameters representing the variation in element coefficients μ ∈ [−∞, 1]; and α ∈ [0, +∞].Compared with the nonlinear function in the fractional-derivative model, which describes the aging characteristics of the model parameters with only one parameter [55, 57], in this paper, the function D used to describe the nonlinear characteristics of model parameters has two parameters, where μ represents the comparison of model parameters before and after creep occurs, and α reflects the trend of the model parameters over time. It is not as simple and easy to calculate as the function in fractional-derivate models, but it can describe the evolution characteristics of model parameters during the creep process more comprehensively.In the Nishihara model, the Kelvin body represents the attenuation and stable creep stages, and the viscoplastic body represents the accelerated creep stage. In the actual creep process, all parameters change with time. Time-varying factorD was introduced to the Kelvin and viscoplastic bodies to describe the relationship between the viscosity coefficient and time.The Kelvin model is composed of an elastic element and a viscous element in parallel, and its constitutive equation is(4)τ=τ1+τ2=GKγ+ηKγ˙,(5)γ=γ1=γ2=τ1GK=τ2ηKt,(6)γt=τGK1−e−GK/ηKt.After introducing time-varying factorD, the Kelvin viscous element with variable parameters is expressed as(7)ηKt=ηK∞1−D1=ηK∞1−μ1e−α1t,where ηKt is the actual viscosity coefficient that changes over time, when t = ∞, ηKt = ηK∞, ηK∞ is the viscosity coefficient of the final moment of rock stable creep; μ1 is the material parameter characterizing the damage/hardening degree of ηK; when μ1 is positive, ηK hardens with time during rock creep; when μ1 ⟶ 1, ηK hardens into a plastomer at t = ∞; when μ1 = 0, ηK is constant, which is the original Kelvin model; when μ1 is negative, ηK deteriorates with time during rock creep; when μ1 ⟶ ∞, ηK softens into an ideal fluid at t = ∞; α1 is a parameter reflecting the trend of ηK over time in the creep process; the larger α1 is, the longer it takes for the ηK to stabilize.Equation (7) indicates that when t = 0, the initial viscous element coefficient ηK0 is(8)ηK0=ηK∞1−D1=ηK∞1−μ1.In accordance with (3)–(7), the constitutive equation of nonlinear Kelvin model with parameter changes is obtained as follows:(9)τ=τ1+τ2=GK+ηKtγ˙,ηKtγ˙=ηK∞1−μ1e−α1tγ˙=τ−GKγ.Combined with the initial conditionst = 0 and γ = 0, the Kelvin nonlinear shear creep model equation with parameter variation is obtained as follows:(10)γkel=τGK1−e−GK/ηK∞lneα1t−μ1−ln1−μ1/α1.Similarly, for the viscoplastic body introduced with time-varying factorD, the viscous element with variable parameters is expressed as(11)η2t=η2∞1−D2=η2∞1−μ2e−α2t,where η2t is the actual viscosity coefficient that changes over time, when t = ∞, η2t = η2∞, η2∞ is the viscosity coefficient of the final moment of rock accelerated creep; when μ2 > 0, η2 hardens with time, and the larger μ2 represents the higher degree of hardening; when μ2 < 0, η2 softens with time, and the smaller μ2 represents the higher degree of softness.The constitutive equation of viscoplastic body of Nishihara model is(12)τ−τs=η2γ˙vp.Combining (11) and (12), the constitutive equation of the viscous element with variable parameters is obtained as follows:(13)τ−τs=η21−D2γ˙vp=η21−μ2e−α2t−tsγ˙vp,(14)γvp=τ−τsα2η2lneα2t−ts−μ2.Combining (10) and (14), the constitutive model of improved nonlinear Nishihara shear creep model with variable parameters can be obtained as follows:(15)γt=τGM+τGK1−e−GK/ηK∞lneα1t−μ1−ln1−μ1/α1,τ<τs,τGM+τGK1−e−GK/ηK∞lneα1t−μ1−ln1−μ1/α1+τ−τsα2η2∞lneα2t−ts−μ2,τ>τs.


## 4.2. Model Verification

To test the correctness of the improved nonlinear shear creep model established in this study, the model was used to fit the test curve with the experimental data of specimens, which was compared with the fitting curve of the original Nishihara model to determine the advantages and disadvantages of the improved model. The parameters of the improved model are given in Table3. The fitting results are shown in Figure 9, and it indicates that the improved Nishihara model was closer to the creep test curve than the original Nishihara model was, thereby reflecting the mechanical behavior of the specimen creep process and verifying the reasonableness of the model.Table 3
Model parameters under various shear stresses conditions.
Normal stress (MPa)
Shear stress (MPa)
G
M(MPa)
G
K(MPa)
η
K
∞(MPa·h)
η
2
∞(MPa·h)
μ
1
α
1
τ
s(MPa)
t
s (h)
μ
2
α
2
Correlation coefficientR2
0.78
1.23
0.8
11.60
44.866
—
0.985
3.861
—
—
—
—
0.998
1.64
0.93
24.35
174.511
—
0.977
0.601
—
—
—
—
0.999
2.05
1.01
22.36
137.981
—
0.996
1.096
—
—
—
—
0.997
2.46
1.06
21.59
105.654
—
0.993
1.324
—
—
—
—
0.996
2.87
1.08
19.67
196.583
—
0.991
0.703
—
—
—
—
0.999
1.56
1.68
1.06
17.58
260.98
—
0.999
0.435
—
—
—
—
0.999
2.24
1.57
116.67
1212.06
—
0.977
0.628
—
—
—
—
0.999
2.8
1.76
39.27
982.89
—
0.971
0.390
—
—
—
—
0.999
3.36
0.74
7.10
3.9693
0.087
0.970
20.32
2.8
0.2634
−14.70
3.64
0.996
2.34
1.935
1.16
17.06
37.564
—
0.984
3.868
—
—
—
—
0.998
2.58
1.26
32.09
1223.15
—
0.983
0.145
—
—
—
—
0.997
3.87
1.48
42.2
603.86
—
0.990
0.787
—
—
—
—
0.999
4.515
1.55
17.94
291.13
—
0.986
0.461
—
—
—
—
0.996
5.16
1.57
9.52
245.49
—
0.977
0.192
—
—
—
—
0.998
3.12
2.199
0.54
22.15
212.33
—
0.993
1.894
—
—
—
—
0.999
2.99
0.65
22.46
370.49
—
0.979
0.246
—
—
—
—
0.998
4.398
0.86
25.86
287.12
—
0.997
0.769
—
—
—
—
0.997
5.131
0.97
21.73
186.73
—
0.999
1.678
—
—
—
—
0.999
5.864
1.06
21.76
290.39
—
0.999
0.821
—
—
—
—
0.998Figure 9
Comparison of the test and calculated values of the creep model under normal stress (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)


## 4.3. Relationship between Creep Model Parameters and Stress

In rock rheology, constitutive model parameters have important physical and mechanical significance. Studying the relationship between model parameters and loading conditions (stress, loading rate, etc.) can effectively reveal the mechanical mechanism of rock rheology [58–60]. In order to better understand the mechanical properties characterized by D, on the basis of Table 3, the relationship between μ1, α1, and shear stress under the normal stress 1.56 MPa is plotted in Figure 10. It could be seen from Figure 10 that μ1 is greater than 0 under each shear stress condition, and the maximum value occurs when the shear stress is 1.68 MPa and then decreases with the increase of shear stress and finally stabilizes to 0.97, indicating that the nonlinear Kelvin viscosity coefficient has hardening under various shear stress conditions, and the hardening degree is the highest under low shear stress condition, and the hardening degree tends to be stable after shear stress increases. As shown in Figure 10, as the shear stress increases, α1 increases first, then decays, and then continues to increase. Combined with (7), it can be seen that the time required for the viscosity coefficient in the nonlinear Kelvin body to reach stability is relatively long under the condition of low shear stress and then fluctuates with the increase of shear stress, but all remain stable at a certain value. Finally, when the shear stress reaches the value beyond the long-term strength, it decreases rapidly, and this phenomenon is caused by the compaction of rock-like specimen under low shear stress conditions. The relationship between ηKt and time in accordance with (7) and Table 4 is plotted in Figure 11, which can better show the change process of ηKt represented by μ1 and α1.Figure 10
The relationship betweenμ1, α1, and shear stress diagram.Table 4
Initial viscosity coefficientη0 under different shear stress conditions.
Shear stress (MPa)
1.68
2.24
2.8
3.36
η
0 (MPa·h)
23.56
67.65
46.84
0.098Figure 11
The relationship betweenηKt and time diagram.
(a)
(b)The relationship between Kelvin model elastic coefficientGK and shear stress is plotted in Figure 12(b), and it shows that when the range of the shear stress was from 1.68 MPa to 2.24 MPa, Kelvin model elastic coefficient GK increased with stress, and this result indicates the obvious internal hardening effect of the specimen and the characteristics of the plastic-elastic deformation of the rock-like specimen that not only hardens under low-stress fast shear but also under low-stress creep. When the shear stress was loaded to 2.8 MPa, the GK was attenuated, indicating that under this shear stress condition, the internal microcracks of the specimen began to develop and damage began to occur. After the shear stress reached 3.36 MPa, the specimen quickly entered the accelerated creep stage. The viscosity parameter of the viscous element in the viscoplastic body that describes the accelerated creep stage was negative, indicating that the internal damage of the specimen was continuously intensified and the creep rate was continuously accelerated until failure.Figure 12
The relationship between viscosity modulusGK-shear stress diagram under normal stress (a) 0.78 MPa, (b) 1.56 MPa, (c) 2.34 MPa, and (d) 3.12 MPa.
(a)
(b)
(c)
(d)


## 4.4. Long-Term Strength Determination Method Based on Improved Nonlinear Model

It can be concluded from (10) that the rock-like specimen creep deformation is controlled by GK; when GK is small, the mass of rock creep deformation is large; that is, the larger GK is, the smaller the creep deformation is. GK characterizes the ability to resist creep deformation. Therefore, in the improved nonlinear model, the author proposes a method to determine the long-term strength by using the curve of GK-shear stress. This method uses the nonlinear Nishihara model proposed in this paper to fit shear displacement-time curve to obtain the curve of GK−τ and then determines the shear stress of GK decay in the curve as the long-term strength. According to the parameters in Table 3, the GK−τ curves of the specimen under different normal stress levels were drawn (Figure 12) and were compared with the long-term shear strength determined by the isochronous curve method, as shown in Table 5.Table 5
Comparison of two methods for determining long-term shear strength.
Specimen number
Normal stress (MPa)
Long-term shear stress (MPa)
Isochronous curve method
Modulus inflection point method
1
0.78
2.04
2.04
2
1.56
2.8
2.8
3
2.34
4.515
4.515
4
3.12
5.131
5.131Table5 shows that the long-term shear strength determined by these two methods is basically the same, which proves the effectiveness of the long-term strength determination method proposed in this paper. And compared with the isochronous curve method, the modulus inflection point method only needs to fit the experimental data to determine the GK under each stress condition. Then, the long-term strength of the rock can be obtained according to the attenuation of GK, which has a simpler and faster operation process.


## 5. Discussion and Conclusion

Research on the mechanical characteristics of plastic-elastic rock deformation under the conditions of fast shear and shear creep indicates that under low shear stress, the mechanical properties of rocks with plastic-elastic deformation are continuously strengthened during creep, thereby continuously enhancing the rock’s resistance to deformation. When the shear stress level increases, the microcracks inside the rock are completely compacted. At this time, the mechanical properties of the rock deteriorate with creep and aging. Meanwhile, the mechanical properties of rocks with plastic-elastic deformation under different normal stress conditions differ. The rock bearing capacity and the damage resistance under high normal stress are better than those under low normal stress. Therefore, for rock mass engineering, the rock with plastic-elastic deformation characteristics should be placed under low stress conditions, and the surrounding area should avoid high-intensity stress field to enhance the stability of the engineering and the project.In this work, a time-varying factorD was introduced to describe the viscous elements in a nonlinear manner, and an improved nonlinear shear creep model based on the Nishihara model was established. After fitting the experimental data, it is found that the model can reflect the shear creep behavior of rock-like specimens effectively, and the model has rationality. In addition, under the normal stress of 0.78 MPa and 1.56 MPa, the ratio of long-term shear strength to fast shear strength of rock samples with plastic-elastic deformation characteristics is 0.5, indicating that the critical damage cumulative value required for rock reaches a failure state. Therefore, in a project, the long-term shear strength at which the cracks inside the rock begin to develop and expand must be determined to prevent rock destruction. At present, relevant scholars often monitor rock damage under creep conditions by calculating the creep rate, the deformation modulus, or the monitoring energy of the acoustic emission equipment. However, in practical engineering, calculating the rate at each moment or providing special equipment is difficult. By fitting the improved nonlinear Nishihara shear creep model to the surrounding rock displacement-time data monitored in the project, the damage can be discriminated according to the identified trend of the viscous modulus Gk of the surrounding rock so as to provide timely support and avoid the structural instability of the project, which is crucial in practical engineering. Furthermore, according to the fact that Gk of the nonlinear model can characterize the creep deformation resistance of rock-like specimen, a method to determine the long-term shear strength is proposed, and compared with the calculation results of isochronous curve method, the rationality and applicability of the method are verified.



---

*Source: 7302141-2020-02-27.xml*