## 3 Numerical Modeling Calculation
### Creep Constitutive Law and Its Parameters  
The widely used classic Burgers creep model primarily accounts for the creep behavior of rocks before the third phase. However, it only considers rock viscoelasticity and does not adequately capture the creep characteristics of soft rocks, which typically exhibit instantaneous plasticity, elasticity, viscoplasticity, and viscoelasticity. To address this limitation, we propose an improved Burgers creep model that combines the Mohr-Coulomb criterion with the classical Burgers creep correction model. By incorporating the Mohr-Coulomb criterion, we can better simulate the complex viscoelastic-plastic properties exhibited by various rock samples. This improved model allows us to accurately determine the creep characteristics of each sample. To establish this new creep model, we conducted creep tests on different rock samples, enabling us to simulate their viscoelastic-plastic behavior. This approach provides a more precise representation of the creep properties exhibited by soft rocks and enhances our understanding of their deformation mechanisms [21; 22; 23]. The classic Burgers creep model is composed of the Kelvin model and the Maxwell model in series, as shown in Figure 3.  
Under constant axial stress \(\sigma_{0}\), the axial strain \(\varepsilon(t)\) is:  
\[\varepsilon(t)=\frac{\sigma_{0}}{9K}+\frac{\sigma_{0}}{3G_{1}}+\frac{\sigma_{ 0}}{3G_{2}}+\frac{\sigma_{0}}{3G_{2}}e^{-E_{2}/\eta_{2}t}+\frac{\sigma_{0}}{3 \eta_{1}}t \tag{1}\]  
In Equation (1), \(K\) represents the bulk modulus of the rock sample, \(G_{1}\) represents the elastic shear modulus of the model, and \(G_{2}\) represents the modulus of control delay elasticity. \(\eta_{1}\) and \(\eta_{2}\) determine the rate of viscous flow and delayed elasticity in the model, respectively. You are correct that weak rocks often exhibit simultaneous properties of instantaneous plasticity, instantaneous elasticity, viscoplasticity, and viscoelasticity. The classical Burgers \begin{table}
\begin{tabular}{c c c c c c} \hline \hline \multirow{2}{*}{**Material Science**} & **Deformation** & **Poisson’s** & \multirow{2}{*}{**Cohesion/MPa**} & **Internal Friction** & \multirow{2}{*}{**Density/kg/m3**} \\  & **Modulus/GPa** & & & & **Angle/\({}^{\circ}\)** & \\ \hline V-grade surrounding rock & 18 & 0.30 & 0.05 & 25 & 2000 \\ Initial expenditure & 1.6 & 0.4 & / & / & 2500 \\ Secondary support & 20 & 0.25 & / & / & 2500 \\ IV-grade surrounding rock & 30 & 0.25 & 0.35 & 30 & 2200 \\ III-grade surrounding rock & 40 & 0.22 & 0.50 & 35 & 2500 \\ \hline \hline \end{tabular}
\end{table} 
Table 1: Basic parameters of simulated materials. Figure 3: Typical Burgers creep model and characteristic curve. creep model has limitations in describing the entire process of rock and soil creep as it only focuses on the viscoelastic behavior of materials before the third stage of creep. To overcome these limitations, a new plastic element called the M-C (Mohr-Coulomb) element is introduced into the Burgers creep model as the basic model. The M-C element is based on the Mohr-Coulomb criterion and provides a representation of the material's behavior. When the stress is below the yield stress determined by the Mohr-Coulomb criterion, the strain of the element remains at 0. However, when the stress exceeds or equals the yield stress, it follows the Mohr-Coulomb plastic flow law. By connecting this M-C element and the viscoelastic Burgers model in series, the modified Burgers creep model can be constructed. This combination allows for the simulation of elastic-plastic volume behavior and captures viscoelastic deviant characteristics. It assumes that the deformation of the viscoelastic and viscoplastic strain rate components is coordinated. To further enhance the model, a Kelvin model and a Maxwell model are connected in series to form a viscoelastic body model.