# Numerical Implementation of the Barcelona Basic Model Based on Return-Mapping Integration
## 4 Results
In the stress-strain curve, the mechanical response of strain increment can be expressed as  
\[\Delta\epsilon=\epsilon_{n+1}^{k}-\epsilon_{n}\] (A2)  
Then, in the new iterative cycle step, the normal tangent modulus can be obtained from the partial derivative of the new stress tensor with respect to the new strain tensor  
\[\mathbf{c}_{n+1}^{k}=\frac{\partial\mathbf{c}_{n+1}^{k}}{\partial\epsilon_{n+1 }^{k}}\] (A3)  
Substitute Equation (A1) into Equation (A3), we get  
\[\mathbf{c}_{n+1}^{k}=\frac{\partial\sigma_{n+1}^{k}}{\partial\epsilon_{n+1}^{ k}}=\mathbf{I}\otimes\frac{\partial p}{\partial\epsilon_{n+1}^{k}}+\sqrt{\frac{2}{3}}q \frac{\partial\mathbf{\hat{n}}}{\partial\epsilon_{n+1}^{k}}+\sqrt{\frac{2}{3} }\mathbf{\hat{n}}\otimes\frac{\partial q}{\partial\epsilon_{n+1}^{k}}\] (A4)  
Corresponding to the partial variables in the above equation, where  
\[\frac{\partial p}{\partial\epsilon_{n+1}^{k}}=\frac{\partial p_{n+1}^{\pi}}{ \partial\epsilon_{n+1}^{k}}-KM^{2}\Bigg{[}\bigg{(}2\frac{\partial p}{ \partial\epsilon_{n+1}^{k}}-\frac{\partial p_{0}}{\partial\epsilon_{n+1}^{k} }\bigg{)}\Delta\phi+(2p+p_{s}-p_{0})\frac{\partial\Delta\phi}{\partial \epsilon_{n+1}^{k}}\Bigg{]}\] (A5)  
\[\frac{\partial p_{0}}{\partial\epsilon_{n+1}^{k}}=p_{0}\theta(s)M^{2}\Bigg{[} \bigg{(}2\frac{\partial p}{\partial\epsilon_{n+1}^{k}}-\frac{\partial p_{0}}{ \partial\epsilon_{n+1}^{k}}\bigg{)}\Delta\phi+(2P+P_{s}-P_{0})\frac{\partial \Delta\phi}{\partial\epsilon_{n+1}^{k}}\Bigg{]}\] (A6)  
\[\frac{\partial q}{\partial\epsilon_{n+1}^{k}}=\frac{\partial q_{n+1}^{\pi}}{ \partial\epsilon_{n+1}^{k}}-6\mu\Delta\phi\frac{\partial q}{\partial\epsilon_ {n+1}^{k}}-6\mu q\frac{\partial\Delta\phi}{\partial\epsilon_{n+1}^{k}}\] (A7)  
In Equations (A5) and (A6), let  
\[B=\Bigg{(}2\frac{\partial p}{\partial\epsilon_{n+1}^{k}}-\frac{\partial p_{0} }{\partial\epsilon_{n+1}^{k}}\Bigg{)}\Delta\phi+(2P+P_{s}-P_{0})\frac{\partial \Delta\phi}{\partial\epsilon_{n+1}^{k}}\] (A8)  
Combined with (A6) and (A8), we can see  
\[B=\frac{1}{P_{0}\vartheta(s)M^{2}}\frac{\partial p_{0}}{\partial\epsilon_{n+1 }^{k}}\] (A9)  
By substituting Equation (A9) into Equation (A5), we can see  
Figure 9:0 (A10)  
By substituting Equation (A10) into Equation (A5), we can see that  
Figure 9:1 (A11)  
Figure 9:2 (A12)  
The above formula can be obtained  
Figure 9:3 (A13)In the above equation  
\[\left\{ \begin{array}{l}a_{1}=\frac{1+p_{0}\theta(s)M^{2}\Delta\phi}{a}\\ a_{2}=\frac{-M^{2}(2P+P_{s}-P_{0})}{a}\\ a=1+2 KM^{2}\Delta\phi+p_{0}\theta(s)M^{2}\Delta\phi\end{array} \right. \tag{101}\]  
By combining Equations (100) and (101), we can see  
\[\frac{\partial p}{\partial e_{n+1}^{k}}=\frac{\partial p_{n+1}^{j}}{\partial e _{n+1}^{k}}-\frac{K\partial p_{0}}{p_{0}\theta(s)\partial e_{n+1}^{k}} \tag{102}\]  
Substitute Equation (101) into Equation (100)  
\[\frac{\partial p_{0}}{\partial e_{n+1}^{k}}=p_{0}\theta(s)M^{2}\Bigg{[}\bigg{(} 2\frac{\partial p_{n+1}^{n}}{\partial e_{n+1}^{k}}-\frac{2K\partial p_{0}}{p_ {0}\theta(s)\partial e_{n+1}^{k}}-\frac{\partial p_{0}}{\partial e_{n+1}^{k} }\bigg{)}\Delta\phi+(2P+P_{s}-P_{0})\frac{\partial\Delta\phi}{\partial e_{n+1 }^{k}}\Bigg{]} \tag{103}\]  
After finishing, we can obtain  
\[\frac{\partial p_{0}}{\partial e_{n+1}^{k}}=a_{3}KI+a_{4}K\frac{\partial \Delta\phi}{\partial e_{n+1}^{k}} \tag{104}\]  
In the above equation  
\[\left\{ \begin{array}{l}a_{3}=\frac{2p_{0}\theta(s)M^{2}\Delta\phi}{a}\\ a_{4}=\frac{p_{0}\theta(s)M^{2}(2P+P_{s}-P_{0})}{a}\end{array} \right. \tag{105}\]  
By sorting out Equation (100), we can get  
\[\frac{\partial q}{\partial e_{n+1}^{k}}=2a_{5}\mu\hat{\mathbf{n}}+2a_{6}\mu \frac{\partial\Delta\phi}{\partial e_{n+1}^{k}} \tag{106}\]  
In the above equation  
\[\left\{ \begin{array}{l}a_{5}=\sqrt{\frac{3}{2}}(1+6\mu\Delta\phi)^{-1}\\ a_{6}=\frac{-3q}{1+6\mu\Delta\phi}\end{array} \right.