Multi-Sized Granular Suspension Transport Modeling for the Control of Lost Circulation and Formation Damage in Fractured Oil and Gas Reservoirs

Abstract

Transport and retention of multi-sized suspended granules are common phenomena in fracture media of oil, gas and geothermal reservoirs. It can lead to severe permeability damage and productivity decline, which has a significant impact on the efficient development of underground resources. However, the granule transport and retention behaviors remain not well understood and quantified. The novel stochastic model is proposed for the multi-sized suspended granule transport in naturally fractured reservoirs accounting for granule retention and fracture clogging kinetics. A percolation fracture network is proposed considering fracture connectivity evolution during suspended granule transport. Granule retention and fracture clogging dynamics equations are proposed to account for incomplete fracture clogging by retained granules. The microscale stochastic model is allowed for upscaling to predict the multi-sized granule transport behavior in naturally fractured reservoirs. The model solution exhibits preferential plugging of fractures with sizes equal to or below the granule size. Multi-sized suspended granule shows great advantages over mono-sized suspended granule in the control of permeability damage induced by granule retention and fracture clogging. The retained granule concentration and permeability damage rate decrease with fracture network connectivity improvement. The experimental investigation on size-exclusion suspended granule flow has been performed. The model-based prediction of the retained granule concentration and permeability variation history shows good agreement with the experimental data, which verifies the developed model.

1. Introduction

Rapid population growth and economic expansion result in a greater requirement for energy consumption. Oil, gas and geothermal energy are the two main underground resources which have a significant influence on the supply of the global energy market [1,2,3]. Transport and retention of suspended granules in porous or fracture media are commonly encountered in oil, gas and geothermal reservoirs [4,5,6,7]. Size exclusion is the main granule retention mechanism during oil and gas production, reservoir well drilling, water storage in geothermal reservoirs, water disposal produced from oilfields, multi-phase heat transfer behaviors, hydraulic fracturing, injectivity impairment in oil, gas reservoirs, granules invasion and fines migration in reservoirs, etc. [8,9,10,11,12,13,14,15]. In the drilling of oil, gas and geothermal reservoirs, granules retention frequently happens during the invasion of granular suspension (e.g., drilling fluid and completion fluid) (Figure 1). It can lead to significant permeability damage and productivity decline, which seriously hinders the efficient development of oil, gas and geothermal resources [16,17,18].

A comprehensive understanding and quantification of suspended granule flow and retention behaviors are of great importance for the forecasting and control of permeability damage and subsequent productivity decline. Large amounts of studies have been carried out on the flow and capture behaviors of particulate suspension in the oil, gas and geothermal reservoirs with porous media properties. For porous media, granule retention by size exclusion is dependent on the sizes of the granule and pore system. The quantification of the transport and capture behaviors of granules in porous media requires the characterization of the size distribution of granules and pores, which is described by stochastic models [19,20]. There are a series of stochastic approaches for the description of the suspended granule flow in porous media, including population balance models, random-walk models, mean-field models, CFD models and Boltzmann models [21,22,23,24]. Pore size distributions, filtration coefficient and suspended and captured granules concentrations have been accounted for in these models [25,26]. Exact upscaling is allowed for the mono-sized suspended granule flow in the porous media [26].

Fracture media can also be commonly found in oil, gas and geothermal reservoirs (Figure 2). It contributes much more permeability and conductivity than normal porous media [27]. Nowadays, naturally fractured oil, gas and geothermal reservoirs in deeper, tighter and more complex conditions have turned into one hotspot for underground energy exploration and exploitation [28]. Naturally fractured reservoirs have an extensive distribution in worldwide, including the Middle East, North America and Central Asia. Sichuan and Tarim basins in China have large areas of naturally fractured reservoir distribution. The widely distributed natural fractures and low-permeability matrix are the most prominent characteristics of naturally fractured tight reservoirs. Developed natural fractures are of great importance for high-efficiency reservoir exploitation. However, they can also result in the invasion of drilling fluid which is one kind of particulate suspension (Figure 1). Severe permeability damage and subsequent productivity decline have been observed in this kind of reservoir because of granule retention and subsequent fracture clogging [29]. Moreover, Granules in the suspension are essentially distributed in size. However, most of the research focuses only on the mono-sized suspended granule flow and retention behaviors in porous media. Stochastic modeling for multi-sized suspended granule flow and retention in naturally fractured reservoirs remains not available in the literature.

In the process of suspended granule transport, granule retention can lead to great variations in the fracture media properties, including fracture density distribution, porosity and permeability. On the other hand, the changes in the properties of fracture media can influence the transport and retention behaviors of suspended granules [30]. Therefore, the evolution of fracture network connectivity with granule retention is of great importance for suspended granule transport modeling in fracture media. Percolation theory belongs to the probability theory used to predict the properties of random media, which usually has an association with network modeling [31]. It has a great advantage in the prediction of medium properties variation caused by granule retention [32]. However, few studies have given an introduction of the percolation model to the modeling of suspension flow and granule retention behaviors in fracture media. To our best knowledge, the behaviors of suspended granule flow and retention in fracture media remain not well understood and quantified.

The current work develops novel stochastic models for multi-sized suspended granule flow and retention in naturally fractured reservoirs. Percolation theory is introduced to develop the fracture network accounting for the connectivity evolution with granule retention and fracture clogging. Granule retention and fracture clogging dynamics equations are proposed considering incomplete fracture clogging by retained granules. The microscale model is allowable to be upscaled and a numerical solution is obtained for the macroscale model system for the prediction and control of permeability damage related to the multi-sized granule transport in fracture network. Laboratory experiments on granule transport in fracture media are performed to validate the proposed model.

The structure of the current paper is as follows. Section 2 derives the governing equations for suspended granule flow and retention in the fractured reservoir, including the models of percolation fracture network, granule retention kinetics and fracture clogging kinetics. The model system is then upscaled and numerically solved for multi-sized granule flow and retention in fracture media. Section 3 presents the sensitivity analysis for the suspended and retained granule concentration profiles, and the fracture network properties evolution (e.g., permeability and fracture density distribution) during suspended granule transport. Section 4 presents the laboratory materials, setup, methodology and model validation by the experimental results. Conclusions of Section 5 finalize the current work.

2. Governing Equations for Suspended Granule Transport

2.1. Percolation Model of Naturally Fractued Reservoir

According to the percolation theory, the geometric model for suspended granule transport and capture in fractured tight reservoirs is developed. Fractures developed in tight reservoirs are interconnected and then form a fracture network (Figure 2). Since the matrix permeability is relatively low (usually less than 0.1 mD), the fluid and granule phases flow through the matrix is neglectable compared with their transport in the fractures in actual operation. The fracture network is adopted to simulate the state of the fractured medium in the model. The main method of operation is to assign the fracture interspace qualities to percolation bonds to develop the percolation fracture network. These bonds whose intersecting surfaces are assumed to be rectangular are mainly formed by mineralization and shearing [4,32,33]. The fracture medium is composed of groups of fractures which are connected to the nodes. The influence of connection nodes on suspended particle flow and retention behaviors is very small and can be neglected.

During the drilling of oil, gas and geothermal reservoirs, the interaction of granules and matrix is repulsive and, therefore, size-exclusion is assumed as the only granule retention mechanism (Figure 3). Higher permeability damage is more likely to be caused by granule size exclusion than other granule retention mechanisms, such as granule segregation or attachment. The main reason is that the fracture aperture and length are greatly decreased after size-exclusion granule retention. Size-exclusion capture of granules in fractures can result in a significant decline in the permeability and productivity of fractured tight reservoirs [4]. During colloidal-suspension flow in fractured tight reservoirs, the transport and capture of granules in fractures occur at the same time. Size-exclusion granule retention occurs in smaller fractures whose aperture is less than granule size ds, while the larger fracture with an aperture greater than the granule diameter allows it to be passed through (Figure 3).

For a percolation fracture network, when a bond is occupied with probability p, the probability that the bond is empty is 1-p. The topological structure of the random networks is not constant, and there is a transition from macroscopic disconnection structure to connection structure. The critical value is referred to as percolation threshold pc. The percolation theory can differentiate between permitted and accessible bonds. Permitted bonds are those that are allowed to be entered by one phase. However, it ignores the probability that the ambient bonds are for example too small so that there is no chance for the invading phase to try. Accessible bonds represent the bonds that are permitted and at the same time, the surrounding bonds will not prevent the invading phase from entering.

Regarding the suspended granule flow in fracture media, a fracture is permitted as long as its aperture is larger than the granule diameter. The fraction p of fractures which are permitted for granules with size ds is expressed by

$$p={\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{d_s}{\overset{\infty }{\int }}f}}(a_f,l_f)d{a}_{f}d{l}_f.$$

(1)

The number of fractures per unit intersecting the surface area of the fracture media is determined by the fracture density distribution function.

$$\tau (x,t)={\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{0}{\overset{\infty }{\int }}T}}(w_f,l_f,x,t)d{w}_{f}d{l}_f,$$

(2)

where T is the fracture density distribution. The probabilistic distribution function for fracture with aperture a_f and length l_f is defined by

$$f(a_f,l_f,x,t)=\frac{T(a_f,l_f,x,t)}{\tau (x,t)}.$$

(3)

Normal and lognormal distributions are the common distributions of fracture aperture and length in the fractured reservoirs. Normal distribution of fracture size is adopted for analysis in this work.

Fracture porosity is equal to the overall intersecting surface area of all fractures per unit intersecting surface area of fracture media and expressed as

$$\mathrm{\Phi }[T]={\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{0}{\overset{\infty }{\int }}A}}(a_f,l_f)T(a_f,l_f,x,t)d{a}_{f}d{l}_f,$$

(4)

where the fracture intersecting surface area A is equal to the fracture aperture a_f multiplied by fracture length l_f. Since fractures can only be entered by these granules with diameter d_s smaller than fracture aperture a_f, the allowed or permitted fraction of overall fracture porosity can be determined according to the integration of the intersecting surface area of fractures with aperture larger than granule diameter

$${\mathrm{\Phi }}_a[t,d_s]={\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{d_s}{\overset{\infty }{\int }}A}}(a_f,l_f)T(a_f,l_f,x,t)d{a}_{f}d{l}_f.$$

(5)

As the above discussion, the allowed fraction of fracture can be directly determined based on the relative sizes of fracture aperture and granule diameter. However, the accessible fraction of fracture also depends on the relative sizes between granule diameter and fracture aperture of the surrounded fractures. Even if one fracture is allowed, it is still not accessible if its surrounded fractures are not allowed. Therefore, the allowed fraction of fracture porosity is not equal to the accessible fraction in most cases. Percolation theory can be used to determine the relationship between the allowed fraction and accessible fraction. The accessible fracture porosity can be expressed by the allowed fracture porosity, fraction of allowed fractures and percolation threshold p_c [26].

$${\varphi }_{ac}=\left\{\begin{array}{ll}\hspace{1em}\hspace{1em}0& p<p_c\\ \frac{p_{sc}}{p}\left(\frac{p-p_c}{p_{sc}-p_c}\right){\varphi }_a& p\le p_c\le p_{sc}\\ \hspace{1em}\hspace{1em}{\varphi }_a& p>p_{sc}\end{array}\right.$$

(6)

where p_sc is equal to 1.65p_c and β is a common constant with value of 0.4.

Besides fracture media porosity, another important parameter for the fracture media is permeability which characterizes the fracture network conductively. First, the fracture conductivity distribution function F(k_p) is defined. The fracture conductivity k_p uniquely depends on the fracture aperture and length. For fluid phase without granules, the relationship between the fracture conductivity distribution and the fraction of allowed fractures can be expressed by

$$F(k_p)d{k}_p=f(a_f,l_f)d{a}_{f}d{l}_f,$$

(7)

where k_p is the single fracture conductivity, which is equal to the permeability of a single fracture with a rectangular intersecting surface.

The effective-medium-approximation method is used to calculate the effective conductivity of the fracture network [19]. For a single fluid phase without granule, all the fracture is conductive, and the fraction of allowed fractures is 1. Application of the effective medium method results in

$${\displaystyle \underset{0}{\overset{\infty }{\int }}F}\left(k_p\right)\frac{k_p-k}{k_p+\left({\gamma }^{-1}-1\right)}d{k}_p=0,$$

(8)

where γ = 2/Z and k is the absolute permeability of the fracture media. Substitution of Equation (7) into Equation (8) leads to

$${\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{0}{\overset{\infty }{\int }}f}}\left(a_f,l_f\right)\frac{k_p-k}{k_p+\left({\gamma }^{-1}-1\right)k}d{a}_{f}d{a}_f=0.$$

(9)

For granular suspension, only part of the fractures is conducive to the granule phase. The fractures are non-conductive and not allowed for granules with diameters larger than the fracture aperture. The fractures are allowed for the granule whose diameter is smaller than the fracture aperture. However, not all the allowed fractures are accessible. The effective medium method takes care of the accessibility issue of fractures automatically [4,30]. For the granule phase, the fracture conductivity distribution function is calculated by the allowed fraction, fracture size distribution function and the relative size between granule and fracture

$$F(k_p)=(1-p)\delta (k_p)+f(a_f,l_f)U(a_f-d_s)\frac{d{a}_{f}d{l}_f}{d{k}_p},$$

(10)

where δ(x) is Dirac delta function. For x equals 0 and nonzero value, δ(x) equals 1 and 0, respectively. U(x) equals 1 for x > 0 and U(x) equals 0 for x < 0. Substitution of Equation (10) into Equation (8) leads to

$$\frac{p-1}{{\gamma }^{-1}-1}+{\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{d_s}{\overset{\infty }{\int }}f(a_f,l_f)\frac{k_p-k_{ac}}{k_p+({\gamma }^{-1}-1)k_{ac}}d{a}_f}d{l}_f}=0,$$

(11)

where k_ac is referred to as accessible permeability of the fracture media.

The overall transport rate is obtained by adding flow rate through the accessible fractures to that through the inaccessible fractures. For any granule with specific size, the accessible transport fraction fac is expressed as

$$f_{ac}=\frac{k_{ac}}{k}.$$

(12)

2.2. Population Balance Model and Granule Retention Kinetics for Size-Exclusion Suspension Transport

Overall, the concentration of suspended and retained granules is determined by the integrations of the relevant concentration distribution functions.

$$c\left(x,t\right)={\displaystyle \underset{0}{\overset{\infty }{\int }}C}\left(d_s,x,t\right)d{d}_s,$$

(13)

$$s\left(x,t\right)={\displaystyle \underset{0}{\overset{x}{\int }}S}\left(d_s,x,t\right)d{d}_s,$$

(14)

$$S(d_s,x,t)={\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{0}{\overset{\infty }{\int }}\underset{¯}{S}(w_f,l_f,d_s,x,t)d{w}_{f}d{l}_f}},$$

(15)

where c and σ are the suspended and retained granule concentration within the fracture media, respectively. C and S is the suspended and retained granule concentration distribution function, respectively. Retained granule concentration can also be distributed according to fracture size. Therefore, S is introduced as the distribution of retained granule concentration over fracture size and granule diameter.

The system of governing equations for suspended granule flow and retention in fracture media includes the mass balance equation for suspended and retained granules, and the kinetic equations for granule capture and fracture clogging rates. The system of three governing equations determines the three unknown distribution functions: C, S and T.

The mass balance equation which neglects the diffusion considers the granules in the suspension and retained in the fracture and is described in the following form:

$$\frac{\partial {\varphi }_{\alpha }\left[p_c,T,d_s\right]C\left(d_s,x,t\right)}{\partial t}+U\frac{\partial f_{ac}\left[Z,T,d_s\right]C\left(d_s,x,t\right)}{\partial x}=-\frac{\partial S\left(d_s,x,t\right)}{\partial t}.$$

(16)

The granule retention kinetics is determined based on the hypothesis that granule retention in the fractures with specific sizes is in direct proportion to the granule flux through these fractures. The retention ratio of d_s granules which are seized by fractures with an aperture less than the granule diameter is equivalent to the total granule flux into the fractures with a smaller aperture.

$$\frac{\partial \underset{¯}{S}\left(a_f,l_f,d_s,x,t\right)}{\partial t}=q\left(a_f,l_f\right)T^{\nu }\left(a_f,l_f,x,t\right)C\left(d_s,x,t\right),$$

(17)

where T^v is referred to as the volumetric fracture density distribution. Besides the granule size, it is also a function of the size of fracture which the granules are retained in. The volumetric fracture density distribution is associated with the surface fracture density distribution by equality $T^{\nu }\left(a_f,l_f,x,t\right)=\frac{1}{l}T\left(a_f,l_f,x,t\right)$.

Integrating Equation (17) with respect to fracture size smaller than the granule diameter leads to the expression for the granule retention kinetics

$$\frac{\partial S\left(d_s,x,t\right)}{\partial t}=\frac{UC\left(d_s,x,t\right)}{lk\left[T\right]}{\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{0}{\overset{d_s}{\int }}{k}_1}}\left(a_f,l_f\right)T\left(a_f,l_f,x,t\right)d{a}_{f}d{l}_f.$$

(18)

In view of the fractional flow function, substitution of Equation (12) into Equation (18), the granule retention kinetics can be described by

$$\frac{\partial S(d_s,x,t)}{\partial t}=\frac{UC(d_s,x,t)}{l}(1-f_{ac}\left[Z,T,d_s\right]).$$

(19)

2.3. Fracture Clogging Kinetics

The dimensionless parameter p_s is proposed to describe the possibility the d_s granule is captured in fractures with a_f aperture and l_f length and the fracture length l_f is subsequently shortened to l_f because of the fracture clogging by the retained d_s granule. By integrating the reduced fracture length, the dimensionless parameter p_s is obtained, which is used to describe the possibility that a granule with diameter d_s is seized within fracture with aperture a_f and l_f length. The parameter has been used to describe incomplete pore plugging by retained granules and subsequent permeability damage [19]. These two introduced dimensionless coefficients are described by

$$\left\{\begin{array}{l}{p}_s(d_s,a_f,l_f\to {l^{\prime }}_f)=0:l_f\le {l^{\prime }}_f\\ {\displaystyle \underset{0}{\overset{\infty }{\int }}{p}_s(d_s,a_f,l_f\to {l^{\prime }}_f)d{l^{\prime }}_f=p_s(d_s,a_f,l_f)}\end{array}\right..$$

(20)

The difference between the “increase” and “decrease” terms can be adopted to characterize the density distribution change of the fractures with aperture a_f and length l_f. The “increase” term means that longer fractures with the same with are plugged by retained granules and their length are reduced to l_f. The “decrease” term means that the fractures with aperture a_f and length l_f change to other fractures with the same aperture but smaller length after they are plugged by retained granules.

$$\frac{\partial T\left({\alpha }_f,l_f,x,t\right)}{\partial t}=I\left({\alpha }_f,l_f,x,t\right)-D\left({\alpha }_f,l_f,x,t\right).$$

(21)

According to the retention possibility defined above by accounting for the dimensionless parameter Ps, the “increase” and “decrease” terms in Equation (21) can be given by

$$I\left({\alpha }_f,l_f,x,t\right)={\displaystyle \underset{l_f}{\overset{\infty }{\int }}d{l^{\prime }}_f{\displaystyle \underset{a_f}{\overset{\infty }{\int }}d{d}_s\left\{p_s(d_s,a_f,{l^{\prime }}_f\to l_f)q(a_f,{l^{\prime }}_f)T(a_f,{l^{\prime }}_f,x,t)C(d_s,x,t)\right\}}},$$

(22)

$$D\left({\alpha }_f,l_f,x,t\right)={\displaystyle \underset{a_f}{\overset{\infty }{\int }}d{d}_s{\displaystyle \underset{0}{\overset{\infty }{\int }}d{l^{\prime }}_f\left\{p_s(d_s,a_f,l_f\to {l^{\prime }}_f)q(a_f,l_f)T(a_f,l_f,x,t)C(d_s,x,t)\right\}}}.$$

(23)

The possibility distribution functions p_s(d_s,a_f,l_f→l_f^’) and p_s(d_s,a_f,l_f^’→l_f) are described by

$$p_s(d_s,a_f,{l^{\prime }}_f\to l_f)=\left\{\begin{array}{ll}0& a_f>d_s\\ \frac{1}{l_{\max }-l_f-L({\alpha }_f,l_f,x,t)}& a_f<d_s\end{array}\right.,$$

(24)

$$p_s(d_s,a_f,l_f\to {l^{\prime }}_f)=\left\{\begin{array}{ll}0& a_f>d_s\\ \frac{1}{l_f-L({\alpha }_f,l_f,x,t)}& a_f<d_s\end{array}\right.,$$

(25)

where L is determined based on the overall length of the retained granules in the fractures with aperture a_f and length l_f:

$$L\left(a_f,l_f,x,t\right)={\displaystyle \underset{0}{\overset{\infty }{\int }}\frac{\underset{¯}{S}\left(a_f,l_f,d_s,x,t\right)d_s}{T^{\nu }\left(a_f,l_f,x,t\right)}dds}.$$

(26)

After substituting Equations (22) as well as (23) to Equation (21), the fracture clogging rate is given by

$$\frac{\partial T(a_f,l_f,x,t)}{\partial t}={\displaystyle \underset{a_f}{\overset{\infty }{\int }}C(d_s,x,t)d{d}_s\left\{\begin{array}{l}\frac{1}{l_{\max }-l_f-L({\alpha }_f,l_f,x,t)}{\displaystyle \underset{l_f}{\overset{\infty }{\int }}q(a_f,{l^{\prime }}_f)T(a_f,{l^{\prime }}_f,x,t)d{l^{\prime }}_f}\\ -q(a_f,l_f)T(a_f,l_f,x,t)\end{array}\right\}}.$$

(27)

Suspended granule concentration C, retained granule concentration S and fracture density distribution T are the three to be solved parameters within the three governing model systems (16), (19) and (27). The behaviors of suspended granule transport and size-exclusion granule retention accounting for incomplete fracture clogging are determined by the stochastic model system. The permeability damage can be effectively predicted by solving the model system. The initial and boundary conditions of the stochastic model system are given by

$$t=0:C\left(d_s,x,0\right)=0,T\left({\alpha }_f,l_f,x,0\right)=T_0({\alpha }_f,l_f,x),$$

(28)

$$x=0:C(d_s,0,t)=C^0(d_s,t),$$

(29)

where T₀ and C⁰ is the initial fracture density distribution and initial granule concentration in the injected suspended granule, respectively.

2.4. Upscaling and Solution for Multi-Sized Suspension Transport in Fracture Media

Dirac delta function δ is used to express the granule concentration distribution under the condition of multi-size suspended granule transport within fracture media.

$$C\left(d_s,x,t\right)=c_1\left(x,t\right)\delta \left(d_s-d_{s1}\right)+c_2\left(x,t\right)\delta \left(d_s-d_{s2}\right)+\cdots \cdots +c_n\left(x,t\right)\delta \left(d_s-d_{sn}\right),$$

(30)

$$c\left(x,t\right)=c_1\left(x,t\right)+c_2\left(x,t\right)+\cdot \cdot \cdot \cdot +c_n\left(x,t\right),$$

(31)

where d_sn is the diameter of the Nth granule. Substitution of c(d_sn,x,t) into system (16), (19) and (27) and integration with regard to d_s leads to the macroscale governing equations as

$$\frac{\partial {\varphi }_{ac}\left[P_c,T,d_s\right]c\left(x,t\right)}{\partial t}+U\frac{\partial f_{\alpha c}\left[Z,T,d_s\right]c\left(x,t\right)}{\partial x}=-\frac{\partial s\left(x,t\right)}{\partial t},$$

(32)

$$\frac{\partial s\left(x,t\right)}{\partial t}=\frac{1}{l}Uc\left(x,t\right)\left(1-f_{ac}\left[p_c,T,d_s\right]\right),$$

(33)

$$\frac{\partial T(a_f,l_f,x,t)}{\partial t}=\frac{Uc(x,t)}{k}\left\{\begin{array}{l}\frac{1}{l_{\max }-l_f-L(a_f,l_f,x,t)}{\displaystyle \underset{l_f}{\overset{\infty }{\int }}{k}_1(a_f,{l^{\prime }}_f)T(a_f,l_f,x,t)d{l^{\prime }}_f}\\ -k_1(a_f,l_f)T(a_f,l_f,x,t)\end{array}\right\}.$$

(34)

The numerical solution of the developed stochastic model is given by finite difference method as follows:

$$\frac{{\varphi }_{acm}^{n+1}{C}_m^{n+1}-{\varphi }_{acm}^n{C}_m^n}{\Delta t}+U\frac{f_{acm+1}^n{c}_{m+1}^n-f_{acm-1}^n{c}_{m-1}^n}{2\Delta x}=-\frac{{\sigma }_m^{n+1}-{\sigma }_m^n}{\Delta t},$$

(35)

$$\frac{s_m^{n+1}-s_m^n}{\Delta t}=\frac{1}{l}U{c}_m^{n+1}(1-f_{acm}^{n+1}),$$

(36)

$$\frac{T{(a_f,l_f)}_m^{n+1}-T{(a_f-l_f)}_m^n}{\Delta t}=\frac{U{c}_m^n}{k_m^n}\left\{\begin{array}{l}\frac{1}{l_{\max }-l_f-L{(a_f,l_f,x,t)}_m^n}{\displaystyle \underset{l_f}{\overset{\infty }{\int }}{k}_p(a_f,{l^{\prime }}_f)T{(a_f-{l^{\prime }}_f)}_m^{n}d{l^{\prime }}_f}\\ -k_1(a_f,l_f)T{(a_f-l_f)}_m^n\end{array}\right\}.$$

(37)

3. Quantitative Analysis of the Solution

In this section, a quantitative analysis of the modeling solution is conducted. The fracture distribution data is adopted from the fractured reservoir in Sichuan Basin, China. It shows that the initial aperture and length of fractures exhibit a normal distribution. The mean value and standard deviation for fracture aperture are 0.525 mm and 0.1 mm, respectively (Figure 4a). Because of its developed natural fractures and ultra-low matrix permeability, the fractured tight reservoir is regarded as fracture media. Multi-sized suspended granule is injected into the medium with granule diameter of 0.3 mm, 0.4 mm, 0.45 mm, 0.5 mm and 0.55 mm. According to the model proposed in the second section, the concentration distribution of capture granules and the evolution process of fracture medium properties are calculated and analyzed in detail.

Figure 4 shows the variation of fracture density distribution because of granule retention and fracture clogging during suspended granule flow in the fracture network. Fractures of aperture smaller than the granule size (a_f < d_s) are plugged during granules capture. The concentration distributions of fractures with smaller aperture decline gradually with granules injection due to fracture clogging by retained granules (Figure 4b,c). For the fractures with the same aperture, the concentration distribution for shorter fractures increases gradually because the longer fractures are plugged into shorter ones during granules capture. With the continuous granule capture within fractures, the fracture length becomes shorter and shorter and the finally vanish. For the fractures with an aperture greater than the granule diameter, no granule retention occurs within them, and their density distribution keeps constant (Figure 4d). The reason is that size exclusion is the fundamental granule retention mechanism considered in this work.

Figure 5 presents the density distribution evolution of the fractures with the same length during multi-sized granules transport in the fracture network. The modeling results show that fracture density distribution declines with time for fractures with aperture a_f less than the granule diameter d_s because of granule size exclusion; while T remains constant for fractures of aperture greater than the maximum granule diameter (a_f > d_smax). The results indicate that greater fractures (a_f > d_smax) remain accessible to the granules all the time.

The fractures whose apertures are equal to the granule size decline faster than those with size smaller than the granule size. The reason is that fracture with greater aperture is characterized by larger suspension flux, which results in higher granule retention rate. Higher retention rate leads to faster decline of fracture density. For different granule diameters in Figure 5, greater granule diameter leads to larger fracture density decrease, which further negatively affects the fracture network conductivity and permeability. Therefore, the fracture clogging effect can be controlled by the rational optimization of granule diameter in the colloidal suspension.

The evolution of retained granule concentration and fractured network permeability during the retention of multi-sized granules of 0.3 mm, 0.4 mm and 0.5 mm diameters are shown in Figure 6 and Figure 7. Figure 6 shows the evolution of retained granule concentration profiles for multi-sized granules. Modeling results in Figure 6a,b indicate that the retained granule concentration raises with time and declines with distance in non-linear form. Figure 7 shows the variation of fractured network permeability with time and distance for multi-sized granules. Because of granule retention and resultant fracture clogging, the conductivity and dimensionless permeability of the fracture network decreases in a non-linear form (Figure 7a). Moreover, the permeability declines more quickly than the permeability far away from the inlet (Figure 7b). The reason is that the retained granule concentration is larger in the near-inlet region based on the retention concentration results in Figure 6.

Figure 8 compares the retained granule concentration between mono-sized and multi-sized granules during suspended granule transport. For the mono-sized granule, larger granule diameter results in larger retained granule concentrations in the near inlet area (Figure 8a). It is because the path tortuosity is much higher for larger granules to pass through the fractured network, which leads to higher capture probability and capture rate for larger granules (Figure 9). However, lower retained granule concentration in the near-inlet region often corresponds to greater granule invasion depth due to the lower fracture clogging by retained granules (Figure 8b). Therefore, the retained granule concentration in the near and far inlet region cannot be adjusted simultaneously for the mono-sized granules. For the multi-sized granules, the retained granule concentration in the inlet region and the granule invasion depth can be controlled effectively at the same time according to the rational optimization of granule size combination and their concentrations (Figure 8b).

Figure 8 compares k/k₀ between mono-sized and multi-sized granules during suspended granule transport. For the mono-sized granule, higher granule size results in lower dimensionless absolute permeability in the near inlet area due to the higher retained granule concentration and fracture clogging effect (Figure 8a and a). However, the permeability decline range is much wider for smaller granule than that of the larger granules (Figure 10b). The permeability decline rate and range cannot be controlled at the same time due to the retained granule concentration profile for the mono-sized granules. However, for the multi-sized granules, the permeability decline rate and range can be controlled effectively and simultaneously by the rational optimization of granule size combination and their concentrations (Figure 10b).

Three key granule sizes are often used to represent the granule size distribution of suspended granules, including D90, D50 and D10. D90 represents the value in which 90% of the granule size in the suspended granule is less than and the same as D50 and D10. Three-sized granules are used to analyze the effect of D90, D50 and D10 on the profile of dimensionless absolute permeability during suspended granule transport (Figure 11). D90 value mainly affects the permeability decline degree in the near inlet area as well as the permeability decline range (Figure 11a). A higher D90 value leads to a larger permeability decline in the inlet area and a smaller permeability decline range. D50 value mainly affects the average permeability decline degree along the radial distance but has little impact on the permeability decline range (Figure 11b). Higher D50 value results in a larger average permeability decline degree. D10 value slightly affects the permeability decline profile (Figure 11c). However, the D10 value can be used to adjust the D90 and D50 values of the suspended granule.

The parameters which characterize the fracture network connectivity mainly include coordination number Z and percolation threshold p_c. Their impacts on the flow and retention behaviors of multi-sized suspended granule are presented in Figure 12, Figure 13, Figure 14 and Figure 15. Larger Z and smaller p_c improve the connectivity of the fracture network. Figure 12 presents the variation of retained granule concentration for multi-sized granules capture in fracture network with different connectivity. Retained granule concentration declines with the network connectivity increase. The reason is that lower fracture network connectivity yields higher granule path tortuosity, which leads to higher probability for the granules to come into the smaller aperture fractures (a_f < d_s) and be retained (Figure 13). Therefore, a higher fracture network connectivity also results in lower decline of fracture permeability as well as fracture density because of lesser retained granule concentration during the transport of multi-sized suspended granule (Figure 14 and Figure 15).

4. Laboratory Experiments

In this part, laboratory experiments of suspension transport and granule retention are carried out. Retained granule concentration, pressure drop and permeability changes during suspended granule flow in fracture media are measured and compared with above modeling results for model validation.

4.1. Experimental Design

Fractured core sample is adopted as the fracture media for suspended granule flow and retention experiments in this paper. The diameter and length for the fractured core sample is 25 mm and 60 mm, respectively. The core sample is obtained by coring operation in the fractured tight reservoirs of Sichuan Basin, China. Compared with fracture permeability, the matrix permeability is lower than 0.01 mD, which can be neglected. Therefore, the flow rate of fluid through pores is negligible compared with the fracture. The fracture aperture distribution of the core sample is determined by three-dimensional laser scanner. The fractured surface morphology images of the fractured core sample obtained from the three-dimensional laser scanner are shown in Figure 16. Fracture aperture distribution density can be determined according to the fracture surface contour images and presented in Figure 17. The contribution rate of fracture permeability and its cumulative distribution in Figure 17 are obtained based on the distribution density of fracture aperture. As shown in Figure 17, when the fracture aperture is 0.7 to 0.8 mm, the maximum contribution rate appears.

Since the repulsion criteria of granule/grain and granule/granule are used in these experiments, no granules can be attached and retained on the fracture surface. Under this condition, we can conclude that size exclusion is the only mechanism of granule retention in the experiment. This conclusion is consistent with the model hypothesis. With regard to a specific granule material and fracture medium, the repulsion condition is realized according to the adjustment of fluid properties, including pH value as well as salinity value [11]. Zero salinity of water, alkalinity of suspension with pH ≈ 11 and granules covered by carboxyl functional groups can result in size exclusion as the only mechanism of granule retention. These criteria are adopted in the granule suspension flow and granule retention experiments in the current paper.

As is shown in Figure 18, we can have an observation about the application of the experimental apparatus in the suspended granule flow and retention experiments. ISCO pump is used to inject granules in suspension through the fractured core sample at a constant rate. Another pump provides confining pressure to the core sample holder. A fraction collector is connected in parallel with a flowmeter. They are used for the constant collection of the outlet flux and measure the output flux. Two pressure sensors are adopted to monitor the pressure at both ends of the core sample, respectively. The outlet pressure is constant and equal to atmospheric pressure. All the pressure values are recorded in real time by a data acquisition system (DAS). All the setup is put in a thermotank in order to keep constant 25 °C. Before we carry out the suspension transport experiments, alcohol, soapy water and deionized water are used to clean all the flow elements of the experimental system so that the effect of residual contaminants can be avoided.

Calcium carbonates are adopted as suspended granules in this work. With the help of the prepared solution, the fractured core sample is saturated in a vacuum condition for 48 h. On this basis, suspensions with granule sizes of 0.3 mm, 0.5 mm and 0.7 mm are, respectively injected into the fractured sample at the same volume flow rate. Mono-sized and multi-sized granules are injected for comparison with the same overall concentration. With regard to each case, outlet flow rate is regularly measured by means of the flow meter. The pressure sensor is used to record the inlet and outlet pressures in real time. After each test, retained granules in the fracture were collected for every 5 mm section to obtain the spatial distribution of retained granule concentration. Retained granules are sieved to calculate the granule number of each size. During the whole experiment, Darcy’s law is used to calculate the permeability changes caused by granule retention and resultant fracture clogging.

4.2. Experimental Results and Discussion

Figure 19 shows spatial concentration distribution of retained granule for mono-sized and multi-sized granule with the same injected granules concentration. For the mono-sized granules, the 0.7 mm granule corresponds to the highest retained granule concentration during the core inlet region and the lowest granule invasion depth. For the multi-sized granules, the spatial concentration distribution of retained granules fall in between the distribution of mono-sized granules. Laboratory experiments results in Figure 19 indicate the same trend of retained granule concentration as those calculated from the proposed stochastic models (Figure 8b).

Figure 20 presents the pressure drop during injection of granule suspensions with the same granule concentration but different granule sizes. The pressure differential is a reflection of fracture clogging which is caused by plugging capture by size exclusion. The pressure drop for injected granule size in the descending order is 0.7 mm, multi-sized granule combination (0.7 mm, 0.5 mm and 0.3 mm), 0.5 mm and 0.3 mm. Corresponding to the pressure drop, the k/k₀ evolutions during the transport of suspended granules with the same concentration but different granule sizes are shown in Figure 20. Dimensionless permeability k/k₀ decline indicates the fracture clogging by retained granules. The dimensionless permeability decline for injected granule size in the descending order is also 0.7 mm, multi-sized granule combination (0.7 mm, 0.5 mm and 0.3 mm), 0.5 mm and 0.3 mm. For mono-sized granules, laboratory results show that 0.7 mm granules can achieve the most efficient plugging of the fractures with maximum permeability contribution rate (0.7–0.8 mm). The dimensionless permeability for multi-sized granules falls in between that of mono-sized granules. The laboratory experiment results in Figure 20 indicate that the decreasing trend of k/k₀ is consistent with the results of the proposed model results in Figure 10a.

Another important finding from the model results is that fractures with apertures close to the injected granule diameter vanish more quickly. In addition, when the fracture aperture is equal to the injected granule diameter, the fracture disappears fastest. Fracture clogging with the largest contribution rate of permeability results in larger declines of k/k₀. The fracture with an aperture of 0.7~0.8 mm makes the largest contribution. From the experimental results in Figure 20, we can see that the dimensionless permeability for 0.7 mm granules decline the most compared with those of 0.5 mm and 0.3 m granules during suspended granule injection. These results indicate that fractures with the largest contribution rate are plugged faster by 0.7 mm granules than the 0.5 mm and 0.3 mm granules. Therefore, the laboratory experiment results indicate that the fractures with an aperture (0.7~0.8) mm close to the injected granule diameter (0.7 mm) disappear more quickly. The results are in good agreement with the model solution and verify the developed model qualitatively.

Because of the lack of fracture length distribution data, the model proposed in the current paper has only been qualitatively verified. In order to realize the quantitative verification of the proposed model, the data set of fracture size density distribution should be determined through seismic response. Moreover, for suspension transport with higher particle concentration, particle bridging can occur and become another particle capture mechanism. In this case, size exclusion is not the only particle capture mechanism. It will be investigated further in future studies. However, in this stage, the prediction of the profile and evolution of the retained granule concentration and dimensionless permeability by the proposed stochastic model is instructive enough for permeability damage control and production enhancement in fractured oil, gas and geothermal reservoirs.

5. Conclusions

Mathematical and laboratory modeling of multi-sized suspended granule transport in naturally fractured reservoir accounting for granule retention and fracture clogging kinetics allows drawing the following conclusions:

Percolation theory is introduced to develop the fracture network model accounting for the connectivity evolution with the granule retention and resultant fracture clogging. Granule retention and fracture clogging dynamics equations are proposed considering the behavior of incomplete fracture clogging by retained granules. The micro-scale stochastic model allows for upscaling and the numerical solution is obtained for the macro-scale model system for the prediction of the multi-sized granule transport behavior in the fracture network.

During suspended granule flow and retention in the percolation fracture network, the granule retention and fracture clogging behaviors can be predicted quantitatively. The model solution exhibits preferential plugging of fractures with an aperture equal to or below the injected granule diameter. The fractures with an aperture equal to the granule diameter vanish fastest. The retained granule concentration and permeability damage rate decrease with the improvement of the network connectivity.

Multi-sized suspended granules show great advantages over mono-sized suspended granules in the control of granule retention and fracture clogging effect. Based on the rational optimization of granule size combination and concentration, the profiles for granule invasion and permeability decline can be controlled effectively and simultaneously. The proposed model shows good agreement with laboratory data.