# Grouting Process Simulation Based on 3D Fracture Network Considering Fluid-Structure Interaction
## 2 Methodology
### Fluid-Structure Interaction Model
#### Computational Fluid Dynamics (CFD) Grouting Numerical Model  
The governing equations of the grouting can be described by the continuity equation, momentum equation, two-phase volume of fluid (VOF) equation and Papanastasiou regularized equation. 1. The two-phase VOF equation In the process of grouting, the grout drives out air or groundwater, which should be treated as a two-phase flow [4]. The accurate description of the interface between two kinds of incompatible and incompressible fluids is one of the most important issues in multi-fluid flow computations [39], this can be solved by the VOF method which is proposed by [PERSON] [40] to track free fluid surfaces under fixed grid condition. Therefore, the VOF method is used to keep track of the grout-air interface in this paper. In this method, a volume fractional variable \(F=F\left(x,y,z,t\right)\) for each phase of the model in the computational domain is introduced. \(F_{g}=1\) indicates that the volume is occupied by grout while \(F_{g}=0\) indicates that the volume contains no grout and is in the air phase, and \(0<F_{g}<1\) stands for the volume that contains both grout and air. Equation (7) is used to describe the motion of the grout-air interface: \[\frac{\partial F_{g}}{\partial t}+\rho\
u\
abla F_{g}=0\] (7) where \(F_{g}\) is the volume fraction of grout; \(\rho\) is the density of fluid in kg/m\({}^{3}\); \(\
u\) is the kinematic viscosity of fluid in m\({}^{2}\)/s. 2. The continuity equation: \[\frac{\partial\rho}{\partial t}+\
abla(\rho u)=0\] (8) where \(\rho\) is the density of fluid in kg/m\({}^{3}\); \(t\) is the time in s; \(u\) is the velocity of the unit section in m/s. 3. The momentum equation: \[\rho\frac{du}{dt}=-\
abla p+\rho g+\
abla\left(\eta\dot{\gamma}\right)+F_{ \mathrm{st}}+S\] (9) where \(p\) is the pressure on the fluid micro-unit in \(Pa\), \(g\) is the acceleration of gravity in m/s\({}^{2}\); \(\eta\) is the apparent viscosity of fluid in Pa-s; the relationship between \(\
u\) and \(\eta\) is \(\
u=\eta\rho\); \(\dot{\gamma}\) is the shear rate in 1/s, \(F_{st}\) is the surface tension force in N/m\({}^{3}\) and is presented as Equation (10) [41]; and \(S\) is the momentum resistance source term in N/m\({}^{3}\), including inertia loss term \(S_{i}\) and viscosity loss term \(Sv\); in this paper, \(S_{i}\) can be neglected because of the low velocity of the grout and \(S=S_{v}\). Equation (11) is the expression of the viscosity loss term \(S_{v}\): \[F_{\text{st}}=-\sigma\
abla\cdot\left(\frac{\
abla F_{i}}{|\
abla F_{i}|} \right)\
abla F_{i}\] (10) where \(F_{i}\) is the volume fraction of phases; \(\sigma\) is the surface tension coefficient in N/m. \[S_{v}=-\frac{\rho\
u}{\alpha}u\] (11) where \(\frac{1}{\alpha}\) is the viscous drag coefficient and its expression is as follows: \[\frac{1}{\alpha}=\frac{g}{K\
u}\] (12) where \(K\) is the permeability coefficient. In order to obtain single set of equations, \(\rho\) and \(\
u\) in Equations (7)-(12) are no longer constants but are variables weighted by the volume fraction of fluid [42]: \[\rho= F_{\text{g}}\rho_{g}+\left(1-F_{\text{g}}\right)\rho_{a}\] (13) \[\
u= F_{\text{g}}\
u_{g}+\left(1-F_{\text{g}}\right)\
u_{a}\] (14) where \(\rho_{g}\), \(\rho_{a}\), \(\
u_{g}\), \(\
u_{a}\) are the density of grout, the density of air, the kinematic viscosity of grout and the kinematic viscosity of air, respectively. 4. The Papanastasiou regularized equation The cement grout with a w/c ratio of less than 1 is usually described by the Bingham model. However, in the Bingham constitutive equation, when the shear rate is close to zero, the apparent viscosity will become infinite, which causes problems in numerical simulation. In order to solve this problem, the Papanastasiou regularized model is used to describe the rheological properties of cement grout [4], as shown in Equation (15): \[\eta=\left\{ \begin{array}{ll}m\tau_{0}&\dot{\gamma}=0\\ \eta_{0}+\frac{\tau_{0}}{\dot{\gamma}}\left[1-e^{-m\dot{\gamma}}\right]&\dot{ \gamma}\
eq 0\end{array}