Simulation of Particle-Fluid Interaction in Fractal Fractures Based on the Immersed Boundary-Lattice Boltzmann Method

Liu, Shenggui; Tang, Songlei; Huang, Jinkuang; Lv, Mindong; Li, Yingjun

doi:https://doi.org/10.1155/2020/6695623

Geofluids

On this page

Abstract Introduction Results and Discussion Conclusion Abbreviations Data Availability Conflicts of Interest Acknowledgments References Copyright Related Articles

Special Issue

Modelling and Simulation of Geofluids Transport in Reservoirs

View this Special Issue

Research Article | Open Access

Volume 2020 | Article ID 6695623 | https://doi.org/10.1155/2020/6695623

Simulation of Particle-Fluid Interaction in Fractal Fractures Based on the Immersed Boundary-Lattice Boltzmann Method

Shenggui Liu,¹Songlei Tang,¹Jinkuang Huang,¹Mindong Lv,¹and Yingjun Li²

Academic Editor: Abdelraheem Mahmoud Aly

Received05 Oct 2020

Revised17 Nov 2020

Accepted26 Nov 2020

Published16 Dec 2020

Abstract

In order to reveal the influence of particles on fluid flow characteristics in rough fractures under fluid-solid coupling, a range of fracture systems with varying roughness were generated using the Weierstrass-Mandelbrot function. Fluid-particle interactions in rough fractal fractures were simulated using the immersed boundary-lattice Boltzmann method. In this paper, the effects of fluid viscosity, particle size, particle quantity, fracture fractal dimension, and particle grading composition are studied. Results illustrate that increasing fluid viscosity hinders the movement of particles, resulting in the decreasing of particle velocity. As particle size and particle quantity increase, the particle occupation of the channel area grows larger, which lead to lower permeability of the channel. Increasing fracture fractal dimension surges the curvature of the fluid channel, but permeability has a negative exponential correlation to fractal dimension. With increasing particle grading composition, the blocking effect of particles on fracture flow also increases with increasing particle proportion.

1. Introduction

Fractures are widespread in rock masses, and fluid flow through fractures affects their mechanical properties by inducing deformation, fracturing, softening, argillization, and corrosion [1]. For oil shale and coal reservoirs, grand and undetected fractures may lead to significant engineering problems because they prominently alter groundwater movement [2, 3]. As fractures increase the transport distance of contaminants, a better understanding of fracture geometry and fracture occurrence will contribute to a more accurate determination of water protection areas [4, 5]. Particle migration and the consequent permeability reduction can cause reservoir damage in numerous petroleum, environmental, and water resource technologies [6, 7]. It is of great scientific significance and engineering application value to accurately and quantitatively describe the structure characteristics of fractured rock and evaluate the influence of internal particle migration.

The permeability performance and seepage direction of rock masses are not only related to the expansion and cutting characteristics of the fracture network [8] but are also closely related to the geometric characteristics, aperture, roughness, connectivity, and contact degree of a single fracture [9–11]. Fracture roughness is the essential element that causes deviation between the actual fracture permeability and theoretical permeability.

Sedimentary transport and adhesion plugging in the fracture space considerably change the permeability characteristics of fractured rock masses, resulting in a change in the groundwater level and the redistribution of pore water pressure within the formation, eventually inducing geological disasters [12–14]. Presently, the porous media capillary bundle model has been widely studied to explore the sediment transport mechanism [15, 16]. The model simplifies the coal medium into a uniform combination of capillary units, and each unit contains certain amount of capillaries. Happel [17] proposed an ideal model that considers hydrodynamic and gravity effects. The method is based on a granular porous medium model, which assumes that the porous medium is composed of numerous spherical particle elements. Yao [18] proposed an analysis method for determining particle motion paths and determined that the control factors of suspended particle motion are diffusion, sedimentation, and interception. For the sake of simplicity, the coarse-grained model ignores the dynamics of fluids and the Van Der Waals forces between molecules. Tufenkji and Elimelech [19] revised the governing equations of diffusion, sedimentation, and interception in the particle motion path analysis method.

James and Chrysiko [20] believe that suspended particles migrate at high speed with central fluid in the pores. Raizah and Ahmed [21] examined the influence of two cross-section heated ellipses on the fluid flow and heat transfer inside porous enclosures. Tien, Sharma, and Mattson et al. [22, 23] established a 3D fracture simulation system. In this method, the Lagrangian method was used to track solid particles, and the influence of the density, size, and flow rate of solid particles on fluid migration velocity was considered. Oseen [24] modified the resistance expression of Stokes formula by comprehensively examining the influence of the inertia term. In porous media, Civan and Shamar [25] divided particle sedimentation into surface precipitation, pore throat blockage, pore filling, and filter cake formation.

In the past decades, there have been numerous direct numerical simulation methods for particle-fluid coupling, such as the finite element-arbitrary Lagrangian Euler method, fictitious domain method, lattice Boltzmann method, immersed boundary method, and smoothed particle hydrodynamics method [26–31]. Feng et al. [29] proposed the immersed boundary–lattice Boltzmann method (IB-LBM) that couples the computational fluid dynamics and immersion boundary method. In this method, fluid flow is simulated using the Navier-Stokes equation, and the hydrodynamic force applied to the solid particles can be calculated using the structure and velocity of the solid particles. Unlike the lattice Boltzmann method, the IB-LBM method uses a Cartesian grid and a set of Lagrangian points to discretize the boundaries of solid obstacles. Discrete Lagrangian lattice points are closer to the curve boundary than those in a staircase fashion. Simultaneously, the IB-LBM solves the Navier-Stokes equations in a rectangular domain, where flexible boundaries can move or change shapes in a complicated fashion. The boundary force is applied to the fluid using a chosen form of the Dirac function. The IB-LBM method has a higher efficiency than body-fitted methods, e.g., finite-element method [32].

In order to understand the flow mechanism through microchannels with varying morphological structures and the influence of sedimentary particles during drainage, the IB-LBM is utilized to analyze and calculate the fluid seepage behavior in a rock mass with a two-dimensional fracture structure in this contribution. Using the Weierstrass-Mandelbrot function, a series of fractal curves with varying roughness and microscale rock fracture structure model are established.

Through the numerical simulations of the steady-state seepage within the fractured structure, the velocity field distribution and the equivalent permeability coefficient of the fluid are calculated. Through the analysis of simulation results, the relationship between the fracture fractal dimension and equivalent permeability coefficient is established. The seepage behavior of fluid in complex fractures with particles is explored, and the influence of particle size, particles quantities, fracture fractal dimension, and particle grading composition on fluid passage is considered.

2. Numerical Methods

2.1. Lattice Boltzmann Method

The Navier-Stokes equations are the fundamental partial differential equations that describe incompressible fluid flow. Using the rate of stress and rate of strain tensors, the velocity-pressure formulas can be written as where is the fluid velocity, is the time, is the kinematic viscosity, is the pressure, is the fluid density, and is the forcing term. In the IB-LB method, the effect of a solid boundary pushing or blocking a fluid is convert into a force term, which is the most essential component of . The other part of is the conservative force, such as gravity or a magnetic force.

The fluid velocity must satisfy the Mach number in order to meet the condition of an incompressible fluid, where is the fluid velocity and is the lattice speed of sound.

In the LB method, the D2Q9 model is the commonly used discrete format for solving Eqs. (1) and (2). The fluid flow can be determined by resolving the particle collision and streaming processes, and the lattice Boltzmann equation is used to determine the streaming and collision processes of fluid particles. In this model, velocity is microscopically discretized into 9 directions and the discrete velocity (Figure 1). The density distribution function is the only quantity of concern, and the discrete lattice Boltzmann equation can be written as where is the density distribution function of moving populations. Here, is the dimensionless mean relaxation coefficient, and is the time step. The equilibrium distribution function is obtained using a Taylor series expansion using the Maxwell-Boltzmann distribution function with velocity up to the second order, which can be written as where the is the weight coefficient associated with velocity , and the nondimensional lattice speed of sound is model-dependent. In the D2Q9 model, commonly equals and .

Applying the Chapman-Enskog multiscale expansion method, the Navier-Stokes equations can be recovered using the LB equation with first-order time precision and second-order space precision. Density and local momentum are determined using the summation and the first order moment of the density distribution function:

Performing the Chapman-Enskog expansion on the LB dynamics with a low Mach number limit gives the following expression for kinematic viscosity: where is the lattice spacing.

In order to implement IB in the LB method, Eq. (3) can be separated into collision and streaming components. Thus, in the IB-LB method, the collision and streaming components are

As mentioned above, is the discrete force term. In this paper, the discrete method for the force term proposed by Guo et al. [33] is used, and the discrete force term is

However, the existence of the external force term changes the density and momentum; hence, a velocity correction is needed:

There are numerous methods for the dispersion of force terms, but the accuracy of this method is relatively high, and a second-order precision in space can be achieved. In addition, the implementation of the method is relatively simple.

2.2. Immersed Boundary Method

Coupling LB and IB requires a procedure to implement the effect of the body force between the Eulerian grids and Lagrangian grids. Currently, there are two methods to couple force terms, the direct force method (DF) [34] and interpolation force method (IF) [35]. The DF method aims at directly obtaining the force term on the Euler grid using the Dirac function under the condition that the solid boundary stress is already known. However, the IF method aims to determine the velocity and force terms on the Euler grid through high-precision interpolation. Dupuis et al. compare the differences between the two methods in detail and found that the difference between the two methods was not notable [36]. The DF is adopted in this paper to exert the force term. The computational scheme for the implementation of the IB-LBM in the particle-fluid interaction problem can be found in the paper by Liu et al. [37].

In order to model the physical surface of a particle, we consider a system with physical particles with discrete points. We assume that the reference point of the particle is stationary, and all boundary points are on a circle. At time , the center of the particle is at , and the instantaneous particle rotational matrix is . The position of a reference point in a Lagrangian coordinate system with respect to the particle may be determined using the following formula: where is the boundary reference point of the particles at timestep . The linear restoring force will generate when the reference point and the boundary point are not in the same position. We suppose that the particle deviation is , where the subscripts indicate the th particle and the th discrete point. The restoration force can be calculated as [37] where is the spring constant.

The internal link force in the simulation of the rigid particle motion is simply the balance force that was mentioned above. Hence, the total force and the total torque exerted on the ith particle are where is the mass of the ith particle.

In order to satisfy the no-slip boundary condition on the particle-fluid interface, the velocity on the particle point must strictly be the same as on its neighboring fluid. We can obtain the velocity of the Lagrangian points from the velocity of the adjacent fluid. where is the displacement of solid particles. We can also use the Dirac function to composite the fluid and particle:

When the IB-LBM is used to simulate particle motions in fluid, the include external force (such as buoyancy), internal force (such as the particle collision force), and restoration force can be calculated. Hence, the total force exert on the particle surface node can be written as [37]

A collision technique is applied if the gap between the particles is less than a given threshold. The functional form of the repulsive force is given by Glowinski et al. [38]:

The force density is obtained from the immersed boundary force density through integration over the immersed boundary.

The Dirac function is a continuous interpolation function. We cannot interpolate all the points one by one in the interpolation process; therefore, the Dirac function must be approximated at the grid points. Numerous methods exist for approximating the Dirac function [39]. We use a more continuous method:

2.3. Characterization of Fracture Roughness

The rock damage process is actually the appearance, growth, bifurcation, and intersection of microcracks in the rock to form a macroscopic fracture. The seepage structure of fractures commonly displays fractal characteristics [40–42].

The three-dimensional fracture structure network of rock masses has complex structural characteristics composed of multiple intersecting and combined two-dimensional fracture surfaces; however, the two-dimensional fracture surface is composed of multiple groups of single fractures. The structure curve of a single fracture obtained by cutting the cross-section along a specific section presents an irregular fractal structure.

A quantitative description method of single fracture can determine the influence of the rough structure fracture shape on fluid seepage and clarify the fluid-particle coupling mechanism in the rough structure. Firstly, we use the Weierstrass-Mandelbrot function to obtain fractal curves for different fractal dimensions , which is where is a parameter greater than 1, and phase is any angle that makes exhibit deterministic or stochastic behavior. is the Hausdorff fractal dimension of the graph of . The fractal governing function is the real part of .

The fractal dimension is equal to 1.0, 1.1, 1.2, 1.3, 1.4, and 1.5, respectively (Figure 2). During the fracture curve generation process, is fixed at 1.4, and the time control parameter is and . The numerical simulation model is 10 mm long by 4 mm wide, and the average gap width of the fracture is 1 mm.

Judging from the shape of the fractal structure curves, with increasing fractal dimension, the crack surface gradually becomes steeper, and wrinkle degree surges. Conversely, the smaller the fractal roughness index, the smoother the fractal surface.

The influence of roughness on fluid flow with particle-fluid interaction is investigated using a two-dimensional incompressible model. The 2D profile of the fracture surface is imported into the IB-LBM in order to generate a fractured model (Figure 3). The length of the fracture model is set to 10 mm to be consistent with a realistic profile. The aperture between the two parallel rough profiles is 1 mm. The pressure is applied at the inlet and outlet of the fracture, and the bounce-back scheme is used to simulate a nonslip boundary at the fluid–solid interface.

3. Results and Discussion

3.1. Benchmark-Single Particle Sedimentation in Fluid

In order to validate the accuracy of the immersed boundary-lattice Boltzmann method, we adopt a model that considers the sedimentation of a circular particle under gravity and buoyancy in a channel filled with a viscous fluid. In 2006, Wan and Turek [43] simulated this single-channel model using FEM, which considered the initial state of water flow in the pipeline as static. In order to verify the correctness of the model in this paper, we compared the results of the immersion boundary-lattice Boltzmann simulation with the FEM results of Wan and Turek. The dimension of the box is 2 cm (in the -direction) by 6 cm (in the -direction). The fluid domain is divided into square lattices with spacing of . The viscosity and density of the fluid are 0.01 Pa·s and 1 g/cm³, respectively. The radius of the particle is 0.125 cm, the particle density is 1250 kg/m³, and the relaxation time in the simulation is 0.53. The four boundaries of the model are stationary walls, which are no-slip boundaries for the fluid. Initially, one stationary solid particle is generated at the position (1 cm, 4 cm). Owing to gravity, a solid particle will move downward gradually.

During the simulation, the particle position coordinates and flow field velocity distribution were recorded. Figure 4 shows the curve of particle velocity in three cases with varying viscosities as time increases. Particle velocity increases continuously in the initial state due to gravity, and at about 0.3 s, the speed reaches a steady state (Figure 4). The particles produce a vertical upward water resistance during acceleration, and the resistance is positively correlated to the velocity. When the velocity reaches a certain value, the gravitational force and the resistance force on the particle remain the same, and the particle reaches the ultimate velocity. The ultimate velocity of the particle will gradually decrease as viscosity increases because the viscosity of the fluid increases the viscous resistance of the particles (Figure 4). Our simulation results are consistent with the FEM results Wan and Turek in 2006, which indicate the correctness of the IB-LBM.

Particle movement and fluid velocity contours at times of 0, 0.2, 0.4, 0.6, 0.8, 0.9, and 1.0 seconds are given in Figure 5. Particle sedimentation impacts the four channel walls. When the particles pass through a certain area, a large impact force is generated on the vertical pipe wall of the area. The particles squeeze the fluid, causing it to spread out on both sides, and the fluid moves in the opposite direction relative to the particle movement. A smaller impact area is formed in front of the particle movement, and a larger vortex area is formed behind the particle. The difference between these two areas will gradually increase as particle speed increases.

3.2. Particle Disturbance in Rough Fractures

This section will focus on the effects of the particulate matter on fluid flow in irregular fractures, including the influence of the particle size on channel permeability, the effect of particle number on fluid flow, the filling behavior of particles in rough fractures with varying fractal dimensions, and the influence of particle gradation on seepage effects.

3.2.1. Effect of the Particle Size on Fluid Flow

Figure 6 shows a schematic diagram of the flow configuration with our coordinate system and the geometrical parameters adopted in this study. In order to examine the effect of the particle size on the seepage in the fractured channel, we set up a sinusoidal uplift wall with a height of , and three circular particles with varying diameters were set up in the trough of the wall.

Refer to the parameters of the benchmark program above, the density and viscosity of the fluid are set to 1000 kg/m³ and 0.01 Pa·s, respectively. The particle density (1250 kg/m³) is slightly heavier than that of the fluid. The inlet and outlet are pressure boundary conditions, the pressure difference is set to 0.01 kPa, and the upper and lower walls are set to complete bound back boundary conditions to ensure that the speed at the side wall is zero.

The spacing of each grid point is 10 μm, the length of the entire field is , the number of grids is , and the relaxation time in the simulation is set to 0.53. In order to examine the effect of the particle diameter on fluid passage, we set the height of the wall protrusion to 0.2 mm. The velocity distribution in the fluid channel under four different cases is analyzed. In the first case, fluid flow through a sinusoidal channel containing no particles was analyzed. The other three cases all examined flow through a sinusoidal channel containing a circular particle with diameters of 0.12, 0.24, and 0.48 mm.

The IB-LBM is used to simulate the fluid flow field and particle motion in this model to obtain the velocity contours of the flow field with particles of different sizes. Figure 6 shows the fluid flow path of different particles in the same channel. When the diameter of the particles is smaller than the height of the corrugated curved channel, there is almost no effect on the seepage velocity in the channel, but as the diameter of the particles increases, the passage of the fluid in the channel is blocked by the particles (Figure 7). The fluid passage deteriorates, thereby increasing the velocity of this part of the fluid and simultaneously reducing the seepage area of the particles as well as the seepage flux. Extending this concept to a complex fracture network indicates that the particle size will hinder the channel only when the particle height reaches the height of the rough protrusions on the fractures, and this hindering effect will increase as particle size increases.

In order to quantitatively determine the influence of particle radius on the channel, the particle drag coefficient and the channel permeability are analyzed. The ratio of the drag force exerted by the fluid on the particle to the product of the projected area of the particle in the direction of movement and the fluid dynamic pressure is where represents the resistance of the particle surface along the incoming flow direction, represents the incident velocity of the fluid, and represents the particle diameter.

Figure 8 shows the drag coefficients of particles with varying diameters. From the stable value of the drag coefficient, we can find that particles with a large radius have greater resistance to the particles because the particles hinder the smooth flow of the fluid and increase the path of movement of the fluid, thereby increasing the velocity of the fluid and further increasing the drag coefficient. The fluid movement in the microscopic seepage channel can also be characterized by permeability, which represents the difficulty for the fluid to pass through complex fractures. Permeability at the microscopic scale is given by where is the flow quantity, which is estimated using a numerical integration algorithm, is the partial path length, and is the cross-sectional area of channel.

Table 1 shows the flow rate and permeability change with particle radius. The existence of large particles changes the dominant path of the flow field, thereby changing the flux and permeability of the flow field (Figure 7). We can make the same observation using the change law of permeability. As particle size increases, the permeability of the flow field decreases. Owing to the obstruction effect of particles, the fluid cannot deliver effectively; so, the pressure drops, and the flux reduces, which presents the decline of permeability.

3.2.2. Effect of the Particle Number on the Flow Channel

In order to study the influence of the number of particles on the fluid movement in complex fractal fractures, we set up a curved and complex fractal fracture keeping the structure of the fractal fracture unchanged. The simulation discussed the influence of the presence of particles on fluid flow when the number of particles is 0, 5, 10, and 15. The fractal dimension of the fracture structure is 1.1, the length of the channel is set to , and the aperture of the fracture is set to 1 mm. The spacing of each grid point is 50 μm; so, the number of grids is . The density and viscosity of the fluid are 1000 kg/m³ and 0.01 Pa·s, respectively. The density of the particles is slightly heavier than that of the fluid (1250 kg/m³), and the radius of the particles is 0.1 mm. The inlet and outlet are pressure boundary conditions, the pressure difference is set to 0.01 kPa, and the upper and lower fracture surfaces are set to bounce back boundary conditions.

Figure 9 shows the seepage process of fractures with a fractal dimension of 1.1 containing different numbers of particles. From numerical analysis of the contour map, it can be inferred that as the number of particles increases, the average seepage speed decreases (Figure 9). The presence of particles decreases the cross-sectional area of the fluid flow and hinders fluid passage. The fluid will form a velocity increase zone around the particles, but this has a lesser impact than the permeability reduction caused by the decrease in the fluid passage cross-sectional area. The more particles, the larger the fluid channel area occupied by the particles, and the stronger the hindrance. As the number of particles increases, the area of the fluid channel occupied by the particles becomes increases, and the blocking effect becomes stronger.

Table 2 shows the seepage velocity, flux, permeability, and permeability reduction rate of complex fractures with varying particle numbers. The presence of more particles will significantly reduce the passage of fluid; that is, reduce the permeability of the fracture space (Table 2). Simultaneously, the permeability reduction coefficient positively correlates with the space occupied by the particles. The larger the space occupied by particles, the greater the permeability reduction coefficient. The fluid-particle interaction flow not only needs to overcome the hindrance of irregular boundaries but also needs to entrain the solid particles to move together, so it causes more energy expenditures. The greater the number of particles, the greater the energy expenditure. This is the main factor of permeability reduction caused by particles.

3.2.3. Particle Motion in Fractures with Varying Fractal Dimensions

Considering the seepage process of complex fractures with varying fractal dimensions, this section constructs six fractal curves with varying fractal dimensions, which correspond to six fracture structures with varying roughness. The fractal dimensions of the six different fractal curves are 1.0, 1.1, 1.2, 1.3, 1.4, and 1.5, and the parameter settings of fluid and particles are exactly the same as those in Section 2.3. Figure 10 shows the contour map of fluid-particle coupled seepage velocity under different fractal dimensions. Fractal dimension has a significant influence on the flow. Increasing the fractal dimension means that the crack roughness increases, the tortuosity of the fluid passing through the crack increases, and the channel permeability decreases.

Figure 11 shows the permeability curves under different fractal dimensions. The larger the fractal dimension of the fracture, the smaller the permeability of the channel. The increase of the fracture fractal dimension means that the structure is more complex, which indicates that the uplift and depression amplitudes are greater. When fluid flows through a fracture with a larger fractal dimension, the boundary impact that needs to be overcome is stronger, resulting in a decrease in permeability. From a numerical point of view, this relationship approximately satisfies the negative exponential correlation between permeability and fractal dimension, which has also been experimentally verified by Jafari and Babadagli [44]. The permeability of complex fractures has a positive correlation with the fractal dimension, and the negative exponential law is also approximately satisfied for fluids containing particles.

3.2.4. Blocking Effect of the Particle Size Distribution

In order to simulate the effect of the particle size distribution on fluid seepage, the following four sets of seepage channel diagrams of size particles with varying proportions are set up (Figure 12). The cavity size of the four models is mm, the radius of the large particles is 9 mm, and the radius of the small particles is 4.5 mm. The large to small particle radius ratio is 2 : 1, and the ratio of the occupied area is 4 : 1. In order to eliminate the error caused by the position of the particles to the fluid channel, we attempt to make the position of the particles as symmetrical and regular as possible, and the area of the total cavity occupied by large and small particles is roughly equal. When the fluid flows through a complex seepage channel, it will autonomously flow through the superior channel; that is, flow in the place where the flow path is the shortest and the seepage width is wide (Figure 12).

Figure 12 shows the velocity contour map with varying the particle size distribution. With increasing amounts of small particles, that is, the ratio of large particles to small particles increases, the fluid seepage channel becomes complex, the advantageous seepage path incorporates more twists and turns, and fluid seepage velocity increases significantly, but the overall permeability decreases. Increasing the particle gradation ratio makes the fluid path more dispersed, which is not enough to form a larger main seepage path and reduces permeability.

4. Conclusion

A series of particle-fluid coupling models are constructed using the IB-LBM and Weierstrass -Mandelbrot function. Numerical simulation results show that the IB-LBM has good numerical stability and convergence in simulating complex fluid seepage flow in fractures with particles, and this method can accurately simulate fluid-particle coupling. In order to explore the seepage behavior of fluid in complex fractures with particles, the effects of the particle size, particle number, fractal dimension, and particle gradation composition on fluid passage are examined. The coupling of fluid particles is analyzed at the microstructure level.

With increasing particle diameter, fluid passage in the channel is blocked by particles, which decreases the seepage area, resulting in a decrease in seepage flux and permeability. With increasing particle number, the average seepage velocity decreases. Fluid velocity will increase around the particles, but this effect is smaller than the overall permeability reduction caused by the decrease in fluid cross-sectional area. With increasing fractal dimension, the roughness of fracture increases, and the complexity of corresponding fluid flow path through the fracture increases. Increasing particle gradation composition makes the fluid channel more dispersed but is not enough to form a large main seepage channel, resulting in a permeability decrease.

Abbreviations

:	Fluid velocity
:	Time
:	Kinematic viscosity
:	Pressure
:	Fluid density
:	Force term
:	Mach number
:	Fluid velocity
:	Lattice speed of sound
:	Discrete velocity
:	Density distribution function
:	Dimensionless mean relaxation coefficient
:	Time step
:	Equilibrium distribution function
:	Discrete force
:	Number of physical particles
:	Number of discrete points in each particle
:	Particle center
:	Instantaneous particle rotational matrix
:	Boundary reference point
:	Particle deviation
:	Restoration force
:	Spring constant
:	Total force
:	Total torque
:	Mass of the ith particle
:	Displacement of solid particles
:	Dirac function
:	Collision force
:	Immersed boundary force density
D:	Fractal dimension
:	Fractal parameter
:	Fractal angle
:	Fractal stochastic index
:	Fractal governing function
:	Resistance coefficient
:	Fluid incident velocity
:	Particle diameter
:	Permeability
:	Flow quantity
:	Partial path length
:	Channel cross-sectional area.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

No potential conflicts of interest were reported by the authors.

Acknowledgments

This project is supported by the National Natural Science Foundation of China (41472130).

References

M. Wang, Y. F. Chen, G. W. Ma, J. Q. Zhou, and C. B. Zhou, “Influence of surface roughness on nonlinear flow behaviors in 3D self-affine rough fractures: lattice Boltzmann simulations,” Advances in Water Resources, vol. 96, pp. 373–388, 2016.
View at: Publisher Site | Google Scholar
M. A. Hariri-Ardebili, S. M. Seyed-Kolbadi, and H. Mirzabozorg, “A smeared crack model for seismic failure analysis of concrete gravity dams considering fracture energy effects,” Structural Engineering and Mechanics, vol. 48, no. 1, pp. 17–39, 2013.
View at: Publisher Site | Google Scholar
Y. Zhang, K. Höeg, W. Wang, and Y. Zhu, “Watertightness, cracking resistance, and self-healing of asphalt concrete used as a water barrier in dams,” Canadian Geotechnical Journal, vol. 50, no. 3, pp. 275–287, 2013.
View at: Publisher Site | Google Scholar
R. E. Jackson, A. W. Gorody, B. Mayer, J. W. Roy, M. C. Ryan, and D. R. Van Stempvoort, “Groundwater protection and unconventional gas extraction: the critical need for field-based hydrogeological research,” Groundwater, vol. 51, no. 4, pp. 488–510, 2013.
View at: Publisher Site | Google Scholar
A. Vengosh, R. B. Jackson, N. Warner, T. H. Darrah, and A. Kondash, “A critical review of the risks to water resources from unconventional shale gas development and hydraulic fracturing in the United States,” Environmental Science & Technology, vol. 48, no. 15, pp. 8334–8348, 2014.
View at: Publisher Site | Google Scholar
M. T. Byrne and S. M. Waggoner, “Fines migration in a high temperature gas reservoir - laboratory simulation and implications for completion design,” in SPE-121897. SPE European formation damage conference held in Scheveningen, the Netherland, 2009.
View at: Google Scholar
F. Civan, Reservoir Formation Damage, Gulf Professional Publishing, Burlington, MA, USA, third ed edition, 2014.
X. Xiong, B. Li, Y. Jiang, T. Koyama, and C. Zhang, “Experimental and numerical study of the geometrical and hydraulic characteristics of a single rock fracture during shear,” International Journal of Rock Mechanics and Mining ences, vol. 48, no. 8, pp. 1292–1302, 2011.
View at: Publisher Site | Google Scholar
T. Petchsingto, Numerical Study of Fracture Aperture Characteristics and their Impact on Single-Phase Flow and Capillary-Dominated Displacement, [Ph.D.thesis], 2008.
A. Kamal, M. Kunieda, N. Ueda, and H. Nakamura, “Evaluation of crack opening performance of a repair material with strain hardening behavior,” Cement & Concrete Composites, vol. 30, no. 10, pp. 863–871, 2008.
View at: Publisher Site | Google Scholar
G. Dimadis, A. Dimadi, and I. Bacasis, “Influence of fracture roughness on aperture fracture surface and in fluid flow on coarse-grained marble, experimental results,” Journal of Geoscience and Environment Protection, vol. 2, no. 5, pp. 59–67, 2014.
View at: Publisher Site | Google Scholar
C. Zou, W. Z. Zhao, C. Z. Jia et al., “Formation and distribution of volcanic hydrocarbon reservoirs in sedimentary basins of China,” Petroleum Exploration and Development, vol. 35, no. 3, pp. 257–271, 2008.
View at: Publisher Site | Google Scholar
J. Becherer, G. Flöser, L. Umlauf, and H. Burchard, “Estuarine circulation versus tidal pumping: sediment transport in a well-mixed tidal inlet,” Journal of Geophysical Research, Oceans, vol. 121, no. 8, pp. 6251–6270, 2016.
View at: Publisher Site | Google Scholar
H. Zhu, D. Wang, and P. Cheng, “Dynamic release process of pollutants during suspended sediment transport in aquatic system,” Appl. Math. Mech.-Engl. Ed., vol. 35, no. 11, pp. 1411–1420, 2014.
View at: Publisher Site | Google Scholar
A. C. Payatakes, C. Tien, and R. M. Turian, “A new model for granular porous media: Part I. model formulation,” American Institute of Chemical Engineers Journal, vol. 19, no. 1, pp. 58–67, 1973.
View at: Publisher Site | Google Scholar
L. N. Reddi and M. V. S. Bonala, “Analytical solution for fine particle accumulation in soil filters,” Journal of Geotechnical And Geoenvironmental Engineering, vol. 123, no. 12, pp. 1143–1152, 1997.
View at: Publisher Site | Google Scholar
J. Happel and H. Brenner, Low Reynolds number hydrodynamics: With special applications to particulate media, Springer Science and Business Media, 2012.
K. Yao, M. T. MJIabibian, and C. R. Omelia, Water and waste through sand [Ph.D. thesis], Department of Civil and Sanitary Engineering, Massachusetts Institute of Technology, Cambridge, 1940.
N. Tufenkji and M. Elimelech, “Deviation from the classical colloid filtration theory in the presence of Repulsive DLVO interactions,” Langmuir, vol. 20, no. 25, pp. 10818–10828, 2004.
View at: Publisher Site | Google Scholar
S. C. James and C. V. Chrysikopoulos, “Effective velocity and effective dispersion coefficient for finite-sized particles flowing in a uniform fracture,” Journal of Colloid and Interface Science, vol. 263, no. 1, pp. 288–295, 2003.
View at: Publisher Site | Google Scholar
Z. A. S. Raizah and S. E. Ahmed, “Convective transport in case of variable properties in porous enclosures with/without two heated ellipsis with rough boundaries,” Alexandria Engineering Journal, vol. 59, no. 5, pp. 3927–3943, 2020.
View at: Publisher Site | Google Scholar
C. Tien and B. V. Ramarao, Granular filtration of aerosols and hydrosols, Elsevier, 2011.
P. Sharma, M. Flury, and E. D. Mattson, “Studying colloid transport in porous media using a geocentrifuge,” Water Resources Research, vol. 44, no. 7, pp. 767-768, 2008.
View at: Publisher Site | Google Scholar
Y. Ghorbani, A. N. Mainza, J. Petersen, M. Becker, J.-P. Franzidis, and J. T. Kalala, “Investigation of particles with high crack density produced by HPGR and its effect on the redistribution of the particle size fraction in heaps,” Minerals Engineering, vol. 43-44, no. 4, pp. 44–51, 2013.
View at: Publisher Site | Google Scholar
M. M. Sharma, H. Chamoun, D. S. H. S. R. Sarma, and R. S. Schechter, “Factors controlling the hydrodynamic detachment of particles from surfaces,” Journal of Colloid and Interface Science, vol. 149, no. 1, pp. 121–134, 1992.
View at: Publisher Site | Google Scholar
S. Ghosh and S. Raju, “R-S adapted arbitrary Lagrangian–Eulerian finite element method for metal-forming problems with strain localization,” International Journal for Numerical Methods in Engineering, vol. 39, no. 19, pp. 3247–3272, 1996.
View at: Publisher Site | Google Scholar
R. Glowinski, T. W. Pan, T. I. Hesla, and D. D. Joseph, “A distributed Lagrange multiplier/fictitious domain method for particulate flows,” International Journal of Multiphase Flow, vol. 25, no. 5, pp. 755–794, 1999.
View at: Publisher Site | Google Scholar
J. Hao, T. W. Pan, R. Glowinski, and D. D. Joseph, “A fictitious domain/distributed Lagrange multiplier method for the particulate flow of Oldroyd-B fluids: a positive definiteness preserving approach,” Journal of Non-Newtonian Fluid Mechanics, vol. 156, no. 1-2, pp. 95–111, 2009.
View at: Publisher Site | Google Scholar
Z. Feng and E. E. Michaelides, “The immersed boundary-lattice Boltzmann method for solving fluid–particles interaction problems,” Journal of Computational Physics, vol. 195, no. 2, pp. 602–628, 2004.
View at: Publisher Site | Google Scholar
A. Main and G. Scovazzi, “The shifted boundary method for embedded domain computations. Part II: Linear advection–diffusion and incompressible Navier–Stokes equations,” Journal of Computational Physics, vol. 372, pp. 996–1026, 2018.
View at: Publisher Site | Google Scholar
M. R. Hashemi, R. Fatehi, and M. T. Manzari, “SPH simulation of interacting solid bodies suspended in a shear flow of an Oldroyd-B fluid,” Journal of Non-Newtonian Fluid Mechanics, vol. 166, no. 21-22, pp. 1239–1252, 2011.
View at: Publisher Site | Google Scholar
M. Wang, Y. Feng, and T. Qu, “On the implicit immersed boundary method in coupled discrete element and lattice Boltzmann method,” International Journal for Numerical and Analytical Methods in Geomechanics, vol. 44, no. 4, pp. 516–532, 2019.
View at: Publisher Site | Google Scholar
Z. Guo and T. S. Zhao, “Lattice Boltzmann model for incompressible flows through porous media,” Physical Review E, vol. 66, no. 3, article 036304, 2002.
View at: Publisher Site | Google Scholar
M. Uhlmann, “An immersed boundary method with direct forcing for the simulation of particulate flows,” Journal of Computational Physics, vol. 209, no. 2, pp. 448–476, 2005.
View at: Publisher Site | Google Scholar
E. A. Fadlun, R. Verzicco, P. Orlandi, and J. Mohd-Yusof, “Combined immersed-boundary finite-difference methods for three-dimensional complex flow Simulations,” Journal of Computational Physics, vol. 161, no. 1, pp. 35–60, 2000.
View at: Publisher Site | Google Scholar
A. Dupuis, P. Chatelain, and P. Koumoutsakos, “An immersed boundary–lattice-Boltzmann method for the simulation of the flow past an impulsively started cylinder,” Journal of Computational Physics, vol. 227, no. 9, pp. 4486–4498, 2008.
View at: Publisher Site | Google Scholar
S. Liu, S. Tang, M. Lv, Y. Zhao, and Y. Li, “Simulating particle sedimentation in a flowing fluid using an immersed boundary–lattice Boltzmann method,” International Journal of Computational Fluid Dynamics, vol. 34, no. 1, pp. 39–49, 2020.
View at: Publisher Site | Google Scholar
R. Glowinski, T. W. Pan, T. I. Hesla, D. D. Joseph, and J. Périaux, “A fictitious domain approach to the direct numerical simulation of incompressible viscous flow past moving rigid bodies: application to particulate Flow,” Journal of Computational Physics, vol. 169, no. 2, pp. 363–426, 2001.
View at: Publisher Site | Google Scholar
C. S. Peskin, “The Immersed Boundary Method,” Acta Numerica, vol. 11, pp. 479–517, 2002.
View at: Publisher Site | Google Scholar
S. Majumdar, H. K. Genant, S. Grampp et al., “Correlation of Trabecular Bone Structure with Age, Bone Mineral Density, and Osteoporotic Status: In Vivo Studies in the Distal Radius Using High Resolution Magnetic Resonance Imaging,” Journal of Bone and Mineral Research, vol. 12, no. 1, pp. 111–118, 1997.
View at: Publisher Site | Google Scholar
J. Lipiec and R. Hatano, “Quantification of compaction effects on soil physical properties and crop growth,” Geoderma, vol. 116, no. 1-2, pp. 107–136, 2003.
View at: Publisher Site | Google Scholar
T. Hirata, “Fractal dimension of fault systems in Japan: Fractal structure in rock fracture geometry at various scales,” Pure and Applied Geophysics PAGEOPH, vol. 131, no. 1-2, pp. 157–170, 1989.
View at: Publisher Site | Google Scholar
D. Wan and S. Turek, “Direct numerical simulation of particulate flow via multigrid FEM techniques and the fictitious boundary method,” International Journal for Numerical Methods in Fluids, vol. 51, no. 5, pp. 531–566, 2006.
View at: Publisher Site | Google Scholar
A. Jafari and T. Babadagli, “Equivalent fracture network permeability of multilayer-complex naturally fractured Reservoirs,” Transport in Porous Media, vol. 91, no. 1, pp. 339–362, 2012.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2020 Shenggui Liu et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

602

Downloads

1029

Citations