# Unsteady Approximate Model of Grouting in Fractured Channels Based on Bingham Fluid

Authors: Bo Ren; Lianguo Wang; Hao Fan; Ke Ding; Kai Wang; Chongyang Jiang
Journal: Geofluids
 (2021)

---

## Abstract

It is of great significance for the improvement of grouting technology and engineering practice to master the flow law of grout between parallel plates. However, the traditional calculation model ignores the influence of the inertia term and only considers the stable flow of slurry, so there is a big error in some cases. It is difficult to solve the motion equation of a Bingham fluid considering the inertial force term directly. Combined with the relationship between the steady-state flow equation of a Bingham fluid and a Newtonian fluid, the approximate unsteady-state flow equation of a Bingham fluid suitable for describing slurry flow is constructed. In addition, according to the unsteady flow equation, the relationship between the time and distance of slurry flow in parallel plate fractures can be obtained, and the simplified conditions of the Bingham fluid unsteady flow model are given. Finally, the accuracy of the flow equations and the simplified conditions are verified by experiments and numerical calculations.


---

## Body

## 1. Introduction

Grouting technology is widely used in many fields such as mining, tunneling, landfilling of municipal solid waste, storage of nuclear waste, and development and utilization of underground space. In coal mining, a soft rock roadway is often strengthened by grouting to prevent roadway damage or large deformation during roadway excavation [1–4]. Bad geological bodies such as lava and faults will have many adverse effects on the construction during the tunneling process. Grouting into the faults can increase the integrity of the rock mass and ensure safe and efficient construction [5–7]. Landfill of municipal solid waste [8, 9] requires grouting and antiseepage treatment in the buried pits that have been dug in advance to avoid water and soil pollution caused by toxic gas and liquid penetration. In the same way, nuclear waste storage [10, 11] also needs to consider the environmental impact of radioactive materials and consider strengthening the surrounding areas of the rock to avoid irreparable damage. In the long run, the scope of human activities is becoming wider and wider, and the use of underground space is gradually being put on the agenda. The use of underground space will certainly use grouting technology. So far, the method of grouting is more environmentally friendly, economical, and effective than other methods.However, in a grouting project, the important construction parameters such as slurry water-cement ratio, grouting pressure, and time are often judged by experience. Studying the flow of slurry in the fractured rock mass can guide on-site grouting and can aid in properly designing the grouting parameters [12]. A grouting project is of great significance. Therefore, it is becoming more and more important to study the characteristics of the slurry itself and the flow of the slurry in the rock mass [13].Fracture is the medium into which the slurry flows in a rock mass, and the simplest type is the smooth parallel plate fracture. By studying the parallel plate fracture, the flow law of the slurry in the rock mass can be well grasped, so some scholars have studied this problem and have obtained a series of results. Among them, the cubic law [14–16] has been widely used because of its accuracy and simplicity. Through theoretical derivation, it is found that the cubic law can be obtained from the Navier-Stokes (N-S) equation [17, 18]. To be precise, the cubic law is the theoretical solution of an incompressible Newtonian fluid doing steady laminar flow in smooth parallel fractures. However, grouting materials are mostly non-Newtonian fluids. Therefore, more flow equations have been proposed. Wereley and Pang [19] proposed a flow equation for the steady laminar flow of an incompressible Bingham fluid in smooth parallel fractures. Yan and Koplik [20] derived the flow equation of the steady laminar flow of an incompressible power law fluid in a smooth parallel fracture. Zou et al. [21] also obtained the flow equation of a Herschel-Bulkley (H-B) fluid which is more general.The above models all have the same disadvantage, that is, they ignore the effect of inertial terms. In some cases, there is a large calculation error. In addition, at present, many grouting materials have a very short solidification time, ranging from a few seconds to a few tens of seconds, and the viscosity changes significantly with time, so it is necessary to take the inertial term into account in the equation derivation. There is a “flow core” when a Bingham fluid flows. Therefore, the unsteady flow equation of a Bingham fluid needs to be divided into two sections. The equation is complicated and difficult to solve. In this paper, the unsteady-state flow model between parallel plates of a Newtonian fluid considering inertial force is first derived. Combined with the relationship between the steady-state solution of a Bingham fluid and a Newtonian fluid, the approximate unsteady-state flow equation of a Bingham fluid suitable for describing slurry flow is constructed and the approximate solution is verified by numerical calculation. Then, based on the unsteady-state flow model, the calculation method of the flow distance of slurry in the fractures of smooth parallel plates is obtained, and the simplified condition of a Bingham fluid unsteady-state flow model is given. Finally, different models are used to calculate the mud diffusion distance under different conditions, and the correctness of the simplified conditions is also proven.


## 2. Theoretical Model

### 2.1. Steady-State Flow Model of Slurry in Fracture

Firstly, the steady-state flow equation is derived based on 2D N-S equations, and the coordinate system is shown in Figure1. The derivation process is based on the following assumptions:
(1)
The flow state of slurry is laminar flow(2)
The slurry is incompressible(3)
The grouting process has reached a stable state(4)
The property of slurry does not change during groutingFigure 1
Diagram of a parallel plate fracture channel.The differential equation of viscous fluid motion is expressed by stress:(1)∂u∂t+u∂u∂x+v∂u∂y=fx+1ρ∂pxx∂x+∂τxy∂y,∂v∂t+u∂v∂x+v∂v∂y=fy+1ρ∂pyy∂y+∂τyx∂x.The constitutive equation of a Newtonian fluid is shown as follows:(2)pxxτxyτyxpyy=μ2∂u∂x∂u∂y+∂v∂x∂v∂x+∂u∂y2∂v∂y−p00p.The continuity equation for incompressible fluids is shown as follows:(3)∂u∂x+∂v∂y=0.Bring equation (2) and equation (3) into equation (1) to obtain the N-S equation [22], and we have
(4)∂u∂t+u∂u∂x+v∂u∂y=fx−1ρ∂p∂x+μρ∂2u∂x2+∂2u∂y2,∂v∂t+u∂v∂x+v∂v∂y=fy+1ρ∂p∂y+μρ∂2v∂x2+∂2v∂y2.According to the four hypotheses mentioned above, we can get(5)v=0,∂v∂x=∂v∂y=0,∂u∂x=0.Therefore, the two-dimensional problem is simplified into a one-dimensional problem. The definite condition that makes it easy to get the problem is shown as follows:(6)μ∂2u∂y2+ΔPL=0,∂u∂yy=0=0,uy=±b/2=0.The velocity distribution equation can be solved as follows:(7)u=b28μΔPL1−2yb2.Integral along theY direction to obtain the average velocity equation, we have
(8)u¯=1b∫−b/2b/2udy=−b212μΔPL.Similarly, equations of the average velocity of a non-Newtonian fluid can be obtained. The average velocity equation of the steady flow of a Bingham fluid is shown as follows:(9)u¯=2b∫y0b/2udy+b0uy=y0=−b212μΔPL1−b0b21+12b0b,where b0=−2τ0L/Δp, which is the width of the “flow core,” like Figure 2.Figure 2
Velocity distribution of a Bingham fluid in a parallel plate fracture channel.Similarly, for a power law fluid, when we calculate the steady-state flow equation of a power law fluid, we only need to replace the constitutive equation of the fluid [4, 20, 23]:
(10)τ=Kγ̇n,γ̇=124∂u∂x2+2∂u∂y+∂v∂x2+4∂v∂y2.And then we can get,(11)K∂∂y∂u∂yn+ΔPL=0,∂u∂yy=0=0,uy=±b/2=0.It is not difficult to acquire the steady-state average velocity equation of a power law fluid:(12)u¯=−n4n+2ΔP2KL1/nbn+1/n.It is not easy to directly calculate the flow equation of an H-B fluid by the N-S equation. Therefore, based on the relationship between this model and the constitutive equation of the model mentioned above, the flow equation of an H-B fluid is obtained by analogy. The average velocity equation obtained is consistent with the research results of other scholars [21, 24].
(13)u¯=−n4n+2ΔP2KL1/nbn+1/n1−b0bn+1/n1+nn+1b0b.


### 2.2. Unsteady Flow Model of Slurry in Fracture

The above formula is based on the assumption of stable flow. Next, we only keep assumptions (1), (2), and (4) to further derive the unsteady flow equation. The influence of inertia force must be taken into account in the establishment of unsteady-state equations. The first is the unsteady flow solution based on a Newtonian fluid:
(14)∂u∂t−μρ∂2u∂y2=ΔPρL,u0,y=u0y,∂u∂yt,0=0,ut,−b2=ut,b2=0.Solving equation (14) is a pure mathematical process. First of all, it needs to be divided into two subproblems, and then the solutions of the two subproblems are superposed to get the original solution.Subproblem 1.
(15)∂u∂t−μρ∂2u∂y2=0,u0,y=u0y,∂u∂yt,0=0,ut,−b2=ut,b2=0.Subproblem 2.
(16)∂u∂t−μρ∂2u∂y2=ΔPρL,u0,y=0,∂u∂yt,0=0,ut,−b2=ut,b2=0.The solution of Subproblem1 is as follows:
(17)u1y,t=∑k=1∞ake−kπ/b2μ/ρtsinkπby+kπ2,k=1,2,⋯,whereak=2/b∫−b/2b/2u0ysinkπ/by+kπ/2dx.The solution of Subproblem2 is as follows:
(18)u2y,t=b2ΔP12μL∑k=1∞fYyfTt,k=1,2,⋯.(19)fYy=241kπ3−1k+1+1sinkπby+kπ2,fTt=1−e−kπ/b2μ/ρt.When the initial velocityu0y=0 in Subproblem 1, the solution of Subproblem 1 is always equal to 0. Therefore, the solution of the original problem is equal to the solution of Subproblem 2.
(20)unewy,t=u2y,t.According to equation (18), the unsteady solution of a Newtonian fluid flowing in parallel plates can be divided into three terms: the first term is the average velocity term representing the magnitude of the velocity, the second term is the velocity distribution term related to the coordinates, and the third is the time-dependent attenuation term.It is assumed that the time for a Newtonian fluid and a Bingham fluid to reach steady flow from 0 is the same. The unsteady flow equation of a Bingham fluid can be obtained by combining the relationship between the steady flow velocity equation of a Bingham fluid [25] and the steady solution of a Newtonian runny nose (equation (7)).
(21)ubiny,t=unewy,t⋅ftrans1y,(22)ftrans1y=1−2y/b2−2b0/b1−2y/b1−2y/b2,b02<y<b2.Through the unsteady solution of a Newtonian fluid and a Bingham fluid, we can get the fluid flow process in a parallel plate, which will go through the entrance stage, transition stage, and stable stage. The entrance stage is located at the entrance of the parallel plate, which is generally very short in length. In this stage, the boundary effect has not affected the whole flow layer, even for a Newtonian fluid; its velocity distribution will show a similar shape of “flow core.” In the transition stage, the boundary effect has affected the whole flow layer, but the influence of inertial force has not disappeared, and the velocity distribution is still changing. When it comes to the stable stage, the influence of the inertial term has disappeared, and the fluid enters a stable flow state (see Figure3).Figure 3
Distribution diagram of velocity in the entrance stage, the transition stage, and the steady stage when a Newtonian fluid and a Bingham fluid flow in fractures.


## 2.1. Steady-State Flow Model of Slurry in Fracture

Firstly, the steady-state flow equation is derived based on 2D N-S equations, and the coordinate system is shown in Figure1. The derivation process is based on the following assumptions:
(1)
The flow state of slurry is laminar flow(2)
The slurry is incompressible(3)
The grouting process has reached a stable state(4)
The property of slurry does not change during groutingFigure 1
Diagram of a parallel plate fracture channel.The differential equation of viscous fluid motion is expressed by stress:(1)∂u∂t+u∂u∂x+v∂u∂y=fx+1ρ∂pxx∂x+∂τxy∂y,∂v∂t+u∂v∂x+v∂v∂y=fy+1ρ∂pyy∂y+∂τyx∂x.The constitutive equation of a Newtonian fluid is shown as follows:(2)pxxτxyτyxpyy=μ2∂u∂x∂u∂y+∂v∂x∂v∂x+∂u∂y2∂v∂y−p00p.The continuity equation for incompressible fluids is shown as follows:(3)∂u∂x+∂v∂y=0.Bring equation (2) and equation (3) into equation (1) to obtain the N-S equation [22], and we have
(4)∂u∂t+u∂u∂x+v∂u∂y=fx−1ρ∂p∂x+μρ∂2u∂x2+∂2u∂y2,∂v∂t+u∂v∂x+v∂v∂y=fy+1ρ∂p∂y+μρ∂2v∂x2+∂2v∂y2.According to the four hypotheses mentioned above, we can get(5)v=0,∂v∂x=∂v∂y=0,∂u∂x=0.Therefore, the two-dimensional problem is simplified into a one-dimensional problem. The definite condition that makes it easy to get the problem is shown as follows:(6)μ∂2u∂y2+ΔPL=0,∂u∂yy=0=0,uy=±b/2=0.The velocity distribution equation can be solved as follows:(7)u=b28μΔPL1−2yb2.Integral along theY direction to obtain the average velocity equation, we have
(8)u¯=1b∫−b/2b/2udy=−b212μΔPL.Similarly, equations of the average velocity of a non-Newtonian fluid can be obtained. The average velocity equation of the steady flow of a Bingham fluid is shown as follows:(9)u¯=2b∫y0b/2udy+b0uy=y0=−b212μΔPL1−b0b21+12b0b,where b0=−2τ0L/Δp, which is the width of the “flow core,” like Figure 2.Figure 2
Velocity distribution of a Bingham fluid in a parallel plate fracture channel.Similarly, for a power law fluid, when we calculate the steady-state flow equation of a power law fluid, we only need to replace the constitutive equation of the fluid [4, 20, 23]:
(10)τ=Kγ̇n,γ̇=124∂u∂x2+2∂u∂y+∂v∂x2+4∂v∂y2.And then we can get,(11)K∂∂y∂u∂yn+ΔPL=0,∂u∂yy=0=0,uy=±b/2=0.It is not difficult to acquire the steady-state average velocity equation of a power law fluid:(12)u¯=−n4n+2ΔP2KL1/nbn+1/n.It is not easy to directly calculate the flow equation of an H-B fluid by the N-S equation. Therefore, based on the relationship between this model and the constitutive equation of the model mentioned above, the flow equation of an H-B fluid is obtained by analogy. The average velocity equation obtained is consistent with the research results of other scholars [21, 24].
(13)u¯=−n4n+2ΔP2KL1/nbn+1/n1−b0bn+1/n1+nn+1b0b.


## 2.2. Unsteady Flow Model of Slurry in Fracture

The above formula is based on the assumption of stable flow. Next, we only keep assumptions (1), (2), and (4) to further derive the unsteady flow equation. The influence of inertia force must be taken into account in the establishment of unsteady-state equations. The first is the unsteady flow solution based on a Newtonian fluid:
(14)∂u∂t−μρ∂2u∂y2=ΔPρL,u0,y=u0y,∂u∂yt,0=0,ut,−b2=ut,b2=0.Solving equation (14) is a pure mathematical process. First of all, it needs to be divided into two subproblems, and then the solutions of the two subproblems are superposed to get the original solution.Subproblem 1.
(15)∂u∂t−μρ∂2u∂y2=0,u0,y=u0y,∂u∂yt,0=0,ut,−b2=ut,b2=0.Subproblem 2.
(16)∂u∂t−μρ∂2u∂y2=ΔPρL,u0,y=0,∂u∂yt,0=0,ut,−b2=ut,b2=0.The solution of Subproblem1 is as follows:
(17)u1y,t=∑k=1∞ake−kπ/b2μ/ρtsinkπby+kπ2,k=1,2,⋯,whereak=2/b∫−b/2b/2u0ysinkπ/by+kπ/2dx.The solution of Subproblem2 is as follows:
(18)u2y,t=b2ΔP12μL∑k=1∞fYyfTt,k=1,2,⋯.(19)fYy=241kπ3−1k+1+1sinkπby+kπ2,fTt=1−e−kπ/b2μ/ρt.When the initial velocityu0y=0 in Subproblem 1, the solution of Subproblem 1 is always equal to 0. Therefore, the solution of the original problem is equal to the solution of Subproblem 2.
(20)unewy,t=u2y,t.According to equation (18), the unsteady solution of a Newtonian fluid flowing in parallel plates can be divided into three terms: the first term is the average velocity term representing the magnitude of the velocity, the second term is the velocity distribution term related to the coordinates, and the third is the time-dependent attenuation term.It is assumed that the time for a Newtonian fluid and a Bingham fluid to reach steady flow from 0 is the same. The unsteady flow equation of a Bingham fluid can be obtained by combining the relationship between the steady flow velocity equation of a Bingham fluid [25] and the steady solution of a Newtonian runny nose (equation (7)).
(21)ubiny,t=unewy,t⋅ftrans1y,(22)ftrans1y=1−2y/b2−2b0/b1−2y/b1−2y/b2,b02<y<b2.Through the unsteady solution of a Newtonian fluid and a Bingham fluid, we can get the fluid flow process in a parallel plate, which will go through the entrance stage, transition stage, and stable stage. The entrance stage is located at the entrance of the parallel plate, which is generally very short in length. In this stage, the boundary effect has not affected the whole flow layer, even for a Newtonian fluid; its velocity distribution will show a similar shape of “flow core.” In the transition stage, the boundary effect has affected the whole flow layer, but the influence of inertial force has not disappeared, and the velocity distribution is still changing. When it comes to the stable stage, the influence of the inertial term has disappeared, and the fluid enters a stable flow state (see Figure3).Figure 3
Distribution diagram of velocity in the entrance stage, the transition stage, and the steady stage when a Newtonian fluid and a Bingham fluid flow in fractures.


## 3. Model Application and Simplification

### 3.1. Model Application

The derivation of the unsteady flow solution of a Bingham fluid is not only to understand the flow pattern of the slurry in the fissure, but more importantly to obtain the parameters or methods that are of great significance to grouting engineering.According tob0=−2τ0L/Δp, when b0=b, the slurry will stop flowing. At this time, the critical pressure gradient is the starting pressure gradient of the slurry, and the critical length is the farthest distance that the slurry can reach in the fracture.
(23)ΔpLcr=−2τ0b,(24)Lmax=−Δpb2τ0.In addition, through further processing of the unsteady-state equation, the corresponding relationship between grouting time and slurry flow distance can be obtained. Also based on Newtonian fluids, equation (18) is integrated to obtain the average velocity:
(25)u¯new=1b∫−b/2b/2uy,tdy=2b2ΔPμI∑k=1∞1kπ4−1k+1+12fTt.Sinceu¯=dI/dt [13], replace the left side of equation (25), which can be obtained after an integral operation:
(26)Inew2=4b2ΔPμ∑k=1∞1kπ4−1k+1+12FTt,(27)FTt=t+bkπ2ρμe−kπ/b2μ/ρt−bkπ2ρμ.According to the previous assumption, the average velocity of a Bingham fluid in the fracture is as follows:(28)u¯bin=ftrans2u¯new,ftrans2=1−b0b21+12b0b.The relationship between time and distance of a Bingham flow in fractures is as follows:(29)Ibin2=ftrans24b2ΔPμ∑k=1∞1kπ4−1k+1+12FTt.In the differential calculation, approximate processing is also adopted. When integrating the flow distanceI, ignoring b0 will also change with I. Therefore, in the actual calculation, it is necessary to correct the value of b0 by using an iterative method.Although the relationship between the slurry flow time and the distance in the fracture has been obtained, it is necessary to distinguish fractures before calculation. Here, fractures are divided into dry fracture and water-saturated fracture according to whether there is water in the fracture before grouting. Different types of fractures use slightly different equations to predict slurry diffusion.First of all, for the dry fracture, it is considered that there is only air inside the fracture before grouting, and the air has high permeability and compressibility, so the influence of air on the slurry flow is not considered here, that is, the gas-liquid two-phase flow is not considered. The traditional method is to use the steady-state flow equation of slurry for iterative calculation, but at the initial stage of slurry diffusion, the diffusion distance tends to zero, that is, the pressure gradient tends to infinity (see Figure4), and the calculated velocity is also very large, so it is difficult to choose the appropriate initial iteration value. One method is to use reverse iteration [26]. Take cubic law as an example, the iteration formula is as follows:
(30)In=In+1−u¯n+1tn+1−tn,u¯n+1=ΔPb212μIn+1.Figure 4
Pressure distribution diagram during grouting in a dry fracture.When the fracture water is saturated, it is assumed that both the slurry and the water are incompressible. Since the pressure ladder will not decrease with the increase of the diffusion distance (see Figure5), the diffusion distance of the slurry can be directly calculated as follows:
(31)Ibin=∫0tu¯binτdτ=ftrans22b2ΔPμL∑k=1∞1kπ4−1k+1+12FTt.Figure 5
Schematic diagram of pressure distribution when grouting in water-saturated fractures.


### 3.2. Simplified Conditions

It can be seen from equation (19) that the attenuation coefficient in the time-dependent attenuation term is inversely proportional to the square of fracture opening, inversely proportional to the density of slurry, and directly proportional to the viscosity of slurry. The development trend of grouting materials is that the viscosity of slurry is decreasing and the density is increasing, which leads to an increase in the time it takes for the slurry to reach steady flow. Especially for the fractures with a large opening, the influence of inertia force cannot be ignored, and the grouting time and flow distance calculated by the steady-state formula are no longer applicable. According to the formula calculation, for the centimeter fracture, it takes nearly 60 s for the slurry to reach the stable state. In the process of coal mining, the grouting time often used for roadway reinforcement is about 200-300 s, which also shows that the influence of inertia force cannot be ignored.When using the Bingham fluid unsteady flow equation to calculate the flow distance of slurry in the fracture, the relative error is ensured to be less than 5%, and the following simplified condition is obtained:(32)μtρb2≥0.3035.The time-dependent terms of the non-steady-state model quickly decay to 0, and they can be replaced by those of the steady-state model:(33)τ0t−3μΔP≤0.062.The size of the “flow core” is negligible with respect to the fracture opening, so the Newtonian fluid model can be used instead of the Bingham fluid model.


## 3.1. Model Application

The derivation of the unsteady flow solution of a Bingham fluid is not only to understand the flow pattern of the slurry in the fissure, but more importantly to obtain the parameters or methods that are of great significance to grouting engineering.According tob0=−2τ0L/Δp, when b0=b, the slurry will stop flowing. At this time, the critical pressure gradient is the starting pressure gradient of the slurry, and the critical length is the farthest distance that the slurry can reach in the fracture.
(23)ΔpLcr=−2τ0b,(24)Lmax=−Δpb2τ0.In addition, through further processing of the unsteady-state equation, the corresponding relationship between grouting time and slurry flow distance can be obtained. Also based on Newtonian fluids, equation (18) is integrated to obtain the average velocity:
(25)u¯new=1b∫−b/2b/2uy,tdy=2b2ΔPμI∑k=1∞1kπ4−1k+1+12fTt.Sinceu¯=dI/dt [13], replace the left side of equation (25), which can be obtained after an integral operation:
(26)Inew2=4b2ΔPμ∑k=1∞1kπ4−1k+1+12FTt,(27)FTt=t+bkπ2ρμe−kπ/b2μ/ρt−bkπ2ρμ.According to the previous assumption, the average velocity of a Bingham fluid in the fracture is as follows:(28)u¯bin=ftrans2u¯new,ftrans2=1−b0b21+12b0b.The relationship between time and distance of a Bingham flow in fractures is as follows:(29)Ibin2=ftrans24b2ΔPμ∑k=1∞1kπ4−1k+1+12FTt.In the differential calculation, approximate processing is also adopted. When integrating the flow distanceI, ignoring b0 will also change with I. Therefore, in the actual calculation, it is necessary to correct the value of b0 by using an iterative method.Although the relationship between the slurry flow time and the distance in the fracture has been obtained, it is necessary to distinguish fractures before calculation. Here, fractures are divided into dry fracture and water-saturated fracture according to whether there is water in the fracture before grouting. Different types of fractures use slightly different equations to predict slurry diffusion.First of all, for the dry fracture, it is considered that there is only air inside the fracture before grouting, and the air has high permeability and compressibility, so the influence of air on the slurry flow is not considered here, that is, the gas-liquid two-phase flow is not considered. The traditional method is to use the steady-state flow equation of slurry for iterative calculation, but at the initial stage of slurry diffusion, the diffusion distance tends to zero, that is, the pressure gradient tends to infinity (see Figure4), and the calculated velocity is also very large, so it is difficult to choose the appropriate initial iteration value. One method is to use reverse iteration [26]. Take cubic law as an example, the iteration formula is as follows:
(30)In=In+1−u¯n+1tn+1−tn,u¯n+1=ΔPb212μIn+1.Figure 4
Pressure distribution diagram during grouting in a dry fracture.When the fracture water is saturated, it is assumed that both the slurry and the water are incompressible. Since the pressure ladder will not decrease with the increase of the diffusion distance (see Figure5), the diffusion distance of the slurry can be directly calculated as follows:
(31)Ibin=∫0tu¯binτdτ=ftrans22b2ΔPμL∑k=1∞1kπ4−1k+1+12FTt.Figure 5
Schematic diagram of pressure distribution when grouting in water-saturated fractures.


## 3.2. Simplified Conditions

It can be seen from equation (19) that the attenuation coefficient in the time-dependent attenuation term is inversely proportional to the square of fracture opening, inversely proportional to the density of slurry, and directly proportional to the viscosity of slurry. The development trend of grouting materials is that the viscosity of slurry is decreasing and the density is increasing, which leads to an increase in the time it takes for the slurry to reach steady flow. Especially for the fractures with a large opening, the influence of inertia force cannot be ignored, and the grouting time and flow distance calculated by the steady-state formula are no longer applicable. According to the formula calculation, for the centimeter fracture, it takes nearly 60 s for the slurry to reach the stable state. In the process of coal mining, the grouting time often used for roadway reinforcement is about 200-300 s, which also shows that the influence of inertia force cannot be ignored.When using the Bingham fluid unsteady flow equation to calculate the flow distance of slurry in the fracture, the relative error is ensured to be less than 5%, and the following simplified condition is obtained:(32)μtρb2≥0.3035.The time-dependent terms of the non-steady-state model quickly decay to 0, and they can be replaced by those of the steady-state model:(33)τ0t−3μΔP≤0.062.The size of the “flow core” is negligible with respect to the fracture opening, so the Newtonian fluid model can be used instead of the Bingham fluid model.


## 4. Results and Verification

### 4.1. Numerical Simulation of Bingham Fluid Flow in Fracture

The numerical calculation results are compared with the calculation results of equation (21). The fracture and slurry parameters selected here [27] are shown in Table 1.Table 1
Fracture parameters and slurry parameters used in calculation.
Injection pressure (kPa)Length (m)Aperture (m)Density (kg/m3)Plastic viscosity (Pa·s)Yield stress (Pa)100.10.00115400.0140.94According to the calculation results (see Figure6), firstly, the approximate solution of the unsteady flow of a Bingham fluid between parallel plates obtained in this paper is highly consistent with the numerical results. Secondly, it also shows that the three stages of velocity distribution mentioned above do exist. It needs to be explained that in the stable stage, the velocity distribution does not observe “flow core,” because the size of “flow core” (b0=1.88e−5m) is too small relative to the fracture opening.Figure 6
Relationship between velocity distribution curve and time when slurry flows in fractures.
(a)
t=0.001s(b)
t=0.002s(c)
t=0.005s(d)
t=0.01s(e)
t=0.02s(f)
t=0.1s


### 4.2. Related Tests

#### 4.2.1. Rheological Property Test of Cement Slurry

In order to study the flow law of a slurry, the first step is to study the flow characteristics of the slurry itself, that is, the rheology of the slurry [28–30]. The commonly used rheological models of slurry include the Newtonian fluid model, the Bingham fluid model, the power law fluid model, and the H-B fluid model. In this paper, a ZNN-D6 rotational viscometer (see Figure 7) is used to test the rheological property of ordinary Portland cement (OPC) slurry, and the rheological curves of the OPC slurry under different water-cement ratios are obtained (see Figure 8). When the water-cement ratio (W/C) of the cement slurry is more than 0.7 but less than 1.0, the rheological curve is basically straight and the intercept between the curve and coordinate Y axis is obviously not zero, that is to say, there is yield stress. Therefore, we think that the water-cement slurry is more suitable for the Bingham fluid model than the OPC slurry with W/C of 0.7~1.0, which is also consistent with the research of other scholars [31]. In the subsequent calculation of the slurry flow distance under different conditions, the slurry parameters of W/C=0.8 in this experiment are selected.Figure 7
ZNN-D6 viscometer physical figure and schematic diagram of key parts.Figure 8
Rheological curves of ordinary Portland cement slurry with different water-cement ratios.


#### 4.2.2. Experiment of Unsteady Fracture Flow of Water

In the derivation of the relationship between the diffusion distance and time of slurry, the differential relationship ofu¯=dI/dt is adopted, so in order to verify equation (26), the experiment of water flowing in the parallel plate fracture is designed. After processing the experimental data, the results obtained and the theoretical calculation results are drawn in Figure 9, which prove the correctness of the differential relationship. At the same time, the correctness of equation (29) can also be guaranteed.Figure 9
The relationship between flow time and flow distance when water flows in a fracture. Model 1 is the result of the iterative calculation of equation (8), and model 2 is the result of direct calculation using equation (26).


### 4.3. Application of Simplified Conditions

Based on formulas (23) and (24), we design two discriminant conditionsXc and Yc:
(34)Xc=μtρb2−0.3035,Yc=τ0t−3μΔP−0.062.According to the discriminant conditions, the discriminant space corresponding to the model can be divided into four quadrants (see Figure10). According to the previous analysis, the entire space can be calculated using model 4, but when the discriminant conditions (Xc, Yc) fall in quadrant II, only model 4 can be used for calculation (to ensure that the model error is within 5%); when the discriminant conditions (Xc, Yc) fall into quadrant IV, all four models can guarantee relatively high accuracy; when the discriminant conditions (Xc, Yc) fall in quadrant II, model 2 can be used for approximate calculation; finally, when the judgment condition (Xc, Yc) is in quadrant I, model 3 can be used for approximate calculation.Figure 10
Schematic diagram of discrimination space.In order to further verify the feasibility of the criteria, using the experimentally obtained cement slurry withW/C=0.8, plastic viscosity is 0.019, yield stress is 2.18, and density is 1600. The slurry flow distance under different conditions is calculated. Among them, model 1 represents the steady-state Newtonian fluid (equation (8)), model 2 represents the unsteady-state Newtonian fluid (equation (26)), model 3 refers to the steady-state Bingham fluid (equation (9)), and model 4 represents the unsteady-state Bingham fluid (equation (29). The calculation results are shown in Table 2.Table 2
Calculation results of four models for slurry flow distance under different conditions.
Aperture,b (mm)Pressure,P (MPa)Time,t (s)Diffusion distance,I (m)Model 1Model 2Model 3Model 40.10.10.010.00940.00930.00930.00930.10.110.09370.09370.09080.09080.10.11000.93660.93660.69900.69900.110.010.02960.02950.02960.02950.1110.29620.29620.29330.29330.111002.96172.96172.68992.68990.1100.010.09370.09330.09360.09320.11010.93660.93650.93370.93370.1101009.36599.36589.08369.083510.10.010.09370.06060.09340.060510.110.93660.93260.90840.904610.11009.36599.36556.9906.9898110.010.29620.19160.29590.19141112.96172.94922.93322.92091110029.617429.616226.899526.89841100.010.93660.60580.93630.605711019.36599.32639.33729.297911010093.658693.654690.835590.8318100.10.010.93660.07740.93370.0774100.119.36596.05789.08365.9391100.110093.658693.263369.900169.68491010.012.96170.24472.95890.2447101129.617419.156629.332019.0370101100296.1744294.9246268.9949267.967310100.019.36590.77399.36300.77391010193.658660.578593.372260.45861010100936.5858932.6333908.3553904.6386Taking model 4 as the standard, first judge its applicability based on the relative errors of the other three models (see Table3), and then calculate the discriminant coordinates under different conditions to get Table 4 through the comparison table. The results in Tables 3 and 4 verify the accuracy of the discrimination method.Table 3
Applicability of the first three models.
Aperture,b (mm)Pressure,P (MPa)Time,t (s)Relative errorModel 1Model 2Model 30.10.10.01√√√0.10.11√√√0.10.1100××√0.110.01√√√0.111√√√0.11100××√0.1100.01√√√0.1101√√√0.110100√√√10.10.01×√×10.11√√√10.1100××√110.01×√×111√√√11100××√1100.01×√×1101√√√110100√√√100.10.01×√×100.11×√×100.1100××√1010.01×√×1011×√×101100××√10100.01×√×10101×√×1010100√√√Note: if the relative error is less than 5%, it is√, otherwise it is ×.Table 4
Discriminant coordinates obtained under different conditions and their quadrants in the discriminant plane.
Aperture,b (mm)Pressure,P (MPa)Time,t (s)Discriminant coordinatesQuadrant0.10.10.01(11.5715, -0.0591)IV0.10.11(1187.2, -0.0331)IV0.10.1100(118750, 0.2267)I0.110.01(11.5715, -0.0611)IV0.111(1187.2, -0.0529)IV0.11100(118750, 0.0293)I0.1100.01(11.5715, -0.0617)IV0.1101(1187.2, -0.0591)IV0.110100(118750, -0.0331)IV10.10.01(-0.1847, -0.0591)III10.11(11.5715, -0.0331)IV10.1100(1187.2, 0.2267)I110.01(-0.1847, -0.0611)III111(11.5715, -0.0529)IV11100(1187.2, 0.0293)I1100.01(-0.1847, -0.0617)III1101(11.5715, -0.0591)IV110100(1187.2, -0.0331)IV100.10.01(-0.3023, -0.0591)III100.11(-0.18475, -0.0331)III100.1100(11.5715, 0.2267)I1010.01(-0.3023, -0.0611)III1011(-0.18475, -0.0529)III101100(11.5715, 0.0293)I10100.01(-0.3023, -0.0617)III10101(-0.18475, -0.0591)III1010100(11.5715, -0.0331)IV


## 4.1. Numerical Simulation of Bingham Fluid Flow in Fracture

The numerical calculation results are compared with the calculation results of equation (21). The fracture and slurry parameters selected here [27] are shown in Table 1.Table 1
Fracture parameters and slurry parameters used in calculation.
Injection pressure (kPa)Length (m)Aperture (m)Density (kg/m3)Plastic viscosity (Pa·s)Yield stress (Pa)100.10.00115400.0140.94According to the calculation results (see Figure6), firstly, the approximate solution of the unsteady flow of a Bingham fluid between parallel plates obtained in this paper is highly consistent with the numerical results. Secondly, it also shows that the three stages of velocity distribution mentioned above do exist. It needs to be explained that in the stable stage, the velocity distribution does not observe “flow core,” because the size of “flow core” (b0=1.88e−5m) is too small relative to the fracture opening.Figure 6
Relationship between velocity distribution curve and time when slurry flows in fractures.
(a)
t=0.001s(b)
t=0.002s(c)
t=0.005s(d)
t=0.01s(e)
t=0.02s(f)
t=0.1s


## 4.2. Related Tests

### 4.2.1. Rheological Property Test of Cement Slurry

In order to study the flow law of a slurry, the first step is to study the flow characteristics of the slurry itself, that is, the rheology of the slurry [28–30]. The commonly used rheological models of slurry include the Newtonian fluid model, the Bingham fluid model, the power law fluid model, and the H-B fluid model. In this paper, a ZNN-D6 rotational viscometer (see Figure 7) is used to test the rheological property of ordinary Portland cement (OPC) slurry, and the rheological curves of the OPC slurry under different water-cement ratios are obtained (see Figure 8). When the water-cement ratio (W/C) of the cement slurry is more than 0.7 but less than 1.0, the rheological curve is basically straight and the intercept between the curve and coordinate Y axis is obviously not zero, that is to say, there is yield stress. Therefore, we think that the water-cement slurry is more suitable for the Bingham fluid model than the OPC slurry with W/C of 0.7~1.0, which is also consistent with the research of other scholars [31]. In the subsequent calculation of the slurry flow distance under different conditions, the slurry parameters of W/C=0.8 in this experiment are selected.Figure 7
ZNN-D6 viscometer physical figure and schematic diagram of key parts.Figure 8
Rheological curves of ordinary Portland cement slurry with different water-cement ratios.


### 4.2.2. Experiment of Unsteady Fracture Flow of Water

In the derivation of the relationship between the diffusion distance and time of slurry, the differential relationship ofu¯=dI/dt is adopted, so in order to verify equation (26), the experiment of water flowing in the parallel plate fracture is designed. After processing the experimental data, the results obtained and the theoretical calculation results are drawn in Figure 9, which prove the correctness of the differential relationship. At the same time, the correctness of equation (29) can also be guaranteed.Figure 9
The relationship between flow time and flow distance when water flows in a fracture. Model 1 is the result of the iterative calculation of equation (8), and model 2 is the result of direct calculation using equation (26).


## 4.2.1. Rheological Property Test of Cement Slurry

In order to study the flow law of a slurry, the first step is to study the flow characteristics of the slurry itself, that is, the rheology of the slurry [28–30]. The commonly used rheological models of slurry include the Newtonian fluid model, the Bingham fluid model, the power law fluid model, and the H-B fluid model. In this paper, a ZNN-D6 rotational viscometer (see Figure 7) is used to test the rheological property of ordinary Portland cement (OPC) slurry, and the rheological curves of the OPC slurry under different water-cement ratios are obtained (see Figure 8). When the water-cement ratio (W/C) of the cement slurry is more than 0.7 but less than 1.0, the rheological curve is basically straight and the intercept between the curve and coordinate Y axis is obviously not zero, that is to say, there is yield stress. Therefore, we think that the water-cement slurry is more suitable for the Bingham fluid model than the OPC slurry with W/C of 0.7~1.0, which is also consistent with the research of other scholars [31]. In the subsequent calculation of the slurry flow distance under different conditions, the slurry parameters of W/C=0.8 in this experiment are selected.Figure 7
ZNN-D6 viscometer physical figure and schematic diagram of key parts.Figure 8
Rheological curves of ordinary Portland cement slurry with different water-cement ratios.


## 4.2.2. Experiment of Unsteady Fracture Flow of Water

In the derivation of the relationship between the diffusion distance and time of slurry, the differential relationship ofu¯=dI/dt is adopted, so in order to verify equation (26), the experiment of water flowing in the parallel plate fracture is designed. After processing the experimental data, the results obtained and the theoretical calculation results are drawn in Figure 9, which prove the correctness of the differential relationship. At the same time, the correctness of equation (29) can also be guaranteed.Figure 9
The relationship between flow time and flow distance when water flows in a fracture. Model 1 is the result of the iterative calculation of equation (8), and model 2 is the result of direct calculation using equation (26).


## 4.3. Application of Simplified Conditions

Based on formulas (23) and (24), we design two discriminant conditionsXc and Yc:
(34)Xc=μtρb2−0.3035,Yc=τ0t−3μΔP−0.062.According to the discriminant conditions, the discriminant space corresponding to the model can be divided into four quadrants (see Figure10). According to the previous analysis, the entire space can be calculated using model 4, but when the discriminant conditions (Xc, Yc) fall in quadrant II, only model 4 can be used for calculation (to ensure that the model error is within 5%); when the discriminant conditions (Xc, Yc) fall into quadrant IV, all four models can guarantee relatively high accuracy; when the discriminant conditions (Xc, Yc) fall in quadrant II, model 2 can be used for approximate calculation; finally, when the judgment condition (Xc, Yc) is in quadrant I, model 3 can be used for approximate calculation.Figure 10
Schematic diagram of discrimination space.In order to further verify the feasibility of the criteria, using the experimentally obtained cement slurry withW/C=0.8, plastic viscosity is 0.019, yield stress is 2.18, and density is 1600. The slurry flow distance under different conditions is calculated. Among them, model 1 represents the steady-state Newtonian fluid (equation (8)), model 2 represents the unsteady-state Newtonian fluid (equation (26)), model 3 refers to the steady-state Bingham fluid (equation (9)), and model 4 represents the unsteady-state Bingham fluid (equation (29). The calculation results are shown in Table 2.Table 2
Calculation results of four models for slurry flow distance under different conditions.
Aperture,b (mm)Pressure,P (MPa)Time,t (s)Diffusion distance,I (m)Model 1Model 2Model 3Model 40.10.10.010.00940.00930.00930.00930.10.110.09370.09370.09080.09080.10.11000.93660.93660.69900.69900.110.010.02960.02950.02960.02950.1110.29620.29620.29330.29330.111002.96172.96172.68992.68990.1100.010.09370.09330.09360.09320.11010.93660.93650.93370.93370.1101009.36599.36589.08369.083510.10.010.09370.06060.09340.060510.110.93660.93260.90840.904610.11009.36599.36556.9906.9898110.010.29620.19160.29590.19141112.96172.94922.93322.92091110029.617429.616226.899526.89841100.010.93660.60580.93630.605711019.36599.32639.33729.297911010093.658693.654690.835590.8318100.10.010.93660.07740.93370.0774100.119.36596.05789.08365.9391100.110093.658693.263369.900169.68491010.012.96170.24472.95890.2447101129.617419.156629.332019.0370101100296.1744294.9246268.9949267.967310100.019.36590.77399.36300.77391010193.658660.578593.372260.45861010100936.5858932.6333908.3553904.6386Taking model 4 as the standard, first judge its applicability based on the relative errors of the other three models (see Table3), and then calculate the discriminant coordinates under different conditions to get Table 4 through the comparison table. The results in Tables 3 and 4 verify the accuracy of the discrimination method.Table 3
Applicability of the first three models.
Aperture,b (mm)Pressure,P (MPa)Time,t (s)Relative errorModel 1Model 2Model 30.10.10.01√√√0.10.11√√√0.10.1100××√0.110.01√√√0.111√√√0.11100××√0.1100.01√√√0.1101√√√0.110100√√√10.10.01×√×10.11√√√10.1100××√110.01×√×111√√√11100××√1100.01×√×1101√√√110100√√√100.10.01×√×100.11×√×100.1100××√1010.01×√×1011×√×101100××√10100.01×√×10101×√×1010100√√√Note: if the relative error is less than 5%, it is√, otherwise it is ×.Table 4
Discriminant coordinates obtained under different conditions and their quadrants in the discriminant plane.
Aperture,b (mm)Pressure,P (MPa)Time,t (s)Discriminant coordinatesQuadrant0.10.10.01(11.5715, -0.0591)IV0.10.11(1187.2, -0.0331)IV0.10.1100(118750, 0.2267)I0.110.01(11.5715, -0.0611)IV0.111(1187.2, -0.0529)IV0.11100(118750, 0.0293)I0.1100.01(11.5715, -0.0617)IV0.1101(1187.2, -0.0591)IV0.110100(118750, -0.0331)IV10.10.01(-0.1847, -0.0591)III10.11(11.5715, -0.0331)IV10.1100(1187.2, 0.2267)I110.01(-0.1847, -0.0611)III111(11.5715, -0.0529)IV11100(1187.2, 0.0293)I1100.01(-0.1847, -0.0617)III1101(11.5715, -0.0591)IV110100(1187.2, -0.0331)IV100.10.01(-0.3023, -0.0591)III100.11(-0.18475, -0.0331)III100.1100(11.5715, 0.2267)I1010.01(-0.3023, -0.0611)III1011(-0.18475, -0.0529)III101100(11.5715, 0.0293)I10100.01(-0.3023, -0.0617)III10101(-0.18475, -0.0591)III1010100(11.5715, -0.0331)IV


## 5. Conclusion

The unsteady flow equation of a Newtonian fluid in a parallel plate fracture is derived in detail based on the N-S equation. Combining the steady-state flow equations of a Newtonian fluid and a Bingham fluid, the unsteady-state flow equation of a Bingham fluid which is more suitable for describing slurry flow is constructed. The application conditions of the new model are analyzed. We can get the following conclusions:(1)
We have deduced the relationship between the flow time and flow distance of a Bingham fluid in the parallel plate fractures, and proved the accuracy of the derivation process through experiments(2)
The unsteady flow equation of a Bingham fluid makes up for the traditional model’s inability to consider the influence of inertial forces on the calculation results, and the slurry flow velocity distribution is divided into the inlet section, the transition section, and the stable section through calculation(3)
We have provided the judgment condition that the unsteady flow equation of a Bingham fluid simplifies the traditional model. After bringing in the calculation criterion of specific construction parameters, we choose whether to use the traditional model and which model to useThe flow equations of power law fluids and H-B fluids can be calculated with reference to the proposed ideas, which can be used for the flow calculation of a high water-cement ratio slurry. In future research, the effect of fracture roughness and tortuosity on slurry flow should also be studied. At the same time, changes in the rheology of the cement slurry over time should be considered.



---

*Source: 6671982-2021-03-28.xml*