## 5 Viscoelastic Artificial Boundary  
For the problems of wave propagation in saturated porous media, a viscoelastic artificial boundary for the _u-p_ equations is adopted. It is expressed as follows:  
\[f_{uBl}=-\mathbf{K}_{Bl}^{\infty}\mathbf{u}_{Bl}-\mathbf{C}_{Bl}^{\infty}\mathbf{u}_{Bl}-\mathbf{ Q}_{Bl}^{\infty}p_{Bl} \tag{21}\]  
\[f_{qBl}=H_{Bl}^{\infty}p_{Bl}+S_{Bl}^{\infty}\dot{p}_{Bl} \tag{22}\]  
where \(\infty\) denotes the exterior domain, \(Bl\) is the node \(l\) on the boundary, \(\mathbf{f}_{uBl}\) and \(f_{qBl}\) are the stress and flow velocity of the node \(l\) on the boundary, respectively, and \(\mathbf{u}_{Bl}\) and \(p_{Bl}\) are the displacement and pore pressure of the node \(l\) on the boundary, respectively. The coefficients in the above equations are expressed as follows:  
\[ \begin{array}{c}\mathbf{K}_{Bl}^{\infty}=\mathbf{W}^{T}L_{l}\left[\begin{array}{cc }K_{BN}^{\infty}&0\\ 0&K_{BT}^{\infty}\end{array}\right]\mathbf{W};\;\mathbf{C}_{Bl}^{\infty}=\mathbf{W}^{T}L_ {l}\left[\begin{array}{cc}C_{BN}^{\infty}&0\\ 0&C_{BT}^{\infty}\end{array}\right]\mathbf{W};\\ \mathbf{Q}_{Bl}^{\infty}=\mathbf{W}^{T}L_{l}\left[\begin{array}{cc}Q_{B}^{\infty}\\ 0\end{array}\right]\\ \end{array}  \tag{23}\]  
\[H_{Bl}^{\infty}=L_{l}H_{B}^{\infty}S_{Bl}^{\infty}=L_{l}S_{B}^{\infty} \tag{24}\]  
where \(L_{l}\) is the length corresponding to node \(l\) on the boundary, and \(\mathbf{W}\) is the coordinate transformation matrix for transforming the local coordinates to the global coordinates. The parameters in the above equations are specifically  
\[\mathbf{W}=\left[ \begin{array}{cc}l_{rx}&l_{ry}\\ l_{\perp x}&l_{\perp y}\end{array} \right] \tag{25}\]  
\[K_{BN}^{\infty}=\frac{\lambda+2G}{2r(1+A)};\;K_{BT}^{\infty}=\frac{G}{2r(1+A)} \tag{26}\]  
\[C_{BN}^{\infty}=B\rho C_{P1}^{r};\;C_{BT}^{\infty}=B\rho C_{S}^{r} \tag{27}\]  
\[Q_{B}^{\infty}=\alpha;\;H_{B}^{\infty}=\frac{\overline{k}}{2r(1+A)};\;S_{B}^{ \infty}=\frac{\overline{k}}{C_{P1}^{r}} \tag{28}\]  
In Equation (25), \(l_{rx}=\cos(r,x)\) is the cosine of the angle between the positive direction of the local \(r\) axis and the positive direction of the global \(x\)-axis. Other parameters in Equation (25) are defined similarly to \(l_{xy}\). A = 0.8 and B = 1.1. The artificial boundary condition as shown in Equations (21) and (22) actually changes only the diagonal values of the coefficient matrix corresponding to boundary nodes in local node system. Assume that there are \(m\) nodes on the boundary in the local node system of node \(k\). Equations (21) and (22) are substituted into Equations (8) and (9) to obtain the finite element equations considering the viscoelastic artificial boundary:  
\[\mathbf{M}_{dk}\bar{\mathbf{u}}_{k}+\sum_{l=1}^{m}\mathbf{C}_{Bl}^{\infty}\dot{\mathbf{u}}_{Bl }+\left(\sum_{L=1}^{n}\sum_{j=1}^{4}\mathbf{K}_{kj}^{L}\mathbf{u}_{lj}^{L}+\sum_{l=1}^ {m}K_{Bl}^{\infty}\mathbf{u}_{Bl}\right)-\left(\sum_{L=1}^{n}\sum_{j=1}^{4}\mathbf{Q}_ {kj}^{L}p_{j}^{L}-\sum_{l=1}^{m}Q_{Bl}^{\infty}p_{Bl}\right)=\sum_{L=1}^{n} \mathbf{f}_{k,u}^{L} \tag{29}\]  
\[\left(S_{dk}\dot{p}_{k}+\sum_{l=1}^{m}S_{Bl}^{\infty}\dot{p}_{Bl}\right)+ \left(\sum_{L=1}^{n}\sum_{j=1}^{4}H_{kj}^{L}p_{j}^{L}+\sum_{l=1}^{m}H_{Bl}^{ \infty}p_{Bl}\right)+\sum_{L=1}^{n}\sum_{j=1}^{4}\left(\mathbf{Q}_{kj}^{L}\right) ^{T}\dot{\mathbf{u}}_{lj}^{L}=\sum_{L=1}^{n}f_{k,l}^{L} \tag{30}\]Compared to Equation (8), Equation (29) has an additional damping term and can be solved by an explicit integration method [47], and the equations of the local node system of node \(k\) at the \(t+1\) time step can be expressed as:  
\[