## 2 Basic Theory and Realization Method of SSD Coupling
### Basic Theory of SSD Coupling
The computational space domain is discretized to obtain the seepage-stress coupling equation:  
\[[\mathbf{K}][\mathbf{U}]=\{\mathbf{F}_{\mathbf{V}}\}+\{\mathbf{F}_{\mathbf{s}}\}+\{\mathbf{F}_{\mathbf{p}}\}+ \{\mathbf{F}_{\mathbf{\sigma}_{\mathbf{0}}}\} \tag{4}\]  
where \([\mathbf{K}]\) is the structural stiffness matrix; \([\mathbf{U}]\) is the nodal displacement matrix; \(\mathbf{F}_{\mathbf{V}}\) and \(\mathbf{F}_{\mathbf{s}}\) are the body and surface loads, respectively; \(\mathbf{F}_{\mathbf{p}}\) is the equivalent load formed by pore pressure; and \(\mathbf{F}_{\mathbf{\sigma}_{\mathbf{0}}}\) is the initial stress load. According to the incremental theory of plasticity in the plastic damage model, the total strain tensor, \(\varepsilon\), is composed of the elastic strain rate, \(\varepsilon^{el}\), and the equivalent plastic strain rate, \(\varepsilon^{pl}\):  
\[\varepsilon=\varepsilon^{el}+\varepsilon^{pl} \tag{5}\]  
When there is no damage to the concrete, the stress-strain relationship of the concrete is as follows:  
\[\mathbf{\sigma}=\mathbf{D}^{el}\Big{(}\varepsilon-\varepsilon^{pl}\Big{)} \tag{6}\]  
where \(\mathbf{\sigma}\) is the total stress, and \(\mathbf{D}^{el}\) is the elastic stiffness matrix. When the concrete material is damaged, according to the theory of continuum damage mechanics, the internal micro-cracks, micro-pores, and other micro-defects under the action of external loads can be described by the damage factor, \(\mathbf{d}\). The damage factor is mainly used to reflect the concrete stiffness degradation under uniaxial or multiaxial loads. Assuming that the damage is isotropic, then, the relationship between the damage and stress of the concrete under the three-dimensional multiaxial state can be expressed by the damage elasticity equation, and the concrete stress, \(\mathbf{\sigma}\), is calculated as follows [36]:  
\[\mathbf{\sigma}=(\mathbf{1}-\mathbf{d})\overline{\mathbf{\sigma}}=(\mathbf{1}-\mathbf{d})\mathbf{D}^{\mathbf{e} \mathbf{I}}\Big{(}\mathbf{\epsilon}-\mathbf{\epsilon}^{\mathbf{p}\mathbf{I}}\Big{)} \tag{7}\]  
where \(\overline{\mathbf{\sigma}}\) is the effective stress, which represents the stress on the net section of the concrete material. The element damage, \(\mathbf{d}\), is expressed by the equivalent plastic strain \(\widetilde{\mathbf{\epsilon}}^{\mathbf{p}\mathbf{I}}\):  
\[\left\{ \begin{array}{l}\mathbf{d_{t}}=\mathbf{d_{t}}\Big{(}\widetilde{\mathbf{\epsilon} }^{\mathbf{p}\mathbf{I}}_{t}\Big{)},\mathbf{0}\leq\mathbf{d_{t}}\leq\mathbf{1}\\ \mathbf{d_{c}}=\mathbf{d_{c}}\Big{(}\widetilde{\mathbf{\epsilon}}^{\mathbf{p}\mathbf{I}}_{c}\Big{)},\mathbf{0}\leq\mathbf{d_{c}}\leq\mathbf{1}\end{array} \right. \tag{8}\]  
where \(\mathbf{d_{t}}\) is the tensile damage factor, \(\mathbf{d_{c}}\) is the compressive damage factor, \(\mathbf{t}\) is the tensile state, and \(\mathbf{c}\) is the compressive state. The equivalent plastic strain \(\widetilde{\mathbf{\epsilon}}^{\mathbf{p}\mathbf{I}}\) is calculated as follows:  
\[\left\{ \begin{array}{l}\widetilde{\mathbf{\epsilon}}^{\mathbf{p}\mathbf{I}}_{t}=\int_{ 0}^{t}\widetilde{\mathbf{\epsilon}}^{\mathbf{p}\mathbf{I}}_{t}\mathbf{d}\mathbf{t}\\ \widetilde{\mathbf{\epsilon}}^{\mathbf{p}\mathbf{I}}_{c}=\int_{0}^{t}\widetilde{\mathbf{ \epsilon}}^{\mathbf{p}\mathbf{I}}_{c}\mathbf{d}\mathbf{t}\end{array} \right. \tag{9}\]  
\[\left\{ \begin{array}{l}\dot{\widetilde{\mathbf{\epsilon}}}^{\mathbf{p}\mathbf{I}}_{t}= \mathbf{r}(\widehat{\mathbf{\sigma}})\dot{\widetilde{\mathbf{\epsilon}}}^{\mathbf{p}\mathbf{I}}_{ max}\\ \dot{\widetilde{\mathbf{\epsilon}}}^{\mathbf{p}\mathbf{I}}_{c}=-\big{(}\mathbf{1}-\mathbf{r}( \widehat{\mathbf{\sigma}})\big{)}\dot{\widetilde{\mathbf{\epsilon}}}^{\mathbf{p}\mathbf{I}}_{ min}\end{array} \right.