RETRACTED: Tracer transport in naturally fractured reservoirs: Analytical solutions for a system of parallel fractures

International Journal of Heat and Mass Transfer 103 (2016) 627–634

Contents lists available at ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Tracer transport in naturally fractured reservoirs: Analytical solutions for a system of parallel fractures Mahdi Abbasi a, Mahdi Hossieni a, Mojtaba Izadmehr b, Mohammad Sharifi a,⇑ a b

Department of Petroleum Engineering, Amirkabir University of Technology (Polytechnic of Tehran), P.O. Box: 15875-4413, Tehran, Iran Department of Chemical and Petroleum Engineering, Sharif University of Technology, Tehran 11365-9465, Iran

a r t i c l e

i n f o

Article history: Received 25 April 2016 Received in revised form 16 July 2016 Accepted 21 July 2016

Keywords: Fractured reservoir Mass transfer Shape factor Advection Diffusion

a b s t r a c t In naturally fractured reservoirs, modeling of mass transfer between matrix blocks and fractures is an important subject during gas injection or contaminant transport. This study focuses on developing an exact analytical solution to transient tracer transport problem along a discrete fracture in a porous rock matrix. Using Gauss-Legendre quadrature, an expression was obtained in the form of a double integral which is considered as the general transient solution. This solution has the ability to account the following phenomena: advective transport in fractures and molecular diffusion from the fracture to the matrix block. Certain assumptions are made which allow the problem to be formulated as two coupled, onedimensional partial differential equations: one for the fracture and one for the porous matrix in a direction perpendicular to the fracture. Using the obtained analytical solution, tracer concentration in matrix block and fracture was calculated. The advective-diffusive equation in matrix and fracture was used for evaluation of the mass transfer shape factor. The derived analytical solution was used for analyzing early and late time periods of mass transport phenomenon in fractured porous media. Finally, validation of the analytical solution was done by comparing the obtained results with laboratory data adapted from a column tracer test conducted on a fractured till. Ó 2016 Elsevier Ltd. All rights reserved.

1. Introduction A discontinuity that divides a rock into several blocks is called a fracture which is formed due to existence of different stresses exposed to formations [1]. Depending on rock type, rocks may be flexible enough to resist fracturing or be brittle. In brittle rocks such as carbonated sedimentary samples, fracturing is caused by stress [2]. Fractured rocks referred to as double porosity systems are composed of two different characteristic media called fracture and porous matrix systems. In a typical double porosity system, fractures are supposed to have high permeability in the range of 0.5–2 Darcy, therefore they could provide the main flow in the reservoir. Matrix blocks however, have low permeability in the range of 0.1–2 millidarcy and provide the main storativity [3]. Different techniques are available in literature in which a fracture network is numerically generated whose characteristics are as close as possible to the real situation. There exist two extreme models to simulate such systems, the first one being totally deterministic and the second one being totally random [4]. In the case of deter⇑ Corresponding author. E-mail address: [email protected] (M. Sharifi). http://dx.doi.org/10.1016/j.ijheatmasstransfer.2016.07.078 0017-9310/Ó 2016 Elsevier Ltd. All rights reserved.

ministic simulation, Warren and Root [5] introduced the dual porosity model for fractured reservoirs, focusing on pressure depletion of fractured system under constant production rate. In dual porosity model, porous matrix is considered to be the source, providing fluid for fracture. Several studies have been reported in the literature for fluid transport in dual porosity model [6–8,9]. In a tracer injection with no adsorption and reaction in an isothermal condition, advective, dispersive and diffusive mechanisms determine the transfer process [10]. Various experiments have been conducted to investigate parameters that affect dispersion, such as fluid velocity through porous media, solute size, porosity and permeability [11–14]. In the case of existence of regions with different concentrations such as fracture and porous media, a mass flux happens because of the diffusion mechanism [10]. Hu and Brusseau [15] reported the transport of different solutes in saturated, constructed structured media. These experiments were carried out in aggregates, stratified and macro-porous systems. Gwo et al. [16] presented the structural information in predicting solute movement within porous medium. Tracer injection is another method that has been used to investigate diffusive mass transport. Maloszewski and Zuber [17]

628

M. Abbasi et al. / International Journal of Heat and Mass Transfer 103 (2016) 627–634

Nomenclature C e C C C D Df h h1/2 Pe qw r rw s

tracer concentration [M L
3] average tracer concentration over the domain [M L
3] tracer concentration in the Laplace domain tracer concentration [M L
3] effective matrix molecular diffusion coefficient [L2 T
1] effective fracture molecular diffusion coefficient [L2 T
1] repetitive element block height [L] half-aperture of the fracture [L] Peclet number injection rate per repetitive element [L
3 T
1] direction of the fracture radius from the wellbore [L] Laplace variable

derived the governing equations of tracer movement for a model of parallel fractures with equal space and width. To develop this model, they investigated and interpreted different conditions including short term, long term and intermediate term experiments with different matrix porosities. Maloszewski and Zuber [18] used multiple tracers with different diffusion coefficients in order to quantify the effect of the matrix diffusion on the overall system dispersion. Using this method, they introduced an additional parameter to match the experimental data. Later, they used a single fracture dispersion model for several tracer tests (Maloszewski and Zuber [19]). Abelin et al. [20] reported that in the case of significant physical no equilibrium processes in subsurface media, tracers with larger molecular diffusion coefficient will be lost to the matrix porosity related to tracers with smaller molecular diffusion coefficients. Jardin et al., [21] provided an improved understanding of processes in heterogeneous matrix-fracture media to resolve inadequate knowledge about the transport processes that control contaminant migration. Analytical and numerical modeling of diffusive mass transport in fractured media is another topic that has been investigated in literature. In order to study the fracture-matrix transport phenomenon, it is extremely simpler to consider a single fracture. This assumption has received considerable attention from researchers. Neretnieks [22] developed an analytical solution for transport in a fracture assuming negligible dispersion and diffusion along the fracture. Grisak and Pickens [23] calculated concentrations in both the fracture and the porous matrix using the finite element method. Their numerical approach could account any arbitrary boundary condition, but the effects of numerical dispersion were difficult to assess in this method. Another work was done by Tang et al. [24], who developed general analytical solutions for the problem of solute transport in fractures in the following cases: longitudinal mechanical dispersion in the fracture, molecular diffusion along the fracture axis, molecular diffusion from the fracture into the porous matrix, adsorption onto the face of the matrix and adsorption within the matrix. Sudicky and Frind [25] developed an exact analytical solution for the situation of transient contaminant transport in discrete parallel fractures. Because the solution was based on analytical inversion of the Laplace transform, numerical inversion problems have been resolved. Feenstra et al. [26] developed an analytical model for contaminant transport in a porous media with a planar horizontal fracture for a radial and semiinfinite system. The introduced model was able to evaluate the influence of matrix diffusion on the contaminant transport away from the injected well. Lee and Teng [27] conducted a theoretical model for the radionuclide transfer phenomenon in a single fracture covering

t z Greek @

r u

time [T] direction normal to the fracture

partial derivative shape factor [L
2] porosity, fraction

Subscripts D dimensionless f fracture i initial m matrix

the entire range of sorption properties of rock in a linear transient case. Chen and Li [28] developed an analytical solution for radionuclide transport in a system of parallel fractures considering the constant flux as the inlet boundary condition. The assumption of negligible longitudinal dispersion along the fractures and steadystate solutions were also considered in this study. Simulation of solute transport is difficult owing to existence of matrix diffusion, which led to some simplifications for formulations reported in literature. Carrera et al., [29] tried to resolve this problem by considering some simplified asymptotic properties for matrix diffusion. For example, matrix diffusion coefficient and block size were neglected for the late time period. Fleming and Haggerty [30] formulated and presented a physical-based mathematical model of diffusion in sediments with variable diffusivity. They investigated the multiple matrix diffusivity along with laboratory studies and concluded that as porosity increases, diffusivity increases. Becker and Shapiro [31] investigated the effect of advective heterogeneity by analysis of the breakthrough during their forced-gradient tracer tests under various hydraulic configurations on a fractured crystalline bedrock. They proposed a theoretical model of transport that predicts tailing behavior in experiments. Liu et al. [32] reported that the effective matrix diffusion coefficient is scale dependent and increases with the test scale. They developed an introductory explanation for this scale dependent treatment. Chu et al. [33] investigated the effects of scale and reactions on lumped mass transfer coefficient. They tried to provide improved knowledge of the nature of a lumped mass transfer coefficient in the non-aqueous liquid source zone. They found out that, mass transfer coefficient stabilizes at some values after a while. Zhou et al. [34] obtained the effective matrix diffusion coefficient. They presented and interpreted the values of the effective matrix diffusion coefficient adapted from various field tests. Mathias et al. [35] used tracer in fractured rocks to investigate the effects of heterogeneity in which the importance of several parameters were presented including the input function for tracer injection, the lateral travel time, characteristic fracture, matrix diffusion times and the Peclet number. Cihan and Tyner [36] developed analytical solutions for mass transfer in macropores assuming that solute transport within the macro-pore is governed by advection while solute transport within the porous matrix is governed only by radial diffusion. They provided approximate solutions for three cases and concluded that in low permeability matrix rocks and short tracer tests, proposed approximations are helpful. In this paper, an analytical solution for describing radial tracer transport in fractured media is presented for the first time where the injected tracer diffuses from the fracture into the adjacent

629

M. Abbasi et al. / International Journal of Heat and Mass Transfer 103 (2016) 627–634

matrix. By transforming advective-diffusive equation into Laplace space, tracer concentration in matrix block and fracture is calculated at different times and locations. Then, dimensionless shape factor for the mass diffusion process is derived by taking into account the advection term in fractures. For this derivation, both advective and diffusive mass transfer processes are considered. By investigation of the results achieved in this derivation, unsteady-state and pseudo-steady-state periods of the matrix– fracture interaction are specified. The governing conditions include one dimensional advective transport in the fracture and onedimensional diffusive transport (normal to the fracture) into the matrix. Finally, validation of the obtained solution is illustrated by comparing it with laboratory tracer test data.

2. Description of the physical system The system considered in this study consists of a thin horizontal fracture placed in porous matrix block (Fig. 1). The following assumptions are considered for this system: 1. The injected fluid velocity is assumed to be constant in the fracture. 2. The fracture’s length is much larger than its width. 3. The horizontal permeability of the fracture is high enough for the purpose of complete mixing at all times. 4. The matrix transport is mainly controlled by molecular diffusion because of matrix low permeability. 5. The transport of fluid along the fracture is much faster than its transport through the matrix porous media. 6. The tracer’s injection rate (qw) is constant. Based on the second and third assumptions, one-dimensional transport is generated along the fracture. Fourth and fifth assumptions infer that the transport within the matrix porous media is perpendicular to the fractures’ fluid transport direction. Considering these assumptions, the problem will be simplified into two one-dimensional systems that are perpendicular to each other, thus, making it easy to solve. A similar approach was used by Tang et al. [24] and Sudicky and Frind [25].

z=0

The processes that should be considered are advective transport in the fracture media, dispersion mass transport along the fracture, and diffusion transport from the fracture to the matrix porous media. Longitudinal mechanical dispersion accounts the effect of mixing on the direction of fracture because of the velocity profile (Taylor [37]) and the fracture wall roughness. On the other hand, the hydrodynamic dispersion term accounts the combination of mechanical dispersion and molecular diffusion (Bear [38]). It should be pointed out that however, adsorption within the matrix and on the fracture wall is considered separately since these two terms are different for various chemical properties of matrix surface. Freeze and Cherry [39] have discussed adsorption processes in detail. The aforementioned assumptions in fractured reservoirs are reasonable according to the geological structure, geomechanical stresses, and position of fractures that are hundreds of kilometers underground. For example, in a typical naturally fractured reservoir, the width of fracture located in the formation is at the scale of micrometer and it is obvious that investigation of fracture concentration changes along the z direction is not an important topic. 3. Mathematical model In order to develop our model, the governing equations for expressing mass transfer phenomenon in a matrix-fracture media are derived, then analytical solutions for fractured rocks system are introduced based on descriptions of the physical system and with the help of presented assumptions. Subsequently, the dimensionless shape factor for an element of the fractured rock system is obtained using the analytical solutions of the coupled partial differential equations for matrix-fracture media. 3.1. Governing equations for advective–diffusive transport The transport processes in the system of Fig. 1 can be described by two coupled, one dimensional equations, one for the fracture and one for the porous matrix. The coupling is provided by the continuity of fluxes and concentrations along the interface. The differential equation for the fracture can be obtained by balancing the total mass of contaminant in the fracture.

No Flux

No Flux

Matrix

Matrix Diffusive Mass Transfer

Diffusive Mass Transfer

Advective Mass Transfer

Advective Mass Transfer

Diffusive Mass Transfer

Fracture

Diffusive Mass Transfer

Matrix

Matrix No Flux

r

No Flux

Fig. 1. Schematic of matrix–fracture system’s elements.

630

M. Abbasi et al. / International Journal of Heat and Mass Transfer 103 (2016) 627–634

According to the assumptions for fluid injection into the fracture reservoir and considering high permeability of fracture system, the injected fluid starts to diffuse radially into the horizontal fractures and then mass transfer takes place linearly between matrix and fracture. Concentration gradient is much higher in vertical direction compared with radial direction. Regarding the fact that fracture permeability is of high value, concentration can move rapidly in the fracture radially, and then the new concentration acts as boundary condition for matrix block. Therefore, radial diffusion can be ignored. Writing material balance for the element shown in Fig. 1 gives the transport equations for matrix and fracture. For a single matrix block, it gives the following equation:

@ 2 C m @C m D ¼ @z2 @t

ð2Þ

where Cf is the concentration of tracer in fracture, Df is the effective fracture molecular diffusion coefficient, r is the radial distance from the wellbore, qw is the injection rate per element, um is the matrix porosity, uf is the fracture porosity, hm/2 is the matrix block half height, and h1/2 is the fracture half thickness. In all formulations, subscripts m and f denotes matrix and fracture, respectively (Sharifi Haddad et al. [40]). In Eq. (1), according to the mentioned assumptions for the model, mass transfer mechanism into the matrix is molecular diffusion control and is perpendicular to the direction of the fraction axis (along z axis). By solving Eq. (1), tracer concentration in the matrix could be calculated in different positions along z direction. Also, according to the mass transfer between matrix and fracture and material balance used in Eq. (2), the amount of mass transfer between matrix and fracture could be expressed through changes of matrix concentration versus time in matrix and fracture contact point (h1/2). The initial and boundary conditions for Eq. (1) are:

C m ðz; r; 0Þ ¼ C mi

C Dm ¼

C mi
C m C mi
C w

ð6Þ

C Df ¼

C mi
C f C mi
C w

ð7Þ

zD ¼

@C m ð0; r; tÞ ¼0 @z

ð3cÞ

The initial and boundary conditions for Eq. (2) are:

C f ðr; 0Þ ¼ C mi

ð4aÞ

C f ðr w ; tÞ ¼ C w

ð4bÞ

C f ðr ! 1; tÞ ¼ C mi

ð4cÞ

The mass transfer rate between the matrix block and fracture per unit volume of the matrix block can be calculated by:

 D @C m  qm ¼ hm =2 @z z¼h=2
h1=2

ðh1=2 Þ

ð5Þ

ð9Þ

2

h h1=2

ð10Þ

rD ¼

r h1=2

ð11Þ

In above formulations, subscripts D denotes dimensionless form of that parameter. Cmi is the initial tracer concentration inside the matrix block, Cw is the tracer concentration in the injection well, and z is the space coordinate in the rock matrix block. rD is the dimensionless radius, and hR is the ratio of the element height to the fracture half thickness. By applying above parameters, dimensionless transport equation in the rock matrix block is achieved as follow:

D

@ 2 C Dm @C Dm ¼ @t D @z2D

ð12Þ

Dimensionless parameters are used for the fracture to give the following equation:

   @C Df Df 1 @ @C Df Pe @C Df um @C m 
¼ rD
r D @r D @tD D r D @r D @r D uf @zD zD ¼ðhR =2Þ
1

ð13Þ

qw where Pe is the Peclet number at the well-bore: Pe ¼ 4ph1=2 u D. f

3.2. Solutions of the governing equations The average matrix block concentration over its height and fracture concentration in the Laplace domain are determined in this section by taking the Laplace transform from Eqs. (12) and (13). For matrix, taking Laplace transform on Eq. (12) and corresponding initial and boundary conditions leads to:

D ð3bÞ

ð8Þ

hR ¼

ð3aÞ

C m ðh=2
h1=2 ; r; tÞ ¼ C f ðr; tÞ

z h1=2 Dt

tD ¼

ð1Þ

where D is the effective matrix molecular diffusion coefficient, Cm is the concentration of tracer in porous matrix, z is the vertical depth, and t is the time. In the derivation procedure of analytical solution, twodimensional mass distribution in the porous matrix, concentration across the thickness of the fracture, mass influx on the top surface of the matrix, and longitudinal dispersivity in the fracture are neglected. Writing material balance for a single fracture gives [25]:

!   @C f Df @ @C f qw C f 1 @ r ¼
r @r 4puf h1=2 @t r @r @r   u hm =2 q j
m uf h1=2 m z¼h1=2

In order to reduce the complexity, dimensionless forms of parameters are used. Governing equations are therefore utilized in dimensionless form by defining parameters as follows:

@ 2 C Dm ¼ sC Dm
0 @z2D

ð14Þ

C Dm ðzD ; rD ; t D Þ ¼ C Df ðrD ; t D Þ t D > 0 zD ¼ hR =2
1

ð15aÞ

@C Dm ðzD ; r D ; t D Þ ¼ 0 t D > 0 zD ¼ 0 @zD

ð15bÞ

where s is the Laplace domain variable. The second boundary condition (Eq. (15b)) gives C1 = 0 and first boundary condition (Eq. (15a)) yields the following equation:

C Dm

pffiffi coshð szD Þ p ffiffi ¼ C Df coshð sððhR =2Þ
1ÞÞ

ð16Þ

The normally high permeability of the fracture leads to an advective dominated piston-like flow; therefore, the first term (diffusion in the fracture) on the right hand side of Eq. (13) is small compared to the second term (advective term), and can be ignored

631

M. Abbasi et al. / International Journal of Heat and Mass Transfer 103 (2016) 627–634

(Maloszewski and Zuber, [17]; Neretnieks, [22]). It is worth noting that far away from the injection point, the velocity decreases and diffusion plays an important role. Since the used Peclet number values are very large (e.g. Pe = 106), diffusion mechanism could be ignored along the fracture length [24–25]. The fracture transport equation is then given by:

 @C Df Pe @C Df um @C Dm  ¼

r D @r D @t D uf @zD zD ¼ðhR =2Þ
1

ð17Þ

sC Df

1 C Df ðrD ; t D Þ ¼ s

0

ð19aÞ ð19bÞ

ð20Þ

By solving the above ordinary differential equation (ODE) and by using the boundary conditions given in Eqs. (19a) and (19b), it could be concluded that:

1 s

"

 exp


pffiffi 1 u pffiffi ðs þ m s tanhððhR =2
1Þ sÞÞðr2D
r 2wD Þ 2Pe uf

#

Thus, the dimensionless concentration in fracture (CDf) is given in terms of its inverse transform L
1 as: T

C Df ¼ 0

pffiffi pffiffi L
1 fexp½
Y s tanhððhR =2
1Þ sÞg ds

ð22Þ

1

p

1 0

2

Z

p

ðr2D
r2wD Þ 2Pe

TP0

1 um 2 ðr
r 2wD Þ 2Pe uf D

ð23Þ ð24Þ

  Z t
a 1 expð
asÞf ðsÞ ¼ L
1 ½f ðsÞ ds t D > a s 0

  1 expð
asÞf ðsÞ ¼ 0 t D < a L s
1

X¼

e expðeR Þ

1 0

1

e

where

0

cosðeI Þds de

0

ð30Þ

s as:

expðeR Þ ½sinðeI ÞjT þ sinðXÞde



ð31Þ

pffiffi  1 X coshð szD Þ p pffiffi ¼ ð
1Þn ð2n þ 1Þ 2 coshð sððhR =2Þ
1Þ ððhR =2Þ
1Þ n¼0 " !# p2 ð2n þ 1Þ2  exp
t 2 4ððhR =2Þ
1Þ   ð2n þ 1ÞpzD  cos 2ððhR =2Þ
1Þ

ð32Þ

ð33Þ

The convolution theorem states that:

Z

t

F 1 ðsÞF 2 ðt
sÞds

ð34Þ

Substituting for C Df from Eq. (31) into Eq. (34) gives:

C Dm ¼

  Z 1 Z tD 1 X ð2n þ 1ÞpzD n ð
1Þ ð2n þ 1Þ  cos 2 2ððhR =2Þ
1Þ 0 0 ððhR =2Þ
1Þ n¼0 2



expðeR Þ sinðeI ÞjT þ sinðXÞ e " ! # p2 ð2n þ 1Þ2  exp
ðt
s Þ dsde D 2 4ððhR =2Þ
1Þ 1

Eq. (35) can be integrated with respect to

ð35Þ

s

 Z 1 1 X ð2n þ 1ÞpzD 1 ð
1Þn ð2n þ 1Þ  cos expðeR Þ 2ððhR =2Þ
1Þ 0 e ððhR =2Þ
1Þ n¼0 2 8 93  2 2 2 > > < p ð2nþ1Þ 2 sinðeI ÞjT
e2 cosðeI ÞjT = 4ððh =2Þ
1Þ 6 7 R i 7 h  2  i h 2  6 e4 p41ð2nþ1Þ4 2 2 2 6 4þ 0 0 p ð2nþ1Þ p ð2nþ1Þ > ;7  6 16ððhR =2Þ
1Þ4 : þ exp
4ððh =2Þ
1Þ2 t D  4ððh =2Þ
1Þ2 sinðX Þ þ e2 cosðX Þ > 7de R R 6 7 h h  2  ii 4 5 4ððhR =2Þ
1Þ2 p ð2nþ1Þ2 þ p2 ð2nþ1Þ2 1
exp
4ððh =2Þ
1Þ2 tD sinðXÞ

C Dm ¼

2

2

R

ð36Þ


1

p

T

Eq. (32) expresses the dimensionless concentration in a fracture with all processes considered. The fracture porosity, uf in the above derivations is the ratio of fracture void space to the fracture bulk volume, which is approximately equal to 1. To complete the derivation procedure, we will now find the solution for the dimensionless concentration in the porous matrix. Noting that the inverse of the hyperbolic cosine function in Eq. (16) becomes (Oberhettinger and Badii [42])

ð25Þ

The inverse of the exponential term appearing in Eq. (25) was evaluated by Skopp and Warrick [41], who presented a model for miscible displacement of solutes in soils containing a significant portion of the immobile liquid phase. We can utilize this inverse term:

pffiffi pffiffi  L exp½
Y s tanhððhR =2
1Þ sÞ Z 1 1 ¼ e expðeR Þ cosðeI Þ de

Z

  Y e sinhððhR =2
1ÞeÞ þ sinððhR =2
1ÞeÞ 2 coshððhR =2
1ÞeÞ þ cosððhR =2
1ÞeÞ



ð26Þ

ð29Þ

where

Eq. (24) is derived by making use of the fact that;

L
1

Z

0

T ¼ tD


ð28Þ

pffiffi pffiffi  L
1 exp½
Y s tanhððhR =2
1Þ sÞ

T

L
1 ½f 1 ðpÞf 2 ðpÞ ¼

where

Y¼

C Df ¼

L
1 ð21Þ

Z

¼

in Eq. (18) gives:

! pffiffi Pe @C Df um pffiffi s tanhððhR =2
1Þ sÞ C Df ¼
sþ r D @r D uf

C Df ¼

  Y e sinhððhR =2
1ÞeÞ þ sinððhR =2
1ÞeÞ 2 coshððhR =2
1ÞeÞ þ cosððhR =2
1ÞeÞ




Eq. (30) can be integrated with respect to

where rwD is the dimensionless wellbore radius. Using Eq. (16) to @C m @zD

2

Z

zD ¼ðhR =2Þ
1

C Df ðrD ; t D Þ ¼ 0 t D > 0 r D ! 1

e2 T



Substitution of Eq. (27) into Eq. (22) and interchanging the order of integration of s and e yields:

ð18Þ

t D > 0 r D ¼ r wD

substitute for the term

eI ¼

Y e sinhððhR =2
1ÞeÞ
sinððhR =2
1ÞeÞ 2 coshððhR =2
1ÞeÞ þ cosððhR =2
1ÞeÞ

C Df ¼

Taking the Laplace transform from Eqs. (17) yields:

 Pe @C Df um @C Dm  ¼

 rD @rD uf @zD 

eR ¼




ð27Þ

where

X0 ¼ X þ

ðr2D
r2wD Þ 2Pe

M. Abbasi et al. / International Journal of Heat and Mass Transfer 103 (2016) 627–634

Eq. (36) is the solution for the dimensionless concentration in the porous matrix. The average matrix concentration can be calculated by taking integration over Eq. (36), which gives:

em ¼ 1 C Vm

Z C m dV m

ð37Þ

Vm

As discussed earlier, during the pseudo steady-state transfer, the mass transfer is proportional to the concentration difference between matrix and fracture and the proportionality constant is related to the shape factor as given in following equation:

e m
Cf Þ qm ¼ rDð C

ð38Þ

In addition, the mass transfer between matrix and fracture can be obtained in terms of matrix concentration by the following equation:

qm ¼

em @C @t

em @C e m
Cf Þ ¼ rDð C @t

ð40Þ

The matrix concentration can be obtained by solving the mass diffusivity equation using its initial and boundary conditions. By using these equations, dimensionless shape factor can be derived as follows:

D 2 h1=2

e m =@t @C e m
Cf Þ ðC

ð41Þ

Substituting for the fracture half thickness (h1/2) in terms of the dimensionless repetitive element height (hR = h/h1/2) gives the dimensionless shape factor (rh2) as:

rh2 ¼ h2R

RD=1-Analytical Model RD=1000-Analytical Model

100 RD=8E+06-Analytical Model

10

1 1.E+04

1.E+05

1.E+06

1.E+07 1.E+08 Dimensionless Time

1.E+09

1.E+10

Fig. 2. Dimensionless shape factor versus dimensionless time for RD = 1, RD = 103, and RD = 8  106.

ð39Þ

Combining Eqs. (38) and (39) yields:

rD ¼

1000 Dimensionless Shape factor

632

e m =@t @C ðC m
C f Þ

ð42Þ

4. Results and discussion In this section, solutions at different times are evaluated using analytical method. Then, for investigation of advective–diffusive mass transfer in the matrix-fracture media, the effects of changing solution parameters on the tracer concentration are shown. The mass transfer inside the porous matrix is considered to be a pure diffusive process for which the depth of penetration inside matrix blocks depends on the matrix block boundary condition. It has been shown that for a short transient period with Dirichlet type boundary condition, there will be a fast penetration and for a long transient period with time dependent boundary condition, there will be a slow penetration. Transient stage (early time) refers to a period of time that the mass transfer shape factor has a transient behavior and changes with time. The fracture width and effective molecular diffusion coefficient in porous media are of small values (e.g. h1/2 = 10
5 m and D = 10
10 m2/s). According to the Peclet number definition in fracture, these parameters are in the denominator of the Peclet number, thus, making the value of Peclet number in the fracture very large (e.g. Pe = 106) [24]. In Fig. 2, the effects of radial distance and injection rate are shown by investigation of dimensionless shape factors for different dimensionless ratios (RD = 1, RD = 103 and RD = 8  106), where hR = 5  104 and RD = (rD2
rwD2)/2Pe is a dimensionless ratio that shows the effect of both the radial distance and rate of injection at the wellbore on the shape factor.

Two different cases are shown here; the first one is for small dimensionless ratios (e.g. RD = 1 and 103), which happen to high injection rates or near wellbore locations and the second one is for large values of dimensionless ratios that happen for low injection rates or locations far from the injection wellbore. Fig. 2 shows two different time regions for the dimensionless shape factor, unsteady-state and pseudo steady-state periods. During the unsteady-state period, the shape factor is not constant and varies with time while in pseudo steady-state period, the dimensionless shape factor remains constant with time. Moreover, the dimensionless shape factor is nearly proportional to t-1/2 during the tranD sient period. It approaches p2 for small dimensionless ratios while for large values of dimensionless ratios, its value is approximately 12 as shown in Fig. 2. In the second case, mass transfer shows longer transient time. It could be concluded that for distances far away from the injection well or at small Peclet numbers, porous matrix blocks have a delayed response to the injected tracer. The reason for a longer transient period in the second case is the fact that the tracer concentration inside the fracture (boundary condition of the porous matrix block) varies with time. For the first case of mass transfer however, the tracer concentration inside fractures is almost constant and therefore the transient period is shorter. Changing the value of Peclet number is proportional with changing the injection rate. Therefore, for investigating the effect of advective flow on the dimensionless shape factor, the injection rate could be changed while keeping radius constant. Here, the magnitude of Peclet number is increased up to 106. Eq. (42) shows the relation derived for the mass transfer dimensionless shape factor showing that as the injection rate at the wellbore increases, the radial distance from the injection point decreases. By increasing the injection rate, the Peclet number value increases and the dimensionless shape factor approaches p2. Fig. 3 shows the average tracer concentration profile in the matrix for different dimensionless time at RD = 1000. As shown in this figure, tracer concentration starts to increase at tD = 106 and has a maximum value at tD = 109, after which the mass transfer between matrix and fracture reaches a minimum value. Figs. 4 and 5 show the effect of the Peclet number on tracer concentration in matrix and fracture for different dimensionless times, respectively. These values have been calculated at a fixed location of (r2D
r2wD ) = 8  109. According to Fig. 4, it could be concluded that for large Peclet numbers, tracer saturation period is shorter compared with that of small Peclet numbers. Figs. 4 and 5 show that at a fixed radius, by increasing the Peclet number (increasing the advective flow), the time required for the arrival of a constant tracer concentration gets shorter.

633

M. Abbasi et al. / International Journal of Heat and Mass Transfer 103 (2016) 627–634

tD=106

tD=107

tD=108

tD=109

0.9

1.2 Pe=1E+06-Analytical Model

1

Pe=1E+05-Analytical Model

0.8

Pe=1E+04-Analytical Model

0.6 0.4 0.2

0.7 0.6 0.5 0.4 0.3 0.2 0.1

0 1.E+00

1.E+03

1.E+06

1.E+09

1.E+12

Dimensionless Time Fig. 4. Dimensionless Matrix tracer concentration versus dimensionless time in different Peclet number values at (r 2D
r 2wD ) = 8  109.

1.4 Dimensionless Tracer Concentration (Cf)

Analytical Solution laboratory-measured Fenstra etal.Model

0.8 Relative Concentration

Dimensionless Tracer Concentration (Cm)

Fig. 3. Matrix Tracer concentration profile in different dimensionless times at RD = 1000.

1.2 1

Pe=1E+06-Analytical Model Pe=1E+05-Analytical Model Pe=1E+04-Analytical Model

0.8 0.6 0.4 0.2 0 1.E+00 1.E+02 1.E+04 1.E+06 1.E+08 1.E+10 1.E+12 1.E+14 Dimensionless Time

Fig. 5. Dimensionless Fracture tracer concentration versus dimensionless time in different Peclet number values at (r 2D
r 2wD ) = 8  109.

Grisak et al. [43] reported an experimental study of tracer injection within a fractured cylindrical porous medium. The injected tracer was chloride and calcium. The diameter and length of the cylinder were 0.67 m and 0.76 m respectively. There were two sets of vertical fractures with about 4 cm spacing in each set. The tracer was injected into the column vertically and the tracer breakthrough time was measured in each of four effluent quadrants. Fig. 6 shows the relative output concentration of chloride for the column effluent. For simulation of the experimental results, reported parameters of Grisak et al. [44] have been considered.

0 0

1

2 Time(Days)

3

4

Fig. 6. Experimental data for chloride breakthrough in the fractured cylindrical porous media, the obtained analytical solution and Feenstra et al. model result [26].

These parameters are: velocity in the fracture V = 29.7 m/day; half-aperture size h1/2 = 20 lm; longitudinal dispersity aL = 0.04 m; matrix (diffusion) porosity u = 0.35; and solute diffusion coefficient D in the matrix = 5.0  10
7 cm2/s. The experimental data could be simulated analytically by applying Eq. (33) and neglecting longitudinal dispersity in the fracture. Fig. 6 shows the obtained results of this simulation in comparison with the experimental data in which a relatively good match is observed; therefore, it can be concluded that in this case, neglecting longitudinal dispersity in the fracture is a reasonable assumption. For the analytical solution, two-dimensional mass distribution in the porous matrix, concentration variation across the thickness of the fracture, mass influx on the surface of the matrix, and longitudinal dispersity in fracture are neglected. In comparison with the experimental results which take into account all the aforementioned factors, the presented analytical model shows a good agreement in the first two days. The maximum error between the experimental data and predicted results from the analytical model is 9% and happens at 2.6 days. It could be concluded that neglecting above assumptions are a reasonable consideration. Also Feenstra et al. model [26] results are included in Fig. 6 which is an analytical model for contaminant transport in a porous media with a planar horizontal fracture for a radial and semi-infinite system. This model underestimates the relative concentration compared with the analytical solution and laboratory-measured data. Thus it can be concluded that for the purpose of mass transfer between matrix blocks and fractures, the assumption of semi-infinite matrix block dimension is not a good approximation under the proposed conditions.

634

M. Abbasi et al. / International Journal of Heat and Mass Transfer 103 (2016) 627–634

5. Conclusions In the present study, for the problem of radial tracer transport in discrete, parallel fractures, a general analytical solution was developed. The obtained solution was then evaluated by Gaussian quadrature in the form of a double integral. Advection along the fractures length and molecular diffusion from the fractures into the porous matrix blocks were accounted by the introduced solution. Although the developed solution is based on limiting assumptions, relatively good agreement between laboratory data and analytical solution results, shows that this solution can capture the tracer diffusion with a good accuracy. Furthermore, it was shown that mass transfer dimensionless shape factor has a transient behavior during the unsteady state period and a constant value during the pseudo steady state period. This constant value is equal to p2 in high Peclet numbers (e.g. Pe = 106) for locations near the injection wellbore. For distances far from the injection point, the Peclet number was decreased and a value of 12 for mass transfer dimensionless shape factor was obtained. These two different values show that the mass transfer dimensionless shape factor depends on boundary condition. Although a number of assumptions were made that resulted in a highly idealized system, the analytical solutions provide a framework for studying the mechanics of transport in fractured porous media and can be regarded as a basic building block for solutions to problems involving more complicated systems. They will also be useful in verifying numerical models developed for this purpose. References [1] T. van Golf-Racht, D. Theodor, Fundamentals of Fractured Reservoir Engineering, vol. 12, Elsevier, 1982. [2] Ali M. Saidi, Reservoir Engineering of Fractured Reservoirs (fundamental and Practical Aspects), Total, 1987. [3] Ashok Kumar Belani, Estimation of Matrix Block Size Distribution in Naturally Fractured Reservoirs, Diss. Stanford University, 1988. [4] Pierre M. Adler, Jean-François Thovert, Valeri V. Mourzenko, Fractured Porous Media, Oxford University Press, 2012. [5] J.E. Warren, P. Jj Root, The behavior of naturally fractured reservoirs, Soc. Petrol. Eng. J. 3 (03) (1963) 245–255. [6] H. Kazemi et al., Numerical simulation of water-oil flow in naturally fractured reservoirs, Soc. Petrol. Eng. J. 16 (06) (1976) 317–326. [7] L.S.-K. Fung, Numerical Simulation of Naturally Fractured Reservoirs. Middle East Oil Show, Society of Petroleum Engineers, 1993. [8] B.L. Litvak, Simulation and characterization of naturally fractured reservoirs, in: Proceedings of the 1985 Reservoir Characterization Technical Conference, Dallas, 1985. [9] James R. Gilman, H. Kazemi, Improved calculations for viscous and gravity displacement in matrix blocks in dual-porosity simulators (includes associated papers 17851, 17921, 18017, 18018, 18939, 19038, 19361 and 20174), J. Petrol. Technol. 40 (01) (1988) 60–70. [10] Peter Dietrich et al., Flow and Transport in Fractured Porous Media, Springer Science & Business Media, 2005. [11] J.B. Passioura, Hydrodynamic dispersion in aggregated media: 1. Theory, Soil Sci. 111 (6) (1971) 339–344. [12] J.B. Passioura, D.A. Rose, Hydrodynamic dispersion in aggregated media: 2. Effects of velocity and aggregate size, Soil Sci. 111 (6) (1971) 345–351. [13] P.S.C. Rao et al., Solute transport in aggregated porous media: theoretical and experimental evaluation, Soil Sci. Soc. Am. J. 44 (6) (1980) 1139–1146. [14] Paul V. Roberts, The influence of mass transfer on solute transport in column experiments with an aggregated soil, J. Contam. Hydrol. 14 (1987) 375–393. [15] Qinhong Hu, Mark L. Brusseau, Effect of solute size on transport in structured porous media, Water Resour. Res. 31 (7) (1995) 1637–1646.

[16] Jin-Ping Gwo, Rachel O’Brien, Philip M. Jardine, Mass transfer in structured porous media: embedding mesoscale structure and microscale hydrodynamics in a two-region model, J. Hydrol. 208 (3) (1998) 204–222. [17] Piotr Małoszewski, Andrzej Zuber, On the theory of tracer experiments in fissured rocks with a porous matrix, J. Hydrol. 79 (3) (1985) 333–358. [18] Piotr Maloszewski, Andrzej Zuber, Mathematical modeling of tracer behavior in short-term experiments in fissured rocks, Water Resour. Res. 26 (7) (1990) 1517–1528. [19] Piotr Maloszewski, Andrzej Zuber, Tracer experiments in fractured rocks: matrix diffusion and the validity of models, Water Resour. Res. 29 (8) (1993) 2723–2735. [20] H. Abelin, et al., 3-D migration experiment-Report 3, part I, in: Performed experiments, results and evaluation. SKB Stripa Project Technical Report, 1987, pp. 87–21. [21] P.M. Jardine et al., Quantifying diffusive mass transfer in fractured shale bedrock, Water Resour. Res. 35 (7) (1999) 2015–2030. [22] Ivars Neretnieks, Diffusion in the rock matrix: an important factor in radionuclide retardation, J. Geophys. Res. 85 (B8) (1980) 4379–4397. [23] G.E. Grisak, J.-F. Pickens, Solute transport through fractured media: 1. the effect of matrix diffusion, Water Resour. Res. 16 (4) (1980) 719–730. [24] D.H. Tang, E.O. Frind, Edward Allan Sudicky, Contaminant transport in fractured porous media: analytical solution for a single fracture, Water Resour. Res. 17 (3) (1981) 555–564. [25] E.A. Sudicky, E.O. Frind, Contaminant transport in fractured porous media: analytical solutions for a system of parallel fractures, Water Resour. Res. 18 (6) (1982) 1634–1642. [26] Stan Feenstra, Matrix diffusion effects on contaminant migration from an injection well in fractured sandstone, Ground Water 22 (3) (1984) 307–316. [27] Ching-Hor Lee, Shi-Ping Teng, An analytical model for radionuclide transport in a single fracture: considering nonequilibrium matrix sorption, Nucl. Technol. 101 (1) (1993) 67–78. [28] Chu-Tien Chen, Shih-Hai Li, Radionuclide transport in fractured porous media—analytical solutions for a system of parallel fractures with a constant inlet flux, Waste Manage. 17 (1) (1997) 53–64. [29] Jesús Carrera, On matrix diffusion: formulations, solution methods and qualitative effects, Hydrogeol. J. 6 (1) (1998) 178–190. [30] Sean W. Fleming, Roy Haggerty, Modeling solute diffusion in the presence of pore-scale heterogeneity: method development and an application to the Culebra dolomite Member of the Rustler Formation, New Mexico, and USA, J. Contam. Hydrol. 48 (3) (2001) 253–276. [31] Matthew W. Becker, Allen M. Shapiro, Interpreting tracer breakthrough tailing from different forced-gradient tracer experiment configurations in fractured bedrock, Water Resour. Res. 39 (1) (2003). [32] H.H. Liu, G.S. Bodvarsson, G. Zhang, Scale dependency of the effective matrix diffusion coefficient, Vadose Zone J. 3 (1) (2004) 312–315. [33] M. Chu, P.K. Kitanidis, P.L. McCarty, Dependence of lumped mass transfer coefficient on scale and reactions kinetics for biologically enhanced NAPL dissolution, Adv. Water Resour. 30 (6) (2007) 1618–1629. [34] Quanlin Zhou, Field-scale effective matrix diffusion coefficient for fractured rock: results from literature survey, J. Contam. Hydrol. 93 (1) (2007) 161–187. [35] S.A. Mathias et al., A parameter sensitivity analysis of two Chalk tracer tests, Q. J. Eng. Geol. Hydrogeol. 42 (2) (2009) 237–244. [36] Abdullah Cihan, John S. Tyner, 2-D radial analytical solutions for solute transport in a dual-porosity medium, Water Resour. Res. 47 (4) (2011). [37] Geoffrey Taylor, Dispersion of soluble matter in solvent flowing slowly through a tube, in: Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, vol. 219, The Royal Society, 1953. No. 1137. [38] Jacob Bear, Dynamics of Fluids in Porous Media, Courier Corporation, 2013. [39] R. Allan Freeze, John A.R.A. Cherry, J.A. Cherry, Groundwater, No. 556.3 FRE, 1979. [40] Amin Sharifi Haddad, Hassanzadeh Hassan, Jalal Abedi, Advective–diffusive mass transfer in fractured porous media with variable rock matrix block size, J. Contam. Hydrol. 133 (2012) 94–107. [41] J. Skopp, A.W. Warrick, A two-phase model for the miscible displacement of reactive solutes in soils, Soil Sci. Soc. Am. J. 38 (4) (1974) 545–550. [42] Fritz Oberhettinger, Larry Badii, Tables of Laplace Transforms, Springer Science & Business Media, 2012. [43] G.E. Grisak, J.-F. Pickens, Solute transport through fractured media: 1. The effect of matrix diffusion, Water Resour. Res. 16 (4) (1980) 719–730. [44] G.E. Grisak, J.F. Pickens, J.A. Cherry, Solute transport through fractured media: 2. Column study of fractured till, Water Resour. Res. 16 (4) (1980) 731–739.