## 2 Materials and Methods
\begin{cases}\frac{-\alpha_{1}\overline{v}_{ij}\psi_{ij}+\beta_{1} \overline{v}_{ij}^{2}}{\overline{v}_{ij}},&\mathbf{v}_{ij}\mathbf{x}_{ij}<0\\ 0,&\mathbf{v}_{ij}\mathbf{x}_{ij}\geq 0\end{cases} , \tag{4}\]where  
\[\phi_{ij}=\frac{\overline{h}_{ij}\varphi_{ij}\mathbf{x}_{ij}}{\left|\mathbf{x}_{ij}\right|^ {2}+0.01\overline{h}_{ij}^{2}}, \tag{5}\]  
\[\overline{c}_{ij}=\frac{1}{2}\big{(}c_{i}+c_{j}\big{)}, \tag{6}\]  
\[\overline{\rho}_{ij}=\frac{1}{2}\big{(}\rho_{i}+\rho_{j}\big{)}, \tag{7}\]  
\[\overline{h}_{ij}=\frac{1}{2}\big{(}h_{i}+h_{j}\big{)}, \tag{8}\]  
\[\mathbf{v}_{ij}=\mathbf{v}_{i}-\mathbf{v}_{j}, \tag{9}\]  
\[\mathbf{x}_{ij}=\mathbf{x}_{i}-\mathbf{x}_{j}, \tag{10}\]  
where \(h\) denotes the smoothing length, and \(c\) represents the speed of sound. Both the parameters \(\alpha_{\Pi}\) and \(\beta_{\Pi}\) of the artificial viscosity are set to 1.0. Ejection involves significant material deformation, resulting in disordered particle distributions. The kernel gradient correction technique [48] is applied to enhance accuracy as follows:  
\[\
abla_{i}W_{ij}^{C}=\mathbf{L}_{i}\
abla_{i}W_{ij}, \tag{11}\]  
where  
\[\mathbf{L}_{i}=\left[\sum_{j=1}^{N}\frac{m_{j}}{\rho_{j}}\big{(}\mathbf{x}_{j}- \mathbf{x}_{i}\big{)}\otimes\
abla_{i}W_{ij}\right]^{-1}. \tag{12}\]  
Ejection involves shock wave loading and unloading as well as tensile fracture, and significant tensile stress develops within the material. Therefore, the artificial stress [49] is applied to remove the tensile instability as follows:  
\[f^{n}R_{ij}^{\alpha\beta}=\left[\frac{W\big{(}\mathbf{x}_{j}-\mathbf{x}_{i},h\big{)}}{ W(\Delta x,h)}\right]^{n}\left(R_{i}^{\alpha\beta}+R_{j}^{\alpha\beta}\right), \tag{13}\]  
where \(\Delta x\) denotes the initial particle distance, \(n\) is set to 4, and \(R_{i}^{\alpha\beta}\) is expressed as follows:  
\[\left\{ \begin{array}{l}R_{i}^{xx}=\cos^{2}\theta_{i}\overline{R}_{i}^{xx}+ \sin^{2}\theta_{i}\overline{R}_{i}^{yy}\\ R_{i}^{yy}=\sin^{2}\theta_{i}\overline{R}_{i}^{xx}+\cos^{2}\theta_{i} \overline{R}_{i}^{yy}\\ R_{i}^{xy}=\sin\theta_{i}\cos\theta_{i}\Big{(}\overline{R}_{i}^{xx}-\overline{ R}_{i}^{yy}\Big{)}\end{array} ,\right. \tag{14}\]  
where \(\overline{R}_{i}^{xx}\) and \(\overline{R}_{i}^{yy}\) are calculated as follows:  
\[\overline{R}_{i}^{xx}= \begin{cases}-\epsilon\frac{\overline{r}_{i}^{xx}}{ \overline{r}^{2}},&\overline{\sigma}_{i}^{xx}>0\\ 0,&\overline{\sigma}_{i}^{xx}\leqslant 0\end{cases} , \tag{15}\]  
\[\overline{R}_{i}^{yy}= \begin{cases}-\epsilon\frac{\overline{r}_{i}^{yy}}{ \overline{r}^{2}},&\overline{\sigma}_{i}^{yy}>0\\ 0,&\overline{\sigma}_{i}^{yy}\leqslant 0\end{cases} , \tag{16}\]  
where \(\epsilon\) is set to 0.3, \(\overline{\sigma}_{i}^{xx}\) and \(\overline{\sigma}_{i}^{yy}\) are expressed as follows:  
\[ \begin{cases}\overline{\sigma}_{ij}^{xx}=\cos^{2}\theta_{i}\sigma_{i}^{xx}+ 2\sin\theta_{i}\cos\theta_{i}\sigma_{i}^{xy}+\sin^{2}\theta_{i}\sigma_{i}^{yy} \\ \overline{\sigma}_{i}^{yy}=\sin^{2}\theta_{i}\sigma_{i}^{xx}-2\sin\theta_{i} \cos\theta_{i}\sigma_{i}^{xy}+\cos^{2}\theta_{i}\sigma_{i}^{yy},\end{cases}  \tag{17}\]where the rotation angle \(\theta_{i}\) is defined as follows:  
\[\tan 2\theta_{i}=\frac{2\sigma_{i}^{xy}}{\sigma_{i}^{xx}-\sigma_{i}^{yy}}. \tag{18}\]  
Ejection induces considerable compression and tension in the material, and maintaining a constant smoothing length would cause notable fluctuations in the number of neighboring particles, ultimately decreasing the precision of the computational outcomes. Hence, the variable smoothing length [50] is used as follows:  
\[\frac{\mathrm{d}h_{i}}{\mathrm{d}t}=-\frac{1}{2}\frac{h_{i}}{\rho_{i}}\frac{ \mathrm{d}\rho_{i}}{\mathrm{d}t}, \tag{19}\]  
In this study, the initial smoothing length \(h_{0}\) is equal to \(1.5\Delta x\), and then, the smoothing length is updated according to Equation (19). The kick-drift-kick method [47] is employed to solve Equation (3).