A new capillary pressure model for fractal porous media using percolation theory

Journal of Natural Gas Science and Engineering 41 (2017) 7e16

Contents lists available at ScienceDirect

Journal of Natural Gas Science and Engineering journal homepage: www.elsevier.com/locate/jngse

A new capillary pressure model for fractal porous media using percolation theory Jun Zheng a, Hongbo Liu b, *, Keke Wang a, Zhenjiang You c a

State Key Laboratory of Oil and Gas Reservoir Geology and Exploitation, Chengdu University of Technology, Chengdu 610059, Sichuan, China Sichuan Water Conservancy Vocational College, Chengdu 611830, Sichuan, China c Australian School of Petroleum, The University of Adelaide, Adelaide, SA, 5005, Australia b

a r t i c l e i n f o

a b s t r a c t

Article history: Received 12 July 2016 Received in revised form 17 January 2017 Accepted 26 February 2017 Available online 2 March 2017

The capillary pressure-saturation relation in traditional models is obtained by representing porous medium as a bundle of capillary tubes. The effect of pore connectivity is ignored, which leads to the evident deviation frequently observed between the measured data and model prediction near the breakthrough pressure region. In the present work, a new fractal-percolation-based capillary pressure model is developed. The proposed model accounts for the characteristics of pore size distribution in fractal media, and is based on the percentage of pore volume rather than that of pore number in conventional percolation theory. In the new model, a simple and yet sufficiently accurate relation between the allowed and accessible saturations is applied near the threshold saturation. It results in a highaccuracy prediction of capillary pressure, especially near the breakthrough pressure region. The new model is validated by experimental data. Parameter sensitivity analysis and comparison with the results from conventional fractal theory are presented. This new model is an extension to the conventional fractal theory, as it is consistent when the threshold saturation is equal to zero. © 2017 Elsevier B.V. All rights reserved.

Keywords: Capillary pressure Fractal porous media Percolation theory Pore connectivity Breakthrough pressure

1. Introduction The capillary pressure in a porous media is the result of the combined effects of the surface and interfacial tensions of the rock and fluids, the pore size and geometry, and the wetting characteristics of the system. The capillary pressure, as a function of saturation, is a fundamental property affecting immiscible multiphase flow in porous media, which makes it of great research significance. The research on the capillary pressure behavior includes both experimental work and theoretical model development. Laboratory methods to obtain the relation between capillary pressure and wetting phase saturation in porous media mainly includes mercury intrusion, semipermeable-membrane, self-absorption and centrifugal method. Among these methods, the mercury intrusion is widely used because this method is simple to conduct and rapid. The data can be used to determine the pore size distribution, to study the behavior of capillary pressure curves, and to infer characteristics of pore geometry. In addition, O'Meara et al.

* Corresponding author. E-mail address: [email protected] (H. Liu). http://dx.doi.org/10.1016/j.jngse.2017.02.033 1875-5100/© 2017 Elsevier B.V. All rights reserved.

(1988) showed that mercury injection capillary pressure data of water-oil systems (normalized using Leverett J-function) are in good agreement with the strongly water-wet capillary pressure curves obtained by other methods. Brown (1951) found that gas-oil capillary pressure data can be made to agree with mercury injection capillary pressure data by using an appropriate scaling factor. For theoretical model development, the capillary pressure model is established mainly through describing the pore structure characteristics of porous media by the simplified capillary bundle model. Researchers use various techniques to represent the pore structure of natural formation, including sandstone, carbonate, shale, coal, etc. (Wong et al., 1986; Friesen and Mikula, 1987; Krohn, 1988; Angulo and Gonzalez, 1992; Zhang et al., 2014; Naik et al., 2015; Xiong et al., 2015; Zhou and Kang, 2016). The experimental results show that, the pore space of porous media formed with randomly stacked natural particles consists with fractal characteristic (Pfeifer and Avnir, 1983; Katz and Thompson, 1985; Hansen and Skjeltorp, 1988; Bu et al., 2015; Lai and Wang, 2015; SakhaeePour and Li, 2016). Hunt (2004) suggested that fractal concepts might be well suited to model fluid movement in a porous medium, because of their simple descriptions of highly ramified spaces. Different fractal

8

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

capillary pressure models have been developed based on different hypothetical fractal models, including de Gennes (1985), Friesen and Mikula (1987), Tyler and Wheatcraft (1990), Rieu and Sposito (1991), Shen and Li (1994), Perrier et al. (1996), Deinert et al.
n and Gonza
lez-Posada (2005), Perfect (2005), Cihan (2005), Milla et al. (2007), Li (2010), Liu et al. (2016), and so on. Among these models, the most widely used ones are de Gennes (1985), Friesen and Mikula (1987), and Rieu and Sposito (1991) models. The last one is usually referred to as the RS model. All the three models agree with the measured data very well, except for the breakthrough pressure neighborhood. For drainage-type capillary pressure curve, the region near the breakthrough pressure corresponds to the large pores. Shen and Li (1994) applied fractal model to fit the capillary pressure curve obtained from mercury intrusion method. Their results indicated that traditional fractal model cannot match the trend of capillary pressure-saturation behavior in the breakthrough pressure neighborhood, as stated above.
n and Gonza
lez-Posada (2005) Based on this phenomenon, Milla proposed the piecewise fractal model by considering that the capillary pressure curve includes both structural and textural pores and assumed two fractal regimes. Nonetheless, Larson and Morrow (1981) indicated that the capillary pressure curve depends not only on the geometrical and wetting properties of individual pores but also upon the pores’ connections to the surface of the sample. Hunt et al. (2014) also argued that some differences between the modeled and measured results must be due to the fractal models having no percolation concepts built in. Models such as RS (Rieu and Sposito, 1991) are in some sense no different from the capillary bundle model: they do not consider pore connectivity, which is one of the most important characteristics of porous media. Sahimi et al. (1990) classified models for transport and reaction in porous media as continuum bundle model and discrete percolation network model. The former is not suitable for describing the phenomena in which the connectivity of the pore space or the fracture network is significant. The latter is free of the limitations of the continuum model. It can be applied to describe phenomena from microscopic to macroscopic scales. Percolation theory is complementary to the network models. The first capillary model using percolation theory was developed by Larson and Morrow (1981), in which the physical concept of accessibility that determines which pores can be invaded by nonwetting phase, and from which pores it can be withdrawn, has been explicitly utilized. Then Heiba et al. (1992) represents a refinement of the model of Larson and Morrow (1981). The basis for Heiba model is that during injection and retraction, the spatial distributions of the pores accessible to and occupied by nonwetting phase, which they refer to as the sub-distributions, are not identical. Consequently, the throat size distribution of the subset of pore space occupied by non-wetting phase differs from the overall throat size distribution. Heiba et al. derived analytical formulae for such sub-distributions. Soll and Celia (1993) defined a set of rules in a network model and developed a computer code (3PSAT) to calculate three-phase capillary-pressure curves. Hunt (2004) presented a capillary pressure model using the percolation theory applied to the fractal of porous media for the first time. The model is different from other models in the following two aspects. First, the particularity of fractal distribution is considered in this model. He presented that using the water content and threshold water content (in the soil field) to replace the allowed fraction p and percolation threshold pc in conventional percolation theory directly. Second, non-equilibrium process is considered in this model rather than the equilibrium one in conventional percolation model. In this paper he points out that the model can be well applied both in drainage and imbibition of soil porous media.

The new fractal capillary pressure model proposed in this paper integrated all the advantages of the methods above. Conceptually, we follow the model developed by Heiba et al. (1992). The difference lies in that the allowed saturation Sp and threshold saturation Sc is used in this paper, instead of the allowed fraction p and percolation threshold pc in conventional percolation theory. The saturation is in analog with the pore fraction ranging between zero and one, while applying the saturation avoids the prominent influence of the minimum pore radius, which is difficult to determine and has little effect on macroscopic properties. The equilibrium method is used in this paper. A simple and yet sufficiently accurate expression is employed to calculate the accessible fraction. The saturation is a function of pore distribution. And so the saturationcapillary pressure model of fractal porous media can be obtained. Finally, parameter sensitivity analysis and comparison with the results from conventional fractal theory are performed. The new model is validated by experimental data.

2. The characteristics of fractal porous media Since the fractal theory was proposed by Mandelbrot (1982), many strict synthetic fractal objects have been proposed by mathematician, which greatly improve our understanding of the fractal media. Among the many synthetic fractal objects, the 2D Sierpinski carpet and its 3D equivalent, the Menger sponge (Fig. 1) have been widely used to model porous materials. The Menger sponge is constructed by starting with a cube of size L (the “initiator length”). For the first iteration, the operational length scale r is (L=b) in which the scaling factor b ¼ 3 for the traditional Menger sponge. Divide every face of the cube into b2 squares, like a Rubik's Cube. This will sub-divide the cube into b3 smaller cubes. Remove n smaller cubes of size r, leaving (b3  n) smaller cubes. Then, repeat the former step for each of the remaining smaller cubes, and continue to iterate ad infinitum. The number of remaining cubes of size r satisfies the following scaling relation (Ghanbarian-Alavijeh et al., 2011):

Nð¼ rÞ ¼ ar D

(1)

where a is a constant coefficient, and D is the fractal dimension. Following Friesen and Mikula (1987), let a ¼ 1 and obtain the fractal dimension of the sponge as:

D¼


 log b3  n log b

(2)

The remaining and the removed cubes in Fig. 1 represent solid particles and pores, respectively. Although the pore phase in the sponge is not geometrically fractal, its number-size distribution is given by a power-law function and it has the same fractal dimensionality as the solid phase (Rieu and Perrier, 1997):

Nð rÞ ¼ kr D

(3)

where Nð rÞ is the number of fractal objects whose size is equal to or greater than r, k is a constant coefficient, and the fractal dimension D typically ranges between 0 and 3 in natural porous materials. So in fact, one fractal dimension scales both solid and pore phases, albeit in different ways. For natural porous media, it only applies within the specified range, rmin to rmax , in which the power-law function (Eq. (3)) is truncated. And natural porous media are more heterogeneous and complex with randomly self-similar instead of the exactly selfsimilar as shown in Fig. 1. The probability density function of fractals-the number of

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

Z Sw ¼ Z

r

Gr 3 f ðrÞdr

rmin rmax

rmin

¼ Gr3 f ðrÞdr

9

l r l  rmin

rl

max

l  rmin

(7)

where l ¼ 3  D, and G is the shape factor. For a cylindrical pore with radius r and length 2r, G ¼ 2p. Note that Snw þ Sw ¼ 1. For the majority of natural porous media, the maximum pore size is much larger than the minimum (rmax [rmin ). Therefore, Eq. (7) can be simplified as:

 Sw y

r

l (8)

rmax

Substituting Eq. (5) into Eq. (8) leads to the De Gennes model (1985) model:

Sw ¼ Fig. 1. Exactly self-similar Menger sponge after three iterations with fractal dimension D ¼ 2:7268 and b ¼ 3.

objects whose size is within the range from r to (r þ dr) is proportional to the first derivative of Eq. (3):

f ðrÞ ¼ cr 1D ; rmin  r  rmax

(4)

where f ðrÞ is the probability density function, and c is a constant coefficient which can be found by taking the integral of Eq. (4) from D  r D Þ. If rmin to rmax and setting it equal to 1: c ¼ D=ðrmin max rmax [rmin , as occurs in many natural porous media, the constant D . coefficient c ¼ Drmin Fig. 2 shows the pore size distribution of fractal sponge with upper limits rmax ¼ 100mm and lower limit rmin ¼ 0:01mm. In the log-log plot, the logarithm of pore size distribution function f ðrÞ linearly decreases with the logarithm of pore size r, as indicated in Eq. (4).

Assuming porous medium is formed by cylindrical tubes with radius r, the required capillary pressure Pc for fluid to enter the pores should obey the Young-Laplace equation:

Pc ¼

lnðSw Þ ¼ l lnðPc Þ þ l lnðPe Þ

(5)

where s is the interfacial tension between the wetting and nonwetting phase. q is the wetting contact angle. In both drainage and imbibition processes, if the capillary bundle model is used, the non-wetting phase can enter and take up all the pores whose radii greater than r corresponding to the displacement pressure Pc . While the rest pores are occupied by wetting phase. So the saturation of non-wetting and wetting phases can be obtained through the individual phase volume divided by the total pore volume:

Zrmax Gr3 f ðrÞdr

rmin

¼ Gr 3 f ðrÞdr

l rmax  rl

l l rmax  rmin

(6)

(10)

Eq. (10) shows a linear relationship between lnðSw Þ and lnðPc Þ. The slope of the straight line lnðPc Þ versus lnðSw Þ is ðD  3Þ. This relationship is the key criterion to judge if a porous medium satisfies fractal distribution, and to obtain the fractal dimension. It will be applied to the derivation of the fractal-percolation-based model for capillary pressure in the next section. For the porous media with narrow pore size distributions, the condition (rmax [rmin ) does not hold. In this case, substituting Eq. (5) into Eq. (7) results in the general form of relation between the saturation and capillary pressure as: l Pcl  Pmax  l l Pe  Pmax

(11)

where Pmax is the maximum capillary pressure, corresponding to the minimum pore size rmin . Eq. (11) can be rearranged into the following form in line with Li model (2010):

Pc ¼ Pmax ð1  hSw Þl

1



2s cos q r

Snw ¼ Z r rmax

(9)

where Pe is the entry pressure, corresponding to the capillary pressure at the maximum pore size rmax according to Eq. (5). Take the logarithm of both sides of Eq. (9) yields

Sw ¼ 3. Conventional fractal model for the capillary pressure

 l Pe Pc

where h ¼ 1 

Pe Pmax

(12)

l .

4. Traditional percolation model Percolation theory aims to describe the morphology of, and transport through, randomly disordered media by statistical means. The theory pertains to network models that consist of branches (bonds or links) and nodes (sites or junctions) (Dullien, 1992; Sahimi, 1993; Stauffer and Aharony, 1994; Yuan et al., 2012). A simplification for flow through these types of networks is a model in which all the properties of the pore space are assigned to the bonds. Bond percolation, applied in the present work, assumes that the porous medium consists of sets of pore throats that are connected to points with zero volume. The number of pore throats connected to the points is called the coordination number Z. The nodes become essentially zero-dimensional features with no role except as markers for branching. A bond is occupied with the probability p and is empty with the probability ð1  pÞ. At some

10

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

Fig. 2. Density function character of pore distribution in fractal porous media.

well-defined value of p, there is a transition in the topological structure of a random network from a macroscopically disconnected structure to a connected one. This value of probability is defined as the bond percolation threshold pc (Stauffer and Aharony, 1994). Fig. 3 is an illustration of a 3D uniformly distributed bond percolation network. The cross points of grids are nodes and lines between nodes are bonds. The volume of nodes is neglected. Fluids of different phases can pass through nodes simultaneously. Bonds are denoted as pore throats. Each pore throat is randomly determined as a bond in the network. Sizes of different pores are represented by lines with different colors and width. Pore throats (bonds) are randomly located following the known distribution. Therefore, a flow network model can be built up. As an example, the pore network in Fig. 3 has the coordination number Z ¼ 6. The percolation threshold can be calculated using the model from Stauffer and Aharony (1994) as

pc ¼

d Zðd  1Þ

(13)

where d is Euclidean dimension. In the example of pore network shown in Fig. 3, d ¼ 3. Therefore, the percolation threshold pc ¼ 0:25 according to Eq. (13). The fractions of allowed pores p and accessible pores c are introduced in percolation theory (Sahimi, 2011). “Allowed” means that we count the bonds that can in principle be invaded by a phase (i.e., the capillary pressure will not prevent this), but ignore the possibility that the surrounding bonds are too small such that the invading fluid will never have the chance to enter. “Accessible” means that the bonds are allowed by the phase, and the surrounding bonds will not inhibit the invading fluid to enter. In Fig. 4, red dots are nodes. Blue lines with various widths are bonds. The radius of bond 4 is larger than the other bonds connecting to it. Considering drainage process, when pressure just reaches the entry pressure value required for entering bond 4, it is not enough for the fluid to enter. Therefore, bond 4 is allowed but inaccessible to the fluid.

Deriving (or estimating) accessibility functions and the associated percolation thresholds is a basic task of percolation theory (Stauffer and Aharony, 1994). An accessibility function, cðpÞ, depends only on network coordination, or topology, and the fraction p of allowed bonds. Moreover, when the fraction of allowed bonds falls below the bond percolation threshold pc , the allowed bonds exist exclusively as isolated clusters that correspond to a hydraulically isolated stationary phase. Hydraulically isolated means that there are no flow paths occupied by that phase and spanning the sample except for thin films and pendula structures retained along pore walls. Above pc , the isolated clusters merge and a samplespanning transport path appears. Finally, as p approaches unity, the probability that any given allowed bond is part of a samplespanning cluster approaches unity. Thus,



c ¼ 0 p < pc c ¼ p p/1

(14)

The equation above is a qualitative expression only. A more accurate quantitative expression is needed in practice. Among all the network models, the simplest is the Bethe lattice network, which has been widely used in statistical mechanics for investigating critical phenomena in the mean-field approximation (Sahimi, 2011). Bethe network, characterized by the coordination number Z, is an endlessly branching structure with no closed loops. It is assumed that the porous space can be represented by interconnected bonds (pores) in Bethe lattice. Based on analytical expressions derived for Bethe network, the accessible fraction cðpÞ as a function of the allowed fraction p can be obtained using numerical procedure as (Fisher and Essam, 1961; Heiba et al., 1992):

8 > ! >   * 2Z2 > Z2 > >

> > > > : 0 p  pc

p > pc

(15)

where p* is the solution to the following transcendental equation:

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

11

Fig. 3. An illustration of a 3D uniformly distributed bond percolation network.

Fig. 4. Schematic of allowed but inaccessible fraction.


Z2 p* 1  p*  pð1  pÞZ2 ¼ 0

(16)

Specifically, p* ¼ 0 as p ¼ 1. The accessible fraction cðpÞ given by Eq. (15) is based on analytical solution using Bethe lattice network, and is inconvenient to implement. Alternatively, a simple formula proposed by Bedrikovetsky and Bruining (1995) has the following explicit form:

8 0 p < pc > > > > b <  cðpÞ ¼ psc p  pc > psc  pc > > > : p p > psc

pc  p  psc

(17)

where psc ¼ 1:65pc . The exponent b ¼ 0:4 is a universal constant for 3D network. Although the pore structure model established by Bethe network is significantly different from the pore structure of actual porous media, as reviewed by Sahimi et al. (1990), the Bethe network approximation to estimate the transport properties of pore networks has been widely used in the statistical physics of disordered media and has succeeded remarkably in explaining basic aspects of multiphase flow in porous media. 5. Fractal-percolation-based model The fractal distribution is monotonically decreasing type with a

span of several orders of magnitude, which contains the upper and lower limits of the pore sizes. The lower limit is the most important parameter of fractal distribution when fractal dimension D is a certain value. But the lower limit is usually difficult to determine. On the other hand, it is the maximum radius rather than the minimum radius that mainly affects the macroscopic property of porous media such as porosity, permeability, capillary pressure and so on. Traditional percolation parameters such as allowed fraction and accessible fraction are based on the percentage of cumulative pore number satisfying the size limit, which is a function of pore size distribution and is also affected strongly by the minimum radius. Hunt (2004) developed a new method to solve the problem. He proposed the water content and threshold water content (in the soil field) to replace the allowed fraction p and percolation threshold pc . The water content as a macroscopic property can avoid the major influence of rmin . In the pore field, the saturation S is more commonly used. In this paper, we will use saturation S as basis parameter to deduce the capillary pressure model. Firstly, some parameters should be defined, Sp means allowed saturation. In drainage process, the entering sequence is from large pores to small pores, therefore Sp is the non-wetting phase saturation, assuming non-wetting phase can enter all the pores whose radius is equal to or greater than the radius limit r:

Zrmax Gr3 f ðrÞdr Sp ¼ Z r rmax rmin

¼ Gr3 f ðrÞdr

l rmax  rl

l l rmax  rmin

y1 

 l Pe Pc

(18)

Similarly, the infinite cluster can be formed when allowed saturation Sp reaches a certain value, which is called the threshold saturation Sc . Sc can be obtained from drainage experiment. When Sp is equal to Sc , solving Eq. (18) results in a critical value of capillary pressure, which is the minimum pressure required for the non-wetting phase to form an infinite cluster through the porous media. This critical value is called breakthrough pressure Pb (Dullien, 1992):

12

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

6. Experimental measurements

Pb ¼ Pe ð1  Sc Þl

1

(19)

Define Sa as the accessible saturation. Sa is the real non-wetting phase saturation that enters the pores, at the moment when Snw ¼ Sa . In continuum percolation, the allowed saturation Sp is the variable analogous to the allowed fraction p, and the accessible saturation Sa is analogous to the accessible fractionc. Consequently, Sa and Sp have the similar relationship as c and p in Eq. (14):



Sa ¼ 0 SP < Sc Sa ¼ SP SP /1

(20)

The formula above shows when Sp /1, Sa ¼ Sp . It indicates that nearly all the allowed pores are accessible to the non-wetting phase fluid. Actually, Sa and Sp only differ near the percolation threshold. Larson et al. (1980) reported that Sa ySp as long as Sp  2Sc , resulting from their simulation. In analog with the expression of the accessible fractioncðpÞ in Eq. (17), the accessible saturation Sa as a function of the allowed saturation Sp can be expressed as:

8 0; SP < Sc > > >  >  < S  Sc b Sa ðSP ; Sc Þ ¼ Ssc P > Ssc  Sc > > > : SP SP > Ssc

Sc  SP  Ssc

(21)

where Ssc ¼ 1:65Sc , b ¼ 0:4. Define Psc as the pressure that when Sp is equal to Ssc . According to Eq. (18), we can get:

Psc ¼ Pe ð1  1:65Sc Þl

1

(22)

Fig. 5 shows the relation between the accessible saturation and the allowed saturation calculated from Eq. (21). When the allowed saturation Sp is lower than the percolation threshold Sc , infinite cluster does not exist. At this moment, accessible saturation Sa is 0; When Sp is just above Sc , a small number of large pores are surrounded by large number of small pores as shown in Fig. 4, which means not all the allowed pores are accessible. At this stage, Sa is an exponential function of ðSp  Sc Þ, as shown in Eq. (21). When Sp is much larger than Sa , majority of the nodes in the network are occupied by non-wetting phase. Almost all the allowed pores are accessible for non-wetting phase, i.e., Sa ¼ Sp . For drainage process, the accessible saturation is the nonwetting phase saturation that enter the porous media, i.e., Snw ¼ Sa . So the wetting phase saturation satisfies:

Sw ¼ 1  Sa

(23)

Substituting Eqs. (18), (19), (22) and (23) into Eq. (21) leads to the new fractal percolation capillary model as:

8 > 1 Pc < Pd > > > >  l 1 > > 0 b > > Pe > 1  S  > c Pc > < B C A Sw ¼ 1  Ssc @ Ssc  Sc > > > > > > >  l > > Pe > > Pc > Psc > : Pc

Pd  Pc  Psc

(24)

It can be seen from Eq. (24) that the saturation-capillary pressure relation of fractal porous media is mainly influenced by the maximum pore radius rmax , fractal dimension D and the threshold saturation Sc .

Mercury intrusion tests were conducted in three sandstone core samples from an oil reservoir. The error of the measurements on pressure and volume was about 1%. The reason for the selection of mercury intrusion tests to measure capillary pressure curves is that the surface tension of mercury and the contact angle are well known and constant during an experiment. With this feature, the measured capillary pressure curves should be a good representative for the pore structure of rock. The measured porosity and permeability data of the three reservoir core samples are listed in Table 1. The three cores were sampled from different depths in the oil reservoir and were expected to have different pore structures and heterogeneity. The surface tension of air/mercury is 480 mN=m and the contact angle through the mercury phase is 140 according to the results reported by Purcell (1949). The mercury intrusion experimental results of the three cores are shown in Fig. 6, which shows similar trend for different cores. It means a minimum mercury pressure is needed to inject into the core, which is the entry pressure Pe . When the pressure is lower than this value, mercury cannot enter the core due to the size limit and the connectivity relation. So, all of the data points at Sw ¼ 1 in Fig. 6 are out of modeling range for a capillary pressure curve. The breakthrough pressure Pb is usually higher than Pe . When the pressure slightly exceeds Pb , a large amount of mercury begins to enter into pores. It yields a significant decrease of the wettingphase saturation over a small variation in capillary pressure. This appears as a plateau-like trend on the curve (Sakhaee-Pour and Bryant, 2015). Afterwards, with the further pressure rise, the speed of mercury injection decreases and the curve becomes steeper (Fig. 6). 7. Results and discussions Fig. 7 shows the relationship between lnðPc Þ and lnðSw Þ for the experimental data of the three oil reservoir sand samples, which should be linear if the pore size distribution is fractal, according to Eq. (10). The linear relationship is established for the three cores except for the data out of modeling scope, i.e., the data points with lnðSw Þ ¼ 0 when the pressure is below the entry pressure (Fig. 7). Another exception is a number of data points at the beginning of mercury entering into the core. These data points distinctly deviate from the line. According to the linear relationship equation, the fractal dimension D and the entry pressure Pe can be calculated (see Table 1), and the saturation function can be obtained from Eq. (9). The linear relationship between lnðPc Þ and lnðSw Þ corresponds to the case of constant fractal dimension in fractal model. This is applicable to a large range of porous media. For some other media which have piecewise linear relationship between lnðPc Þ and lnðSw Þ such as coal bed (Friesen and Mikula, 1987; Liu et al., 2015), the proposed fractal-percolation-based model can be extended for this case, by introducing piecewise constant fractal dimensions into the model. The comparison between the experimental data and the conventional fractal model is shown in Fig. 8 with dots and dashed lines. It shows that the calculated entry pressure Pe and its neighborhood values evidently deviate from the measured data points. But when the pressure goes up, the modeled results agree with the measured data. In comparison, the calculated results from the new fractal-percolation-based model are also shown in Fig. 8 (the solid curves). The new model employs Pe and D calculated by the conventional fractal model. Threshold saturation Sc is calculated using trial and error method. The values of fitting parameters are shown in Table 1. Fig. 8 indicates that the new model can be applied to

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

13

Fig. 5. Relation between allowed fraction and accessible fraction.

Table 1 Measured data and modeling results. Sample

fð% Þ

kðm d Þ

S1 S2 S3

13.22 19.95 30.10

0.185 20.900 490.000

Conventional model

New model

D

Pe

D

Sc

Pe

Pb

2.7924 2.6886 2.7918

0.5000 0.0625 0.0068

2.7924 2.6886 2.7918

0.19 0.30 0.36

0.5000 0.0625 0.0068

1.3751 0.1965 0.0580

calculate the breakthrough pressure Pb according to Eq. (19) efficiently, and can better fit the low-pressure data by accounting for the percolation threshold Sc . The successful fitting result shows that pore connectivity has great influence on capillary pressure curve, especially near the threshold saturation. In contrast, the conventional fractal model ignores the important character, which is the reason why obvious deviation occurs at the right end of the curve. So this deviation part

Fig. 6. Capillary pressure curves of the three core samples from an oil reservoir.

14

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

Fig. 7. The relation between lnðPc Þ and lnðSw Þ of three cores.

of pores still obey the same fractal distribution rule, rather than belong to Euclidean space, or obey another fractal distribution rule. The new model overcomes the limitations of the conventional model by introducing percolation concepts. There are two main parameters in conventional fractal theory to form a capillary pressure curve, one is the maximum pore radius rmax , which determines the entry pressure Pe . The other is the fractal dimension D, which affects the slope and shape of the capillary pressure curve. While in fractal percolation theory, besides rmax and D, the most important parameter is the threshold

saturation Sc , which is the function of size distribution and connectivity of pores. Here we emphasize that we ignore the impact of the minimum pores rmin because the inequality (rmax [rmin ) is satisfied in most natural porous media. The influence of Sc on the capillary pressure curve is shown in Fig. 9. The increase of Sc leads to the larger breakthrough pressure and the longer flat line. The red line represents the case Sc ¼ 0, which is equivalent to the conventional fractal model coupled with capillary bundle model without considering the pore connectivity. The capillary pressure curve calculated by the new model is

Fig. 8. Comparison between new and conventional fractal models based on the experimental data obtained from mercury injection. The solid line is the fitting line of new model, and the dashed line is that of conventional fractal model.

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

15

Fig. 9. The capillary pressure curves of fractal percolation model at different threshold saturation Sc (D ¼ 2:7rmax ¼ 10mm). The red line represents the case Sc ¼ 0, which is identical to the conventional fractal model.

consistent with the experimental results, i.e., there exists a breakthrough pressure. When the pressure reaches this value, a large number of pores start to connect. This shows a sudden mercury intrusion amount in the mercury intrusion experiment. In twophase drainage experiment, the phenomenon occurs as the driving phase begins to form the infinite cluster and the waterflooding begins, i.e., the driving phase enters the porous media and starts to flow. The large pore distribution of porous media has great influence on the parameters such as porosity and permeability. It demonstrates that accurate description of large pore structure is very important. Through the model proposed in this paper, the fractal model can be extended from special form to general case. The thorough understanding of the relation between pore structure and two-phase flow will greatly expand the application of fractal theory in the description of pore structure.

8. Conclusions A new fractal-percolation-based capillary pressure model is developed. It enables quantify the percolation threshold of heterogenous porous media thus better matching experimental capillary pressure data. The capillary pressure curve resulting from the conventional fractal theory, which corresponds to zero threshold, is a particular case of the new model proposed. The effect of percolation threshold reflects the impacts of pore connectivity to the shape of capillary pressure curve, which highly improves the prediction accuracy at the breakthrough pressure neighborhood. The new fractal percolation model in this paper combines the fractal theory with percolation model, accounting for the pore connectivity, which is one of the most important characteristics of porous media. By applying the analytical model to drainage process, the entire capillary pressure curve can be more accurately predicted. This new model is a general extension to the conventional fractal theory.

Acknowledgement Supports from the projects National Natural Science Foundation of China (41302116), Educational Commission of Sichuan Province of China (14ZB0438), and Open Fund of State Key Laboratory of Oil and Gas Reservoir Geology and Exploitation (Chengdu University of Technology) (PLC201410) are gratefully acknowledged. Useful suggestions given by Prof. P. Bedrikovetsky (University of Adelaide) are also acknowledged. Nomenclature a b c d D f ðrÞ G h k L n N Nð rÞ p pc Pb Pc Pe Psc r S Sa Sc

constant scaling factor constant Euclidean dimension fractal dimensionality of pore space distribution function of radius shape factor constant constant cube length number of cubes number of fractal object (pores or solid) number of fractal objects with sizes equal to or greater than r fraction of allowed pores percolation threshold breakthrough pressure capillary pressure entry pressure universal constant radius saturation accessible saturation threshold saturation for percolation

16

Sp Ssc Z

J. Zheng et al. / Journal of Natural Gas Science and Engineering 41 (2017) 7e16

allowed saturation universal constant coordination number

Greek

b q l m s c

universal constant contact angle

constant viscosity interfacial tension fraction of accessible pores

Subscripts max maximum min minimum nw non-wetting phase w wetting phase x x coordinate y y coordinate z z coordinate References Angulo, R.F., Gonzalez, H., 1992. Fractal Dimensions from Mercury Intrusion Capillary Tests. SPE 23695, Presented at the 2nd Latin American Petroleum Engineering Conf. of SPE Held in Caracas, Venezuela, March 8e11. Bedrikovetsky, P., Bruining, J., 1995. A percolation based upscaling technique for viscous force dominated waterflooding in uncorrelated heterogeneous reservoirs. In: Proceedings of the 8th European Symposium on Improved Oil Recovery. European Association of Petroleum Geoscientist & Engineers (EAPG), Vienna, Austria, pp. 15e17. Brown, H.W., 1951. Capillary pressure investigations. J. Petroleum Technol. 3 (03), 67e74. Bu, H., Ju, Y., Tan, J., Wang, G., Li, X., 2015. Fractal characteristics of pores in nonmarine shales from the Huainan coalfield, eastern China. J. Nat. Gas. Sci. Eng. 24, 166e177. Cihan, A., Perfect, E., Tyner, J.S., 2007. Water retention models for scale-variant and scale-invariant drainage of mass prefractal porous media. Vadose Zone J. 6, 786e792. de Gennes, P.G., 1985. Partial filling of a fractal structure by a wetting fluid. In: Adler, D., Fritzsche, H., Ovshinsky, S.R. (Eds.), Physics of Disordered Materials. Plenum, New York, pp. 227e241. Deinert, M.R., Parlange, J.-Y., Cady, K.B., 2005. Simplified thermodynamic model for equilibrium capillary pressure in a fractal porous medium. Phys. Rev. E 72, 041203. Dullien, F.A.L., 1992. Porous Media: Fluid Transport and Pore Structure, second ed. ACADEMIC PRESS, INC. Fisher, M.E., Essam, J.W., 1961. Some cluster size and percolation problems. J. Math. Phys. 2, 609e619. Friesen, W.I., Mikula, R.J., 1987. Fractal dimensions of coal particles. J. Colloid Interf. Sci. 120, 263e271.
n, H., Huang, G., 2011. A review of fractal, prefractal Ghanbarian-Alavijeh, B., Milla and poresolid-fractal models for parameterizing the soil water retention curve. Can. J. Soil Sci. 91, 1e14. Hansen, J.P., Skjeltorp, A.T., 1988. Fractal pore space and rock permeability implications. Phys. Rev. B 38 (4), 2635e2638. Heiba, A.A., Sahimi, M., Scriven, L.E., Davis, H.T., 1992. Percolation theory of twophase relative permeability. SPE Reserv. Eng. 7 (1), 123e132. Hunt, A.G., 2004. Continuum percolation theory for water retention and hydraulic conductivity of fractal soils: 2. Extension to non-equilibrium. Adv. Water Resour. 27, 245e257. Hunt, A., Ewing, R., Ghanbarian, B., 2014. Percolation Theory for Flow in Porous Media, 3nd ed. SpringerLink. Katz, A.J., Thompson, A.H., 1985. Fractal sandstone pores implication for conductivity and formation. Phys. Rev. Lett. 54 (3), 1325e1328. Krohn, C.E., 1988. Sandstone fractal and Euelidean pore volume distribution.

Geophys. Res. 93 (B4), 328e329. Lai, J., Wang, G., 2015. Fractal analysis of tight gas sandstones using high-pressure mercury intrusion techniques. J. Nat. Gas. Sci. Eng. 24, 185e196. Larson, R.G., Striven, L.E., Davis, H.T., 1980. Percolation theory of two phase flow in porous media. Chem. Engng Sci. 36, 57e73. Larson, R.G., Morrow, N.R., 1981. Effects of sample size on capillary pressures in porous media. Powder Technol. 30, 123e138. Li, K., 2010. More general capillary pressure and relative permeability models from fractal geometry. J. Contam. Hydrol. 111 (1), 13e24. Liu, C.J., Wang, G.X., Sang, S.X., Gilani, W., Rudolph, V., 2015. Fractal analysis in pore structure of coal under conditions of CO2 sequestration process. Fuel 139, 125e132. Liu, P.C., Yuan, Z., Li, K., 2016. An improved capillary pressure model using fractal geometry for coal rock. J. Pet. Sci. Eng. 145, 473e481. Mandelbrot, B.B., 1982. The Fractal Geometry of Nature. W.H. Freeman, San Francisco.
n, H., Gonz
Milla alez-Posada, M., 2005. Modelling soil water retention scaling: comparison of a classical fractal model with a piecewise approach. Geoderma 125, 25e38. Naik, S., You, Z., Bedrikovetsky, P., 2015. Rate enhancement in unconventional gas reservoirs by wettability alteration. J. Nat. Gas Sci. Eng. 26, 1573e1584. O'Meara Jr., D.J., Hirasaki, G.J., Rohan, J.A., 1988. Centrifuge Measurements of Capillary Pressure: Part1-outflow Boundary Condition. SPE18296, 63rd Annual conf., Houston, TX, Oct. 2-5, p. 14. Perfect, E., 2005. Modeling the primary drainage curve of prefractal porous media. Vadose Zone J. 4, 959e966. Perrier, E., Rieu, M., Sposito, G., de Marsily, G., 1996. Models of the water retention curve for soils with a fractal pore size distribution. Water Resour. Res. 32, 3025e3031. Pfeifer, P., Avnir, D., 1983. Chemistry non-integral dimensions between two and three. J.Chem Phys. 7 (7), 3369e3558. Purcell, W.R., 1949. Capillary pressures-their measurement using mercury and the calculation of permeability. Trans. AIME 186, 39. Rieu, M., Sposito, G., 1991. Fractal fragmentation, soil porosity, and soil water properties: I. theory ,II. applications. Soil Sci. Soc. Am. J. 55, 1231e1244. Rieu, M., Perrier, E., 1997. Fractal models of fragmented and aggregated soils. In: Baveye, P., et al. (Eds.), Fractals in Soil Science. CRC Press, Boca Raton, pp. 169e201. Sahimi, M., 1993. Flow phenomena in rocks: from continuum models to fractals, percolation, cellular automata, and simulated annealing. Rev. Mod. Phys. 65 (4), 1393e1534. Sahimi, M., 2011. Flow and Transport in Porous, Media and Fractured Rock, second ed. WILEY-VCH Verlag GmbH & Co. KGaA. Sahimi, M., Oavalas, O.R., Tsotsis, T.T., 1990. Statistical and continuum models of fluid-solid reactions in porous media. Chern. Eng.Sci. 45, 1443e1502. Sakhaee-Pour, A., Bryant, S.L., 2015. Pore structure of shale. Fuel 143, 467e475. Sakhaee-Pour, A., Li, W., 2016. Fractal dimensions of shale. J. Nat. Gas. Sci. Eng. 30, 578e582. Shen, P., Li, K., 1994. A New Method for Determining the Fractal Dimension of Pore Structures and its Application. Proceedings of the 10th Offshore South East Asia Conference, Singapore, December 6e9. Soll, W.E., Celia, M.A., 1993. A modified percolation approach to simulating threefluid capillary pressure-saturation relationships. Adv. water Resour. 16 (2), 107e126. Stauffer, D., Aharony, A., 1994. Introduction to Percolation Theory, 2nd revised edn. Taylor and Francis, London. Tyler, S.W., Wheatcraft, S.W., 1990. Fractal processes in soil water retention. Water Resour. Res. 26, 1047e1054. Wong, P., Howard, J., Lin, J.S., 1986. Surface roughening and the fractal nature of rocks. Phys. Rev. Lett. 57, 637e640. Xiong, J., Liu, X.J., Liang, L.X., 2015. Experimental study on the pore structure characteristics of the upper Ordovician Wufeng Formation shale in the southwest portion of the Sichuan Basin, China. J. Nat. Gas. Sci. Eng. 22, 530e539. Yuan, H., Shapiro, A., You, Z., Badalyan, A., 2012. Estimating filtration coefficients for straining from percolation and random walk theories. Chem. Eng. J. 210, 63e73. Zhang, S., Tang, S., Tang, D., Huang, W., Pan, Z., 2014. Determining fractal dimensions of coal pores by FHH model: problems and effects. J. Nat. Gas. Sci. Eng. 21, 929e939. Zhou, L., Kang, Z.H., 2016. Fractal characterization of pores in shales using NMR: a case study from the lower Cambrian Niutitang formation in the Middle Yangtze platform, southwest China. J. Nat. Gas. Sci. Eng. 35, 860e872.