# A Strain Hardening and Softening Constitutive Model for Hard Brittle Rocks
## 4 Establishment and Verification of the Constitutive Model  
### Establishment of the Constitutive Model  
In FLAC\({}^{3\text{D}}\), the strain hardening and softening constitutive model considering confining pressures is coded by the C++ language. The new stress vector \(\sigma^{N}\) of the current step can be calculated from the stress vector \(\sigma^{p}\) of the previous step and the strain increment vector \(\mathbf{\Delta\varepsilon}\) of the current step. First, the elastic guess stress \(\sigma_{i}^{I}\) is defined as follows:  
\[\sigma_{i}^{I}=\sigma_{i}^{p}+S_{i}\Delta\varepsilon_{i},i=1,2,3 \tag{26}\]  
\[S_{i}= \begin{bmatrix}K+\frac{4}{3}G&K-\frac{2}{3}G&K-\frac{2}{3}G\\ K-\frac{2}{3}G&K+\frac{3}{3}G&K-\frac{2}{3}G\\ K-\frac{2}{3}G&K-\frac{2}{3}G&K+\frac{4}{3}G\end{bmatrix}  \tag{27}\]  
where \(S_{i}\) is the elastic stiffness matrix, \(K\) is the bulk modulus, and \(G\) is the shear modulus. Figure 7: Comparison of the overall fit using the cohesion and internal friction angle and the individual fit of two hard brittle rocks at different confining pressures: (**a**,**b**) are for Beishan granite [26], and (**c**,**d**) are for Danjiangkou diabase [37]. \(\Delta\varepsilon_{i}\) is the total principal strain increment, which can be calculated as the sum of the principal elastic strain increment \(\Delta\varepsilon_{i}^{e}\) and the principal plastic strain increment \(\Delta\varepsilon_{i}^{ps}\), as follows:  
\[\Delta\varepsilon_{i}=\Delta\varepsilon_{i}^{e}+\Delta\varepsilon_{i}^{ps} \tag{28}\]  
For the elastic part, the principal elastic strain increment \(\Delta\varepsilon_{i}^{e}\) and the principal stress increment \(\Delta\sigma_{i}\) follow [PERSON]'s law, as follows:  
\[\Delta\sigma_{i}=S_{i}(\Delta\varepsilon_{i}^{e}),i=1,2,3 \tag{29}\]  
Second, judge whether the simulation element is yielded. The maximum principal guess stress \(\sigma_{1}^{I}\) and the minimum principal guess stress \(\sigma_{3}^{N}\) are substituted into the yield criterion to judge whether the element is yielded. In FLAC\({}^{3\text{D}}\), the Mohr-Coulomb yield criterion is given as follows:  
\[f^{s}=\sigma_{1}-\sigma_{3}\frac{1+\sin\varphi}{1-\sin\varphi}+2c\sqrt{\frac{1 +\sin\varphi}{1-\sin\varphi}} \tag{30}\]  
Third, if the element is yielded, then its stress should be updated. Substituting Equations (28) and (29) into Equation (8), the principal stress increment \(\Delta\sigma_{i}\) would be recalculated as follows:  
\[\Delta\sigma_{i}=S_{i}(\Delta\varepsilon_{i}^{e})=S_{i}(\Delta\varepsilon_{i })-\lambda^{s}S_{i}\bigg{(}\frac{\partial g^{s}}{\partial\sigma_{i}}\bigg{)},i=1,2,3 \tag{31}\]  
The principal stress \(\sigma_{i}^{N}\) of the current step can be expressed as follows:  
\[\sigma_{i}^{N}=\sigma_{i}^{p}+\Delta\sigma_{i},i=1,2,3 \tag{32}\]  
Substituting Equations (26) and (31) into Equation (32), the principal stress \(\sigma_{i}^{N}\) of the current step can be written as follows:  
\[\sigma_{i}^{N}=\sigma_{i}^{I}-\lambda^{s}S_{i}\bigg{(}\frac{\partial g^{s}}{ \partial\sigma_{i}}\bigg{)},i=1,2,3 \tag{33}\]  
In FLAC\({}^{3\text{D}}\), when the element yields, the position of the principal stress point of the current step is on the yield surface, which can be expressed as follows:  
\[f^{s}\Big{(}\sigma_{1}^{N},\sigma_{3}^{N}\Big{)}=0 \tag{34}\]  
where \(\sigma_{i}^{N}\) and \(\sigma_{3}^{N}\) are the maximum and minimum principal stresses of the current step, respectively. Substituting Equations (8), (26)-(30) and (33) into Equation (34), a formula for plastic multiplier \(\lambda^{s}\) can be written as follows:  
\[\lambda^{s}=\frac{f^{s}\big{(}\sigma_{1}^{I},\sigma_{3}^{I}\big{)}}{\big{(} \alpha_{1}-\alpha_{2}N_{\psi}\big{)}-\big{(}-\alpha_{1}N_{\psi}+\alpha_{2} \big{)}N_{\psi}} \tag{35}\]  
\[\alpha_{1}=K+\frac{4}{3}G \tag{36}\]  
\[\Delta\varepsilon_{i}=\Delta\varepsilon_{i}^{e}+\Delta\varepsilon_{i}^{ps} \tag{28}\]0  
According to Equation (12), the maximum principal plastic strain increment \(\Delta\varepsilon_{1}^{ps}\) can be calculated.