Study on the Microscopic Mechanism of Grouting in Saturated Water-Bearing Sand Stratum Based on VOF-DEM Method

Abstract

During the construction of tunnels in saturated water-bearing sand stratum, water and mud inrush disasters often occur. Grouting is the most convenient and effective method to improve the mechanical properties of the soil and prevent groundwater seepage. The reasonableness of the Discrete Element Method (DEM) contact parameters is verified by comparing the repose angle test with that obtained by simulations. The grouting model of saturated water-bearing sand stratum was established based on the Volume of Fluid-Discrete Element Method (VOF-DEM). Then, the effects of sand stratum porosity and grouting pressure on grouting were discussed. The results show that (1) in dense sands, the permeation diffusion of the slurry dominates, whereas in loose sands, the compaction zone is well developed. (2) Loose sand has a wider propagation of stress and dense sand has a larger increase in the stress state. (3) In loose sand the slurry diffusion is significantly affected by gravity, whereas in dense sand the slurry diffusion along the dominant path of the grouting pipe wall boundary is obvious. (4) The range of slurry diffusion and compaction zone is positively related to the porosity and grouting pressure. (5) In loose sand with grouting pressure below 200 kPa, no compaction diffusion occurs during the grouting process, whereas in dense sand, when the grouting pressure is below 500 kPa no compaction diffusion occurs during the grouting process.

1. Introduction

China’s urban rail transit is growing rapidly, with a total of 8708 km of subways nationwide by 2022. During the construction of urban subways, water-rich sand stratum is often encountered [1,2]. Due to its high porosity and permeability, the sand stratum can serve as an essential groundwater reservoir. Water and mud inrush often occur when tunneling in water-bearing sand stratum. Among the technical techniques for water plugging and reinforcement, grouting is the most convenient and effective [3,4]. Grouting is a technique that uses a grouting pipe to inject slurry into the ground. During this operation, the slurry displaces the air and groundwater between the soil particles. After slurry solidification and hardening, the soil is modified by the slurry into reinforced soil with exceptional waterproof ability, integrity, and strength [5,6].

Due to the hidden and complex features of grouting projects, grouting theory lags far behind engineering practice. In recent decades, researchers have studied the theory of grouting by various methods. Liu J et al. proposed guidelines for the selection of composite materials in different flowing water environments through site surveys [7]. Li Z et al. investigated the effects of the water–cement ratio and curing time on the reinforcement effect of sandy soils through laboratory tests [8]. Wang D et al. established a grouting test model for deep saturated sand layers and studied the diffusion of urea–formaldehyde resin in saturated sand layers [9]. Wang Q et al. established a grouting model for loose sandy soil to study the diffusion law of slurry in sandy soil from factors such as grouting pressure, water–cement ratio, and sandy soil saturation, respectively [10].

However, site surveys and laboratory tests were unable to reveal the diffusion mechanism of the slurry. Therefore, many researchers have used numerical simulation methods to study the grouting mechanism. Ren B et al. studied the diffusion pattern of fracture grouting under dynamic and hydrostatic water conditions by using the VOF model [11]. Cui W et al. introduced the viscosity time-varying equation into the VOF model and investigated the effect of viscosity time-varying velocity on the washout resistance of slurry [12]. Li S et al. introduced the Sequential Diffusion and Solidification (SDS) method into VOF to represent the spatial variation of viscosity by sequential injection. [13]. Liu X et al. used the DEM method to simulate the grouting process by using the “tube domain” to reveal the effect of grouting pressure on the form of slurry diffusion [14]. Boschi K et al. used the DEM method to discrete the slurry into particles. They investigated the micromechanics of grouting under different grouting pressures and different soil properties [15]. Yang C et al. used the Computational Fluid Dynamics–Discrete Element Method (CFD–DEM) algorithm to study the effects of different grouting pressures on slurry diffusion and soil mechanical mechanisms [16]. Guo W et al. used the two-phase Darcy law for grouting simulation to study the effect of different grouting methods and well spacing on the effect of slurry penetration on the leading-edge lap [17].

However, the above numerical studies mainly focus on the field of fracture grouting and unsaturated soil grouting. The microscopic mechanism of grouting in saturated water-bearing sand stratum is still unclear. The grouting process of saturated water-bearing sand stratum involves the three-phase coupling of slurry–water–sand, and the process is highly complex. The VOF method can accurately simulate the intersection of slurry and water, and the DEM method can simulate the real soil state. Therefore, the grouting process is simulated by coupling the VOF method and the DEM method. First, the reasonableness of the sandy soil DEM contact parameters is verified by comparing the repose angle test with the repose angle simulation test. Then a grouting model based on the VOF–DEM method is established to study the microscopic mechanism of grouting in saturated water-bearing sand stratum. The study results are of great significance for improving the theory of grouting in water-bearing sand stratum.

2. Methods

The slurry and water are divided by a distinct interface, known as the slurry–water intersection, during the grouting process in saturated water-bearing sand stratum. The sandy soil particles can be located in water or slurry, and their motion is governed by their interactions with the sandy soil particles and fluid around them. By coupling the DEM and VOF methods, the grouting process in saturated water-bearing sand stratum can be simulated.

2.1. Governing Equations for the Particle Phase

Cundall et al. first proposed the DEM for simulating the motion of granular solids [18]. In the DEM model, the particles are governed by Newton’s laws of motion. The conservation of momentum and angular momentum of each particle is given by

$$m_A\frac{d{v}_A}{dt}={\displaystyle \sum }_j{F}_{AB}^c+F_A^f+F_A^g$$

(1)

$$I_A\frac{d{\omega }_A}{dt}={\displaystyle \sum }_B{M}_{AB}$$

(2)

where $m_A$ is the mass of particle A, $v_A$ is the translational velocity of particle A,$F_{AB}^c$ is the contact force between particle A and particle B, $F_A^f$ is the fluid–particle interaction force acting on particle A, $F_A^g$ is the gravity acting on particle A, $I_A$ is the moment of inertia of particle A, ${\omega }_A$ is angular velocities of particle A, and $M_{AB}$ is contact moment acting on particle A by particle B.

The interaction model of the particles has a significant impact on the final computation results in DEM. As a result, choosing an acceptable contact model can improve the accuracy of the final simulation. Because sandy soil has low cohesion and large dispersion, the Hertz–Mindlin (no-slip) model was chosen as the particle contact model in this paper. Figure 1 illustrates its contact schematic diagram.

In the DEM, the contact force includes normal contact force Fn and tangential contact force Ft. The normal force is based on the contact theory of Hertzian [19], and the tangential force [20,21] depends on the tangential overlap and tangential stiffness, and their expressions are as follows.

$$F_n=\frac{4}{3}{E}^{*}\sqrt{R^{*}}{\delta }_n^{\frac{3}{2}}$$

(3)

$$F_t=-k_t{\delta }_t$$

(4)

where $E^{*}$ is the equivalent Young’s modulus, $R^{*}$ is the equivalent radius, ${\delta }_n$ is the normal overlap of particle A and particle B, $k_t$ is the tangential stiffness, and ${\delta }_t$ is the tangential overlap of particle A and particle B. According to the theory of Cundall [18], the tangential force cannot be greater than ${\mu }_s{F}_n$, where ${\mu }_s$ is the static friction coefficient.

$E^{*}$, $R^{*},$and $k_t$ are defined as

$$\frac{1}{E^{*}}=\frac{\left(1-{\nu }_A^2\right)}{E_A}+\frac{\left(1-{\nu }_B^2\right)}{E_B}$$

(5)

$$\frac{1}{R^{*}}=\frac{1}{R_A}+\frac{1}{R_B}$$

(6)

$$k_t=8G^{*}\sqrt{R^{*}{\delta }_n}$$

(7)

where$E_{A }$and$E_B$are the Young’s modulus of particle A and particle B, respectively; ${\nu }_A$and${\nu }_{B }$are the Poisson’s ratio of particle A and particle B, respectively$; R_A$and$R_{B }$are the radius of particle A and particle B, respectively; and$G^{*}$is the equivalent shear modulus.

In addition, there are normal damping forces$C_n$and tangential damping forces $C_s$which are given by [22,23]

$$C_n=-2\sqrt{\frac{5}{6}}\beta \sqrt{k_n{m}^{*}}{v}_n^{rel}$$

(8)

$$C_s=-2\sqrt{\frac{5}{6}}\beta \sqrt{k_t{m}^{*}}{v}_t^{rel}$$

(9)

where$m^{*}$is the equivalent mass, $v_n^{rel}$ is the normal relative velocity,$S_n$ is the normal stiffness,$\beta$is the equation for the restitution coefficient, and$v_t^{rel}$is the tangential relative velocity. $m^{*}$,$\beta ,$and$k_n$are defined as

$$m^{*}={\left(\frac{1}{m_A}+\frac{1}{m_B}\right)}^{-1}$$

(10)

$$\beta =\frac{-\ln e}{\sqrt{{\ln }^{2}e+{\pi }^2}}$$

(11)

$$k_n=2E^{*}\sqrt{R^{*}{\delta }_n}$$

(12)

where$m_A$and$m_B$are mass of particle A and particle B, respectively, and$e$is the restitution coefficient.

Rolling friction can be used to simulate sand particle shape, which can be illustrated in commercial software EDEM by applying torque to the contact surface [24].

$${\tau }_i=-{\mu }_r{F}_n{R}_i{\omega }_i$$

(13)

where${\mu }_{r }$is the rolling friction coefficient,$R_{i }$ is the distance between the contact point and the center of mass, and${\omega }_i$is the unit angular velocity vector of the particle at the contact point.

2.2. Governing Equations for the Fluid Phase

The VOF model can simulate two or more unmixed-phase fluids by solving a set of momentum equations and tracking the volume fraction of each fluid over the whole region. When dealing with the slurry–water two-phase flow, assume that the volume fraction of the slurry phase in each cell is${\alpha }_{s }$, and the volume fraction of the water phase is${\alpha }_w$. If${\alpha }_s=0$, then the cell contains no slurry. If${\alpha }_s=1$, then the unit contains only slurry. If$0<{\alpha }_s<1$, then there is a slurry–water interface in this cell. The dynamic behavior of the phase interface is solved by the following equation

$$\frac{\partial {\alpha }_s}{\partial t}+\left(\overrightarrow{v}\cdot 𝛻\right){\alpha }_s=0$$

(14)

where$\overrightarrow{v}$is the density and velocity of the fluid. The water phase volume fraction is found by the following equation.

$${\alpha }_w+{\alpha }_s=1$$

(15)

The momentum equations are related to the material physical parameters of density, viscosity, and volume fraction. The phase density equation and the phase viscosity equation are as follows.

$$\rho ={\rho }_w{\alpha }_w+{\rho }_s{\alpha }_s$$

(16)

$$\mu ={\mu }_w{\alpha }_w+{\mu }_s{\alpha }_{s }$$

(17)

The well-known volume-averaged control equation was used to describe the phase motion of the fluid [25]. The continuity equation is as follows.

$$\frac{\partial \left({\epsilon }_f\rho \right)}{\partial t}+𝛻\left({\epsilon }_f\rho \overrightarrow{v}\right)=0$$

(18)

The slurry flows slowly in porous media, so the laminar flow model is used. The momentum equation is given by

$$\frac{\partial \left({\epsilon }_f\rho \overrightarrow{v}\right)}{\partial t}+𝛻\left({\epsilon }_f\rho \overrightarrow{v}\overrightarrow{v}\right)={\epsilon }_f\left(-𝛻p+𝛻\mu \left(𝛻\overrightarrow{v}+𝛻{\overrightarrow{v}}^T\right)+\rho \overrightarrow{g}\right)+{\overrightarrow{F}}_{pf}$$

(19)

where${\epsilon }_f$is the fluid phase volume fraction, $\rho$is the density of the fluid, and${\overrightarrow{F}}_{pf}$is the fluid–particle interaction force acting on the fluid.

2.3. Fluid–Particle Interaction Force

Particle–fluid interaction forces include pressure gradient force, drag force, Magnus force, Saffman lift force, Bassett force, and virtual mass force [26]. In the grouting process, the interaction force between the fluid and the particles is mainly the drag force. The drag force occurs because of the different velocities of the particle and the fluid, resulting in the shear of the fluid on the solid particles, which acts along the relative velocity direction between the fluid and the particle. For dense particle populations, the Di–Felice drag force model [27] is generally chosen to calculate the particle–fluid interaction forces. Di–Felice’s drag force model is shown below

$$\begin{array}{c}{F}_d=\frac{\pi {\left(d_p\right)}^2}{2}{C}_D{\rho }_f\left|v_f-v_s\right|\left(v_f-v_s\right){\epsilon }_s^{-\chi +1}\\ \chi =3.7-0.65exp\left[-\frac{{\left(1.5-{\log }_{10}R{e}_s\right)}^2}{2}\right]\\ C_D={\left(0.63+4.8\sqrt{R{e}_s}\right)}^2\end{array}$$

(20)

where$d_p$is the particle diameter, $C_D$is the drag coefficient, ${\rho }_f$is the phase density of the fluid phase, ${\epsilon }_s$is the volume fraction of the particle phase, $v_f$is the velocity of the fluid phase at the center of the particle, $v_{s }$is the velocity of the particle phase, and$R{e}_s$is the Reynolds number of particles. In this paper, a recently developed virtual dual grid method [28] was used to calculate the volume fraction of the particle phase, which eliminates the limitation that the fluid grid volume must be larger than the particle volume in the conventional method.

2.4. Implementation of Coupled VOF–DEM Model

The VOF–DEM coupling process was completed on commercial software FLUENT and EDEM. Figure 2 illustrates the coupled simulation flow chart.

At the beginning of the calculation, FLUENT uses the user-defined initialization function, which uses the EDEM Application Programming Interface (API) to pass the particle position information and calculate the grid porosity and fluid–particle interaction force. EDEM starts the calculation and calculates the velocity and position of particles with Equations (1) and (2). After the particle state is updated, FLUENT starts to calculate the particle–fluid interaction forces with Equations (18) and (19). Once the iterations converge, the fluid states and interaction forces are updated. If the final computation time is reached, the computation is finished. If the final computation time is not reached, the computation cycle returns to the EDEM side and a new computation cycle begins.

If the time step is unreasonable, the calculation will be unstable. Therefore, the time step for the coupled VOF–DEM calculation should be determined. The information exchange between DEM and VOF is performed at each step of the VOF simulation. For single-phase fluid–particle coupling, the time step of the fluid calculation is usually set to N times the time step of particle calculation. However, to accurately understand the effect of particles on the motion of the slurry–water interface, VOF and DEM take the same time step. After much research and many simulations, the time step was finally set to 1 × 10⁻⁵ s.

3. Simulation Parameters of Sandy Soil

Figure 3 depicts the particle size distribution curve of the sandy soil. The sandy soil particles with particle sizes ranging from 0.075 to 5 mm accounted for 93.30 percent of the total mass, and the average particle diameter of the sandy soil was 1.43 mm.

The repose angle is a macroscopic parameter that determines granular material movement and friction. The deposition shape of sandy soil is generally conical, and the sand particles are deposited in a certain direction, forming a certain repose angle of deposition. It is a standard method for calibrating contact parameters of particles in DEM. In addition, Ghazavi M et al. determined the correlation between the repose angle and the internal friction angle of sandy soils through a series of laboratory tests, which further demonstrated that the repose angle can be used to calibrate sandy soil contact parameters [30]. The comparison of the repose angle test and the repose angle simulation test is commonly used in EDEM to ensure that the DEM contact parameters are suitable.

The lifting cylinder test was used to determine the repose angle of the sandy soil. A hollow cylinder is placed vertically on a base of known roughness. Then, as in the concrete slump test, the cylinder is lifted upward from the base at a specific speed to cause the sand to form a cone at the base directly below. The test setup and test results are shown in Figure 4a, with a cylinder diameter of 5 cm, a height of 15 cm, and a cylinder lifting speed of 5 mm/s. After filling the cylinder with sandy soil, the cylinder was lifted. After the sand cone was stabilized, vertical photographs were taken. The photographs were imported into CAD to measure repose angles. The test was repeated 10 times, and the final result of the average repose angle was 30.12°.

The sand stratum model requires a large number of sand particles. The computational speed is restricted by the computer’s performance and the number of particles. To improve the computing speed of DEM simulations, researchers typically increase particle size, decrease particle type, and then calibrate the contact coefficient. Finally, the sandy soil model replicates the contact and friction effects of real sandy soil. In DEM simulations, there is currently no strict standard for the magnification ratio of particle diameter. According to previous studies [31,32,33], the sandy soil particle diameters are usually magnified less than 10 times. In this paper, the particle diameter of all sandy soil particles was finally confirmed to be 2.86 mm. To ensure that the simulation was not distorted, the numerical model was built according to a ratio of 1:1 based on the measured rest angle test setup, as shown in Figure 4b. The repose angle simulation test is the same as the repose angle test. First, the cylinder is filled with particles. The cylinder is then propelled vertically upward at a velocity of 5 mm/s, allowing the sand particles to accumulate naturally due to gravity. Finally, the simulation concludes when all particles cease to move. The repose angle of sandy soil was calculated using EDEM’s protractor function. Figure 4b illustrates the repose angle simulation process. Following the completion of a significant number of repose angle simulation tests and a comparison of simulation and actual testing, a contact coefficient with an error of 1.03% was finally obtained, as shown in

Table 1. Contact parameters for sandy soil with a repose angle of 30.43°.

Contact Model	Hertz–Mindlin (No Slip)
Particle diameter (mm)	2.86
Density (kg/m3)	2600
Poisson’s ratio	0.3
Shear modulus (MPa)	10
Coefficient of restitution	0.5
Coefficient of static friction	0.5
Coefficient of rolling friction	0.3

.

4. Calculation Model and Parameters

As shown in Figure 5, a 3D geometry model for the grouting test was divided into three functional sections: the grout-inlet section, the test section, and the outlet section. The sandbox was set to 10 × 10 × 10 cm. The diameter and length of the grouting pipe were set to 1 cm and 5 cm, respectively. The 3D geometry model was created in the Geometry module of ANSYS. Boundaries and meshes were set via the Meshing module of ANSYS.

The simulation used EDEM to model the sandy soil and FLUENT to model the groundwater and slurry. Before VOF–DEM coupling, sandy soil particles were generated in the sandbox in EDEM. A dense sand model with a porosity of 0.3 was filled with 56,929 particles and a loose sand model with a porosity of 0.4 was filled with 48,795 particles. Figure 6 depicts their porosity contours, where the porosity is fairly uniform and the influence of rigid walls on the porosity is obvious. During the coupling calculation, the slurry flowed vertically from the slurry inlet with a certain pressure to displace the groundwater in the sand layer. When a two-phase flow is transported in the computational domain using FLUENT, the flow field information is continually updated throughout the coupling with EDEM; hence, the transient calculation was used. The water and slurry were incompressible, so a pressure-based solver was chosen. The slurry and water flowed relatively slowly in the sand stratum, so the laminar flow model was chosen. The slurry flowed in the y direction, with the inlet (y = 0) set to be the pressure inlet and the outlet (y = 5 cm) set to be the pressure outlet. Cement slurry with a water–cement ratio of 1 was used as the grouting material. The pore space of the sand stratum was filled with groundwater. The parameter settings in FLUENT are shown in

Table 2. Slurry and water parameters in FLUENT.

Fluid	ρ (kg/m³)	μ (Pa·s)
Cement slurry	1510	0.0123
Water	1000	0.001

.

To study the effects of grouting pressure and sand porosity on the diffusion mechanism of grouting in saturated water-bearing sand stratum, numerical simulation schemes were designed, as shown in

Table 3. Numerical simulation scheme.

Scheme	nsand (-)	Pinj (kPa)
1	0.3	50
2	0.3	200
3	0.3	500
4	0.4	50
5	0.4	200
6	0.4	500

.

5. Analysis of Calculation Results

The grouting time was 0.05 s in all of the tests discussed in this section. Figure 7 shows a comparison of the final arrangement of the slurry and sand particles under different grouting conditions. The spatial distribution of the slurry along the X- and Z-axes is depicted in Figure 7 in the central portion, where the thickness equals the pipe diameter. The qualitative analysis of all injection simulations demonstrates that two distinct diffusion mechanisms occurred within the sandy soil.

• Compaction diffusion mechanism: Due to the dislocation of sandy soil particles and the formation of only slurry, it did not contain a sandy soil compaction zone.
• Permeation diffusion mechanism: The slurry flows in the pores of the sandy soil did not seriously change the microstructure of the sandy soil.

It is obvious from Figure 7 that when the initial density was large enough (Figure 7c), permeation diffusion dominated over compaction diffusion. The size of the compaction zone was small and most of the slurry permeated into the pores of the sand stratum. In contrast, in the loose sand (Figure 7b), compaction diffusion was prevalent and the compaction zone was well developed. Relative to Figure 7c (dense state), the displacement vector amplitude of each sand particle in Figure 7b (loose state) was larger. In fact, in loose sands, the wider void space allowed the rearrangement of sand particles to be less constrained, whereas in dense sands, the dense tissue structure limited the displacement of sand particles. This constraint resulted in a higher increase in contact force due to slurry injection relative to the loose sandy soil (Figure 8). The displacement vector amplitude of each sandy soil particle in Figure 7b (high grouting pressure) was larger compared to Figure 7d (low grouting pressure).

Figure 8 illustrates the variation in the force chain during grouting. At the early stage of grouting, both dense sandy soil and loose sandy soil formed strong chains, at which time compaction diffusion prevailed over permeation diffusion. The propagation of force chains in loose sand was more extensive than in dense sand. In loose sands (Figure 8a), the wider pore space made the sand more mobile and more conducive to the propagation of force chains. In contrast, in dense sands (Figure 8b), the dense skeleton, although leading to larger contact forces, also limited the propagation of force chains. With time, the strong force chain gradually disappeared, at which time permeation diffusion was dominant over compaction diffusion. After the grouting was completed, a small number of medium force chains also appeared in the dense sand, whereas only weak force chains were present in the loose sand. From Figure 8d,e, it can be seen that the boundary of the compaction zone of loose sand was far away from the high-flow-rate zone of the slurry and the particles were subjected to less drag force compared to the dense sand. From the degree of velocity decay of the slurry in Figure 8d,e, it can be seen that the dense sandy soil had a high viscous resistance to the slurry, which was not conducive to the diffusion of the slurry. Finally, it is worth noting that the grouting process only caused a large increase in contact force near the grouting hole and did not cause a large increase in contact force throughout the specimen, as the stresses generated by the slurry during permeation were mainly borne by the pore water.

The pore space distribution changed after the slurry was injected, and some of the pores were filled with slurry. Because the sandbox’s boundary was fixed during the injection process, the development of the sand particle distribution was a disruption of the material’s homogeneity rather than a change in its volume average. Since the slurry was connected to the sand particles by drag, the slurry did not exist as a physical entity in the DEM. Here the variation in pore space structure was considered only for sandy soil particles. In Figure 9, the local variations in porosity of loose and dense sands were compared. In particular, the text showed clouds of porosity changes associated with slices of the specimens perpendicular to the injection direction (Figure 9a,c) and central slices parallel to the injection direction (Figure 9b,d). The thickness of the two slices was the diameter of the grouting pipe.

From the observed rearrangement of sand particles, the skeletal damage of dense sand at the corresponding area of injection was small, whereas the skeletal damage of loose sand was relatively large. Regardless of dense sand or loose sand, the arrangement of the sand particles’ position in the slurry permeation region was unchanged. The vast majority of the stresses generated by the slurry on the saturated sandy soil during permeation diffusion were borne by the pore water. For saturated soils, after external loading, the corresponding stress was generated in the soil, which was borne by both soil particles and pore water. However, only the stress borne by the soil particles deformed the soil, whereas the stress transmitted through the groundwater in the pore space did not affect the strength or deformation of the soil.

Figure 10 illustrates the slurry diffusion boundaries of the six grouting tests after 0.05 s of grouting (slurry saturation equal to 10%). It can be seen that the diffusion range of the slurry was positively related to the grouting pressure and porosity. The slurry in both loose and dense sand was symmetrically diffused in the horizontal plane. The dominant path along the boundary of the grouting pipe wall, where the porosity was higher, was evident (Figure 10a). The slurry in dense sand was symmetrically diffused in the vertical plane, whereas it was asymmetrically diffused in loose sand, and the downside was larger than the upside, for two reasons:

• The larger compaction zone made the diffusion center of the slurry move downward.
• The slurry diffusion was influenced by gravity due to the low viscous resistance of loose sand to the slurry.

The boundaries of the compaction zone for the six grouting tests are given in Figure 11. It can be seen that the extent of the compaction zone was positively correlated with the grouting pressure and porosity. Figure 12 illustrates the forces on the sand particles during the flow of the slurry. In dense sand, the contact forces on the sand particles were more and greater, the drag force on the sand particles during the flow of the slurry was not enough to resist the contact forces, and the effect on the particle skeleton structure was smaller. In loose sand, the contact force on sand and soil particles was less and smaller, and the drag force on sand and soil particles during the flow of slurry was enough to resist the contact force, which had more influence on the particle skeleton structure. At this time, the diffusion of slurry not only repelled the pore water but also repelled part of the sandy soil particles. In particular, the slurry flow rate in the loose specimen was greater than that in the loose specimen under the same grouting pressure condition (Figure 8d,e), which made the above phenomenon more obvious.

The compaction direction was in the negative direction of the Y-axis, so the compaction rate CR was proposed to quantify the degree of compaction of the slurry under each condition, which is defined as

$$CR=\frac{l_c}{l}$$

(21)

where$l_{c }$is the maximum compaction distance in the negative direction of the Y-axis, and $l$is the maximum diffusion distance in the negative direction of the Y-axis.

Table 4. Results of the compaction ratio.

Scheme	l [mm]	lc [mm]	CR [%]
1	22.6	0	0
2	30	0	0
3	44.1	3.2	7.25
4	23.8	0	0
5	37.1	3.3	8.89
6	50.0	12.9	25.8

shows the calculated results of CR after completion of the grouting. In the case of loose sand, the CR was 25.8% when the grouting pressure was 500 kPa. When its grouting pressure was reduced to 200 kPa, its CR dropped to 8.89%. When its grouting pressure was reduced to 50 kPa, no compaction occurred. In the case of dense sand, when the grouting pressure was 500 kPa, the CR was 7.25%. When the grouting pressure was lower than 500 kPa, no compaction occurred. It can be seen that in loose sand when the grouting pressure was less than 200 kPa, no compaction occurred during the grouting process. In dense sand, when the grouting pressure was less than 500 kPa, no compaction occurred during the grouting process.

Figure 13 depicts the amount of grout injected with the same cement slurry (with a w/c ratio of 1) at increasing grouting pressure with degrees of porosities of 0.3 and 0.4. At a grouting pressure of 50 kPa, 7.4 mL of slurry was injected into sandy soil with a porosity of 0.3, whereas 11.5 mL of slurry was injected into sandy soil with a porosity of 0.4. As the grouting pressure increased, the injected volume increased almost linearly for both porosity of 0.3 and 0.4, reaching 63.7 and 81.3 mL, respectively. Figure 13 also clearly shows that for a porosity of 0.3, the injection volume was much lower than for a porosity of 0.4.

In the two cases, the increase rate of grouting volume with the increase in grouting pressure was different, and the increase rate of grouting volume in sandy soil with a porosity of 0.4 was larger.

6. Conclusions

To study the microscopic mechanism of grouting in saturated water-bearing sand layers, this paper adopted the coupled VOF–DEM method to simulate the grouting process in saturated water-bearing sand layers. Through simulations, some conclusions were drawn, as explained below.

(1)

In dense sands, the permeation diffusion of the slurry dominated, whereas in loose sands, the compaction zone was well developed.

(2)

The evolution of the force chain state during grouting indicates the loose sand had a wider propagation of stress and the dense sand had a larger increase in the stress state. The diffusion model of the slurry changed from compaction diffusion dominated by permeation diffusion to permeation diffusion dominated by compaction diffusion.

(3)

The slurry diffusion range and injected grout volume show that the larger the grouting pressure and the larger the sand layer porosity the more favorable the slurry diffusion. When the sand layer porosity was small, the slurry diffusion phenomenon along the dominant path of the grouting pipe wall boundary was obvious. When the sand layer porosity was larger, the permeation diffusion of the slurry was significantly affected by gravity.

(4)

In loose sand, when the grouting pressure was less than 200 kPa, no compaction occurred during the grouting process. In dense sand, when the grouting pressure was less than 500 kPa, no compaction occurred during the grouting process.