A Simple Statistical Model for Transmissivity ... - Springer Link

Transp Porous Med (2015) 108:649–657 DOI 10.1007/s11242-015-0493-x

A Simple Statistical Model for Transmissivity Characteristics Curve for Fluid Flow Through Rough-Walled Fractures Cheng Yu1,2

Received: 9 December 2014 / Accepted: 26 March 2015 / Published online: 2 April 2015 © Springer Science+Business Media Dordrecht 2015

Abstract The fracture transmissivity characteristics curve (Witherspoon et al. in Water Resour. Res. 16(6):1016–1024, 1980) is found to deviate from cubic law as aperture decreases and still have residual transmissivity when aperture is very small. The existing models can partly explain the deviation from cubic law (e.g., Sisavath et al. in PAGEOPH 160:1009– 1022, 2003), or the residual transmissivity due to irreducible flow (e.g., Nolte et al. in PAGEOPH 131(1/2):111–138, 1989). In order to predict the transmissivity curve with both the above characteristics, in this study, a simple statistical model is employed with the following assumptions: (1) fracture boundaries are assumed parallel flat at global scale, but with normally distributed aperture variations at local scale (like frosted glasses); and (2) in this case, the flow field is assumed regular with straight head-contours and flow-lines. Then the equivalent transmissivity can be approximated as a series of parallel-connected local transmissivities. The transmissivity curve can be fitted very well with both the above characteristics. It is suggested that the reason for the deviation from cubic low is possibly due to the variations of local apertures which induce redistribution of hydraulic gradients, and the residual foot is because of residual open apertures or micro-fractures in the fracture surfaces. Keywords

Fracture transmissivity · Statistical method · Cubic law

1 Introduction Fracture flow through cracked rocks is a major concern in geosciences (e.g., Garven 1995), environmental engineering (e.g., Long and Ewing 2004; Myers 2012), unconventional petro-

B

Cheng Yu [email protected]

1

Key Laboratory of Hydraulic and Waterway Engineering of the Ministry of Education, Chongqing Jiaotong University, Chongqing, China

2

Center for Water Research, Peking University, Beijing, China

123

650

C. Yu

leum engineering (e.g., Curtis 2002), and geologic hazards sciences (e.g., Wintsch et al. 1995). However, so far the fluid flow and transport processes through single fractures are still not yet fully understood. Previously, single fractures were usually considered with uniform average apertures and in this case the fracture transmissivity (T ) followed the “cubic law.” Since single fractures are rough-walled with stagnant surface contacts (e.g., Tsang 1984; Brown 1987; Tsang and Tsang 1987; Mourzenko et al. 1995), the hydraulic gradient always induces 3-D, non-uniform, and tortuous flow field. The impact of fracture roughness has been observed in various experiments, for example the transmissivity characteristics curve by Witherspoon et al. (1980) is still calling attentions. This curve suggests that as nominal aperture (H ) decreases, the relationship T ∼ H experiences three regimes. Regime III: when H is large, the impact of roughness is negligible, and T follows cubic law; Regime II: when H is comparable with aperture deviation, T tends to decrease faster, and the transmissivity curve tends to deviate from cubic law. The log–log exponent can even be as large as 8–10 (Pyrak-Nolte et al. 1987); Regime I: when H is very small, T becomes very small but still tends to have residual foot as H decreases. Researchers have supposed that the deviation from cubic law (regime II) is due to: (1) stagnant zones where fracture surfaces contact (Witherspoon et al. 1980), which reduce the effective flowing area; (2) tortuosity of flow-lines induced by fracture roughness and contacts, which reduces the effective hydraulic gradient (Nolte et al. 1989); (3) the impact of tortuous boundaries which affects the flow pattern along flow-lines (Sisavath et al. 2003). Nolte et al. (1989) employed fractal geometry to describe the impact of contact areas, providing an idea that fracture transmissivity can be examined statistically, regardless of flow field details. Sisavath et al. (2003) predicted an accurate transmissivity for 2-D Navier– Stokes (N–S) flow between symmetric sinusoidal boundaries. However, real fractures are always more complex with asymmetric boundaries (e.g., Velde et al. 1990); therefore, their model can hardly represent real 3-D fracture flows. In this study in order to predict the transmissivity curve with both the above characteristics, a simple statistical model is employed with the following assumptions: (1) fracture boundaries are assumed parallel flat at global scale, but with normally distributed aperture variations at local scale (like frosted glasses); and (2) in this case, the flow field is assumed regular with straight head-contours and flow-lines. Then the equivalent transmissivity can be approximated as a series of parallel-connected local transmissivities. The transmissivity curve can be fitted very well with both the above characteristics. It is suggested that the reason for the deviation from cubic low is possibly due to the variations of local apertures which induce redistribution of hydraulic gradients, and the residual foot is because of residual open apertures or micro-fractures in the fracture surfaces.

2 About the Model by Sisavath et al. (2003) Sisavath et al. (2003) gave an analytical solution to 2-D N–S flow between symmetric sinusoidal curves: 2π x 1 y(x) = ± H0 + δ sin ; x ∈ [0, λ] (1) 2 λ where λ and δ are the wavelength and amplitude of sinusoidal curves, H0 is the average aperture (or nominal aperture), and the local fracture apertures are found as: 2π x (2) h(x) = H0 + δ sin λ

123

A Simple Statistical Model for Transmissivity Characteristics…

651

Fig. 1 The plot T ∼ H in terms of Eq. (5) when δ = 0.4 and λ = 10. When α = δ, T ∼ H looks very close to cubic law

Then the fracture transmissivity is found by Sisavath et al. (2003): −1 g H03 (1 − δˆ2 )2.5 36 2 δˆ2 1 − δˆ2 T (H0 ) = 1+ π 12ν 1 + 0.5δˆ2 15 λˆ 2 1 + 0.5δˆ2

(3)

2λ where δˆ = Hδ0 , λˆ = H , and υ is dynamic viscosity. The T decreases faster and deviates 0 from cubic law as H0 decreases, as the log–log exponent of H0 ∼ T decreases from around 1/3 to 0. This deviation is due to the impact of boundary roughness, but if supposing this impact can be some correction α to the measure of H0 :

H = H0 − α where H is the corrected nominal aperture, then the corrected Tc is: 2.5 −1 3 1 − δˆ 2 (H + α) g 36 2 δˆ2 1 − δˆ2 Tc = 1+ π 12ν 15 λˆ 2 1 + 0.5δˆ2 1 + 0.5δˆ2

(4)

(5)

As α increases from 0, the trend of Tc ∼ H varies as Fig. 1. Equation (5) is the same with (3) as solution for the N–S fracture flow, giving other descriptions of transmissivity with respect to other nominal apertures. It is worth discussing when α = δ, the Tc ∼ H still approximately follows power law with two regimes and an inflection at H ∝ δ (Fig. 1). The regime II is cubic law, but regime I shows a slight deviation. This is because when H δ, if supposing all hydraulic drawdown is applied between equivalent throat neighborhoods, it turns into fracture flows at smaller scale only around the throats (Fig. 2). As H decreases, the equivalent neighborhoods shrink, and thus Tc shows the slight deviation.

123

652

C. Yu

Fig. 2 Demonstration of the throat neighborhood where most of the hydraulic drawdowns are concentrated when α = δ

The above discussion implies the inaccuracy of local transmissivities of “Local Cubic Law” (LCL). In fact LCL considers local boundaries parallel flat and applies cubic law to local transmissivities: g T = (6) H3 12ν 0 However as key issues, (1) LCL is a solution of Reynolds equation rather than N–S equations; and (2) in fact local boundaries are definitely not flat but still rough-walled even at local scales (e.g., Brush and Thomson 2003; Koyama et al. 2008, 2009; Li et al. 2014). However, it is implied that the impact of local roughness to N–S flow transmissivity can be some correction to the nominal aperture: g Tc = (H0 − β)3 (7) 12ν where β is similar to α. It is somewhat confusing but Tc is for real N–S flow while T is only for LCL. Moreover, it is implied that the β should be around the variance of H0 within local neighborhood. Usually, the variance is much smaller at most locations (β H0 ), which suggests the rationality of LCL: g H3 Tc ≈ T = (8) 12ν 0 As a brief summary, the purpose of this section is to provide a general idea how fracture roughness affects the fracture transmissivity and to suggest the rationality of LCL employed in the following section.

3 Statistical Model for 3-D Fracture Transmissivity Distribution of local apertures is important to the overall transmissivity, but is difficult to be figured out because of randomness and complexity of fracture geometries. However, it is shown that apertures almost always tend to follow normal distributions (e.g., Koyama et al. 2009). Therefore in this study, the fracture geometry is assumed that: at local scale, the variation of local apertures follows a normal distribution, and the flow field is non-uniform; however, at global scale, the flow field is assumed statistically uniformly regular, with straight headcontours and flow-lines. In other words, fractures are like flat-frosted glasses with rough

123


653

Fig. 3 Conceptual model for the statistical model. Although at local scale the flow field is tortuous, the flow field is assumed to be regular at global scale with straight head-contours and flow-lines

Fig. 4 The LCL model is discretized into M(x) × N(y) cells, and the fracture flow is considered through a series of equivalent cells

surfaces. In this case, the roughness and stagnant contact areas are both assumed to have impacts only locally. It is noteworthy although fracture has been widely agreed fractal, the above assumption does not follow any fractal characteristics. Because theoretically lateral scales of fractal profiles can extend even to infinity and it’s not very true at least at global scale for samples prepared under controlled conditions (e.g., Brazilian fracturing). Another concern is that local aperture is in fact the difference between upper and lower boundaries, so apertures don’t have to be fractal even if boundaries are. LCL is then employed to calculate the overall transmissivity. As the flow field is assumed regular, typically the conceptual model can be set as Fig. 3. Flow-lines flow along x-axis and head-contours are along y-axis. If discretizing the flow field into M(x) × N(y) lattice square cells, then the flow can be considered 1D through a series of M equivalent lattice cells T (i) (Fig. 4). Each T (i) is equivalent to a parallel-connected column of N lattice cells, which has transmissivity T j ( j = 1, 2, . . . N ) . Then Ti , (i = 1, 2, . . . M) turns out as:

T (i) =

N N 3 1 g 1

· h j (i) T j (i) = N 12ν N j=1

(9)

j=1

123

654

C. Yu

where h j (i) are local apertures of the N cells in column i. If cells are subdivided to infinite small (N → ∞), then h j (i) becomes h(i) and T (i) turns into the expectation: g (10) T (i) = · E h 3 (i) 12ν Similarly, if cells are divided to finite small in x, the overall transmissivity is found as: T =

1 1 M→∞ −−−−→ 1 M 1 1 E M i=1 T (i) Tx

(11)

where Tx is the equivalent transmissivity at coordinate x. Equations (10) and (11) give a semiquantitative transmissivity. The “semi” indicates its inaccuracy, mainly due to: (1) the assumption about the statistical distribution of geometries; and (2) the inherent inaccuracy of LCL. During the past decades, many studies have introduced using “series” or “parallel” conductance bounds to estimate the effective transmissivity of 3-D heterogeneous fracture flow (e.g., Dagan 1979; Neuzil and Tracy 1981; de Marsily 1986; Silliman 1989; Dagan 1993; Zimmerman and Bodvarsson 1996); however, these discussions are more semiqualitative, and none of them has represented the transmissivity characteristic curve. More specifically, as aperture variation is assumed normally distributed, the h (a function of y) at location x can be noted as: h = h0 + (12) where h 0 is the expectation, is the aperture variation in y, which follows distribution f ( ): f ( ) = √

1 2πσ

e

−

2 2σ 2

(13)

where σ is the standard deviation. Then mathematically Tx is the expectation: Tx =

g 12ν

∞ (h 0 + )3 f ( )d

(14)

−∞

The h 0 , , and f are all functions of x. But in fact stagnant areas at closed apertures (h 0 + ≤ 0) do not contribute to Tx , thus Tx is revised: ∞ g Tx = (h 0 + )3 f ( )d = F(h 0 , σ ) (15) 12ν −h 0 − h 20 2 2 g 3 h0 2 (16) h 0 σ + 2σ 3 e 2σ 2 h 0 + 3σ h 0 1 + erf √ F(h 0 , σ ) = + 24ν π 2σ If supposing h 0 can be regarded as a random variable as well: h0 = H + ε

(17)

where H is the expectation (average aperture of the whole fracture), and ε is the variation of h 0 which follows distribution g(ε): g(ε) = √

123

1 2πτ

e

−

ε2 2τ 2

(18)


655

Fig. 5 The transmissivity characteristic curve fitted to the data from Witherspoon et al. (1980). Good fits are achieved in regime III and I, while errors mostly appear in regime II. The errors can be due to the assumptions of this model and the inherent inaccuracy of LCL as well

where τ is the standard deviation. If further assuming the deviation σ (x) is constant (σ (x) = σ , or called “average σ ”), the T can be estimated: ⎛ T =⎝

∞

−∞

⎞−1 ε2 1 1 − e 2τ 2 dε ⎠ √ F(H + ε, σ ) 2π τ

(19)

Hopefully, the plot T ∼ H can be carried out numerically. Different from Eq. (15), we do not exclude the ε where H +ε < 0, because in fact F is always positive at all H +ε ∈ (−∞, ∞), which suggests the fracture is permeable at all x locations even if there are widely distributed stagnant areas. Note that the variations of local apertures are described in two independent random variables and ε, which in fact quantitatively represent local heterogeneities of the planar flow field. But if supposing aperture variation is described only in a single random variable, then no ε will leave after averaging “in parallel”; or even if and ε are independent but the averaging “in parallel” is over finite field width w, the ε will depend on w and σ . In both case the heterogeneity along x is not independent, and the final transmissivity (19) cannot be achieved. Therefore in fact, and ε are introduced more like purpose-driven as equivalent coefficients in order to fully represent the local uniformities of flow field which have been stated in the assumptions but have not been considered in the “series” or “parallel” connected local cells.

4 About the Experiment by Witherspoon et al. (1980) By adjusting σ, τ and taking H as the argument, good fits are achieved to the transmissivity curve (Fig. 5). It makes sense for rock samples the aperture deviations are of magnitude of µm. Quantitatively, τ mainly determines the deviation of T ∼ H from cubic law (regime II). As τ increases the plot deviates more from cubic law; while the σ mainly enhances the residual foot (regime I). As σ increases regime I moves rightward.

123

656

C. Yu

It was argued the residual foot was due to irreducible flow, or preferential channels (e.g., Pyrak-Nolte et al. 1988). In this study, the residual foot seems due to the pattern of F(H + ε, σ ). As H decreases, the slope of F decreases rapidly and F tends to approach 0. Thus the T (Eq. 19) tends to approach the foot as well. In fact, if supposing aperture variation follows other distributions such as a uniform distribution here, the Tx should be revised as: ⎧ ⎪ |h 0 | ≥ σ ⎨0 σ T (x) = F(h 0 , σ ) = g (h 0 + )3 ⎪ d |h 0 | < σ ⎩ 12ν −h 0 2σ ⎧ ⎨0 |h 0 | ≥ σ σ3 2 2σ 3 4 (20) = g hσ 3h h h ⎩ |h 0 | < σ + + + + 12ν 8 2 4 2 8σ where the “over line” indicates parameters for the uniform distribution. Then the overall transmissivity turns out as: ⎞−1 ⎛ ∞ 1 1 dε ⎠ T =⎝ (21) F(H + ε, σ ) 2τ −∞

The plot T tends to approach 0 as H decreases without any foot (Fig. 5). Therefore, the foot of Eq. (19) must be due to its normally distributed aperture variation, which always remains open apertures. Moreover, other paths might contribute to the irreducible flow as well, such as micro-fractures inside fracture surfaces under significant normal stresses (e.g., Dixon et al. 2009), which are usually of magnitude of µm, too.

5 Conclusions and Discussion In this study, in order to predict the fracture transmissivity curve (Witherspoon et al. 1980) with both the deviation and residual foot characteristics, a statistical model is employed with the following assumptions: (1) fracture boundaries are assumed parallel flat at global scale, but with normally distributed aperture variations at local scale (like frosted glasses) and (2) in this case, the flow field is assumed regular with straight head-contours and flow-lines. Then the equivalent transmissivity can be approximated as a series of parallel-connected local transmissivities. Although previously many studies have introduced similar approaches using “series” or “parallel” conductance bounds to estimate effective transmissivity, none of them have quantitatively represented the transmissivity characteristic curve. In this study, this curve is fitted well with both the deviation and residual foot. It is suggested that the reason for deviation from cubic low is possibly due to aperture variations which induce redistribution of hydraulic gradients, and the residual foot is due to residual open apertures or micro-fractures in the fracture surfaces. The limitation of this model is also about its assumption that: the flow field is regular at global scale. However, hopefully it makes sense for most of laboratory tests and the random variables in this model can to some extent represent the flow field heterogeneities as corrections to retain local flow uniformities. Acknowledgments This study is financially supported by National Natural Science Foundation of China (NSFC) (No. 51409028).

123


657

References Brown, S.R.: Fluid flow through rock joints: the effect of surface roughness. J. Geophys. Res. 92(B2), 1337– 1347 (1987) Brush, D.J., Thomson, N.R.: Fluid flow in synthetic rough-walled fractures: Navier–Stokes, Stokes, and local cubic law simulations. Water Resour. Res. 39(4), 1085 (2003) Curtis, J.B.: Fractured shale-gas systems. AAPG Bull. 86(11), 1921–1938 (2002) Dagan, G.: Models of groundwater flow in statistically homogeneous porous formations. Water Resour. Res. 15, 47–63 (1979) Dagan, G.: Higher-order correction of effective permeability of heterogeneous isotropic formations of lognormal conductivity distribution. Transp. Porous Media 12, 279–290 (1993) de Marsily, G.: Quantitative Hydrogeology. Academic Press, San Diego (1986) Dixon, N.A., Durham, W.B., Suzuki, A.M., Mei, S.: Measurements of grain size sensitivity in olivine deformed at high pressure in the deformation DIA. Resented at AGU 2009 Fall Meeting, San Francisco. MR41A1849 Poster (2009) Garven, G.: Continental-scale groundwater flow and geologic processes. Ann. Rev. Earth Planet. Sci. 23, 89–118 (1995) Koyama, T., Li, B., Jiang, Y., Jing, L.: Numerical modelling of fluid flow tests in a rock fracture with a special algorithm for contact areas. Comput. Geotech. 36, 291–303 (2009) Koyama, T., Neretnieks, I., Jing, L.: A numerical study on differences in using Navier–Stokes and Reynolds equations for modeling the fluid flow and particle transport in single rock fractures with shear. Int J. Rock Mech. Mining Sci. 45, 1082–1101 (2008) Li, B., Wong, R.C.K., Milnes, T.: Anisotropy in capillary invasion and flow fluid through induced sandstone and shale fractures. Int. J. Rock Mech. Mining Sci. 65, 129–140 (2014) Long, J.C.S., Ewing, R.C.: Yucca mountain: earth-science issues at a geologic repository for high-level nuclear waste. Ann. Rev. Earth Planet. Sci. 32, 363–401 (2004) Mourzenko, V.V., Thovert, J.F., Adler, P.M.: Permeability of a single fracture: validity of the Reynolds equation. J. Phys. II Paris 5(3), 465–482 (1995) Myers, T.: Potential contaminant pathways from hydraulically fractured shale to aquifers. Groundwater 50(6), 872–882 (2012) Neuzil, C.E., Tracy, J.V.: Flow through fractures. Water Resour. Res. 17, 191–199 (1981) Nolte, D.D., Pyrak-Nolte, L.J., Cook, N.G.W.: The fractal geometry of flow paths in natural fractures in rock and the approach to percolation. PAGEOPH 131(1/2), 111–138 (1989) Pyrak-Nolte, L.J., Cook, N.G.W., Nolte, D.D.: Fluid percolation through single fractures. Geophys. Res. Lett. 15(11), 1247–1250 (1988) Pyrak-Nolte, L.J., Myer, L.R., Cook, G.W., Witherspoon, P.A.: Hydraulic and mechanical properties of natural fractures in low permeability rock. In: Proceedings of the Sixth International Congress on Rock Mechanics, Montreal, Canada, August 1987, G. HergetS. Vongpaisal225-231, Pubs. A.A. Balkema, Rotterdam (1987) Silliman, S.E.: An interpretation of the difference between aperture estimates derived from hydraulic and tracer tests in a single fracture. Water Resour. Res. 25, 2275–2283 (1989) Sisavath, S., Al-Yaaruby, A., Pain, C.C., Zimmerman, W.R.: A simple model for deviations from the cubic law for a fracture undergoing dilation or closure. PAGEOPH 160, 1009–1022 (2003) Tsang, Y.W.: The effect of tortuosity on fluid flow through a single fracture. Water Resour. Res. 20(9), 1209– 1215 (1984) Tsang, Y.W., Tsang, C.F.: Channel model of flow through fractured media. Water Resour. Res. 23(3), 467–479 (1987) Velde, B., Dubois, J., Touchard, G., Badri, A.: Fractal analysis of fractures in rocks: the Cantor’s Dust method. Tectonophysics 179(3–4), 345–352 (1990) Wintsch, P., Christoffersen, R., Kronenberg, A.K.: Fluid-rock reaction weakening of fault zones. J. Geophys. Res. Solid Earth (1978–2012) 100(B7), 13021–13032 (1995) Witherspoon, P.A., Wang, J.S.Y., Iwai, K., Gale, J.E.: Validity of cubic law for fluid flow in a deformable rock fracture. Water Resour. Res. 16(6), 1016–1024 (1980) Zimmerman, R.W., Bodvarsson, G.S.: Hydraulic conductivity of rock fractures. Transp. Porous Media 23, 1–30 (1996)

123