Elastic wave attenuation and dispersion induced by mesoscopic flow in double-porosity rocks

Elastic wave attenuation and dispersion induced by mesoscopic flow in double-porosity rocks

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111 Contents lists available at ScienceDirect International Journal of Rock ...

796KB Sizes 0 Downloads 33 Views

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111

Contents lists available at ScienceDirect

International Journal of Rock Mechanics & Mining Sciences journal homepage: www.elsevier.com/locate/ijrmms

Elastic wave attenuation and dispersion induced by mesoscopic flow in double-porosity rocks

MARK



Pei Zhenga, , Boyang Dingb, Xiuting Suna a b

School of Mechanical Engineering, University of Shanghai for Science and Technology, Shanghai 200093, China Department of Civil Engineering, Zhejiang University of City College, Hangzhou 310015, China

A R T I C L E I N F O

A BS T RAC T

Keywords: Double-porosity medium Mesoscopic flow Wave attenuation

The double-porosity, dual-permeability theory is employed to predict the wave attenuation and phase velocity dispersion induced by wave-induced mesocopic fluid flow. Instead of using an up-scaled, single-porosity approximation scheme proposed by previous researchers, we develop an analytical method to exactly solve the wave equations for double-porosity materials. We first propose a new form of wave equations formulated in terms of displacements. This new form of wave equations enables us to decouple the field equations into two second-order symmetric dynamic systems, namely, the P-system for compressional waves and the S-system for shear wave. We then implement Newton iteration for solving the cubic dispersion equation for compressional waves. Finally, to understand the loss mechanism caused by mesoscopic flow, we compare the attenuation curves of the first (P1-wave), the second (P2-wave), and the third compressional waves (P3-wave), as well as the shear wave (S-wave), with the mesoscopic flow present to those with the mesoscopic flow absent. Furthermore, the effects of matrix porosity, pore-fluid viscosity, and values of fluid transport coefficient on wave attenuation are also investigated in numerical examples.

1. Introduction Most earth materials such as rocks and sediments are generally heterogeneous and often fractured or cracked. In such materials, the pore and crack/fracture space may be filled with water, oil, or gas. When elastic waves propagate through such porous materials, fluidrelated wave attenuation and dispersion occur. It is commonly accepted that the relative fluid-solid movement before the establishment of pore pressure equilibrium, known as the wave-induced fluid flow (WIFF), is the main cause of wave energy losses. Biot's theory of poroelasticity1,2 may be the most commonly used theoretical model to estimate the frequency-dependent attenuation and dispersion caused by WIFF. Nevertheless, during the past decades, a number of studies3–5 have shown that the level of attenuation predicted by using the Biot's theory is drastically underestimated in the seismic band and can only be significant at frequency above 10 kHz, well outside the seismic exploration band. The principal reason for this is that in Biot's theory the porous material is assumed to be spatially homogeneous (macroscopic homogeneous) and fluid-saturated, and, as a result, only the macroscopic flow, resulting from pore pressure gradients established between the peaks and troughs of the wave, occurs, which does not produce enough loss. This mechanism is broadly known as the Biot loss and, roughly speaking, the intrinsic



loss, quantified using the inverse quality factor Q−1, attains Biot-loss maximum when the wavelength-scale pressure gradients just equilibrate in a wave cycle. A possible alternative to Biot's loss is the socalled “mesoscopic-loss mechanism”. Many geological materials such as reservoir rocks are generally heterogeneous and often fractured or cracked. When a compressional wave stress a material containing such mesoscopic heterogeneities (i.e., the heterogeneities exist on length scales smaller than the wavelength but greater than the pore scale), the pore pressure developed on the mesoscale drives the fluid flow between the more compliant parts of the material, e.g., cracks/fractures, and the stiffer portions, e.g., the background pores. Such wave-induced mesoscopic flow is increasingly believed to be a dominant mechanism of fluid-related attenuation in the seismic band. This mechanism is the so-called “mesoscopic-loss mechanism”. Several theoretical models for wave attenuation and dispersion in porous rocks containing mesoscopic heterogeneities have been proposed during the past decades3,6–17 (a comprehensive review of different models is provided by Müller et al.18). Among them, the dual-porosity, dual-permeability model may be the simplest and ideal approach to modeling fractures in porous media based on Biot's theory of poroelasticity.9,12,13 In such models, the mesoscopic cracks/fractures are treated as a fracture network embedded in porous matrix blocks, and it is generally assumed that the fracture permeability is

Corresponding author.

http://dx.doi.org/10.1016/j.ijrmms.2016.11.018 Received 12 April 2016; Received in revised form 30 September 2016; Accepted 15 November 2016 Available online 24 November 2016 1365-1609/ © 2016 Elsevier Ltd. All rights reserved.

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111

P. Zheng et al.

⎛ e ⎞ ⎛ a11 − a12 − a13 ⎞ ⎛ σ ⎞ ⎜ ζ1 ⎟ = ⎜− a12 a22 a23 ⎟⎟ ⎜ pf 1 ⎟ , ⎜ ⎜ ⎟ ⎜− a a33 ⎠ ⎝ pf 2 ⎟⎠ ⎝ ζ2 ⎠ ⎝ 13 a23

much greater than the matrix permeability, while the volume fraction of the fracture phase is much smaller than that of the matrix phase. Early double-porosity models were initially proposed for large-scale transport problems in fissured rocks in hydrology.19 More sophisticated models9,12,13,17,20–29 were developed later to account for the coupling between rock deformation and flow. In particular, the doubleporosity, dual-permeability models developed by Berryman and Wang9 and Pride and Berryman12,13 have aroused much attenuation in recent years. Rigorously following Biot's approach, Berryman and Wang9 make a generalization of Biot's single-porosity poroelasticity to double-porosity, dual-permeability theory to incorporate fractures. In their model, the fractures/cracks are assumed to be throughgoing joints, namely the fracture volume is all void space. The model of Berryman and Wang predicts three distinct types of compressional waves exist in a doubleporosity medium. Nevertheless, the fluid flow between the matrix and fracture phases, that is the so-called “mesoscopic flow”, are neglected which effect is believed to be a major cause of fluid-related attenuation in the seismic band. As a result, only Biot loss is present in their model. Thus, more recent work of Pride and Berryman12,13 using the volume averaging technique to derive the governing equations for double-porosity media has incorporated mesoscopic flow mechanism by allowing fluid transport between two poroelastic phases. In their model, the mesoscopic heterogeneities, which are not confined to fractures/cracks, are envisioned to be a porous continuum (e.g., a less-well consolidated sandstone), embedded within a stiffer porous host (e.g., a consolidated sandstone). But instead of solving the governing equations for the double-porosity model directly, they propose an ingenious scheme to reduce their model to an effective single-porosity Biot theory having complex frequency-dependent coefficients. Using the opproposed scheme, they demonstrate that the predicted attenuation level of the fast compressible wave (P1-wave) is consistent with that measured by Sams et al.30 Furthermore, they extend the double-porosity framework to incorporate two other theoretical models, namely the squirt-flow model and patchy-saturation model.5 Although we can apply the proposed approximation scheme to the study of P1-wave attenuation and dispersion, it is impossible to employ their method to investigate wave phenomena associated with the second slow wave (P3-wave). Thus a need remains for exactly solving the wave equations for a heterogeneous doubleporosity model, which incorporates the mesoscopic-loss mechanism. The purpose of the present paper is first to provide an analytical method to exactly solve the field equations for double-porosity materials, and then investigate the effect of mesoscopic flow on the attenuation and phase velocity dispersion of compressional and shear waves. To this end, using the solid displacement u and the relative fluid-solid displacements w1, w2 as independent variables, we first reformulate the governing equations of Pride and Berryman in terms of displacements. This new form of wave equations enables us to decouple the field equations into two second-order symmetric dynamic systems, each of them consisting of three coupled wave equations. Next, we substitute plane wave solutions into the obtained symmetric dynamic systems to find the dispersion relations for compressional and shear waves, and subsequently, we implement Newton iteration for solving the cubic characteristic equation for compressional waves. Finally, to understand the mesoscopic loss mechanism in double-porosity, dualpermeability models, we give examples of P1-, P2-, P3-, and S-waves attenuation and phase velocity dispersion.

(1)

where e ≡ ∂V / V is the volumetric strain, σ ≡ (σxx + σyy + σzz )/3 is the mean normal stress, pf 1 and pf 2 are the pore pressures in the phases 1 and 2, and the other two variables ζ1 and ζ2 are respectively the increments of fluid content in the two phases. The six independent coefficients aij appearing in the matrix on the right hand side are material constants. If phase 2 is assumed to be a fracture phase with through-going joints, the constants aij are given by31–33

a11 =

1 , K

(2)

a12 = 1/ Ks − v1/ K1,

(3)

a13 = v1/ K1 − 1/ K ,

(4)

a22 = v1 α1/ B1 K1,

(5)

a23 = v1/(Ks )1 − 1/ Ks,

(6)

a33 = v2 / Kf + 1/ K − v1/ K1,

(7)

where K and K1 are the drained bulk moduli of the whole and the matrix (phase 1), Ks and (Ks )1 are the unjacketed bulk moduli of the whole and the matrix, α1 = 1 − K1/(Ks )1 is the corresponding Biot-Willis coefficient, Kf is the pore fluid bulk modulus, B1 is the Skempton's porepressure buildup coefficient for the matrix, and v2 = 1 − v1 is the total volume fraction of the fractures in the whole. If, on the other hand, phase 2 is assumed to be a porous continuum, aij are determined by5,12,32

a11 =

1 , K

(8)

a12 = −

v1 Q1 α1, K1

(9)

a13 = −

v2 Q2 α2, K2

(10)

a22 =

v1 α1 ⎛ 1 α (1 − Q1) ⎞ − 1 ⎟, ⎜ K1 ⎝ B1 1 − K1/ K2 ⎠

(11)

a23 =

α1 α2 K1 K2 ⎛ v1 v 1⎞ + 2 − ⎟, ⎜ (K2 − K1)2 ⎝ K1 K2 K⎠

(12)

a33 =

v2 α2 ⎛ 1 α (1 − Q2 ) ⎞ − 2 ⎟, ⎜ K2 ⎝ B2 1 − K2 / K1 ⎠

(13)

where the Qi are auxiliary constants given by

v1 Q1 =

1 − K2 / K , 1 − K2 / K1

v2 Q2 =

1 − K1/ K . 1 − K1/ K2

(14)

In Eq. (1), the variables e, ζ1, ζ2 are expressed in terms of the three independent variables σ , pf 1, pf 2 , and thus, this form of constitutive equations is termed the pure compliance formulation.33 On the other hand, the variables e , ζ1, ζ2 can be also chosen to be independent variables and, by solving e , ζ1, ζ2 in terms of σ , pf 1, pf 2 , the inverse relation of (1) (termed the pure stiffness formulation) can be found to be33

2. Field equations

⎛ σ ⎞ ⎛ H + 2/3μ − C1 − C2 ⎞ ⎛ e ⎞ ⎟ ⎜ pf 1 ⎟ = ⎜ − C M1 N ⎟ ⎜⎜ ζ1 ⎟⎟ , 1 ⎜p ⎟ ⎜ ⎝ f 2 ⎠ ⎝ − C2 N M2 ⎠ ⎝ ζ2 ⎠

2.1. Constitutive equations The linear constitutive relations among strain, fluid content, and pressure for isotropic materials with double-porosity can be written in the form9,12,13,31–33

(15)

where the constant μ is the shear modulus of the drained media and other constants H , C1, C2 , M1, M2 , N can be written in terms of aij : 105

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111

P. Zheng et al.

H + 2/3μ =

2 a22 a33 − a 23 , Δ

a13 a23 − a12 a33 , Δ

(17)

a a − a13 a22 C2 = 12 23 , Δ

(18)

2 a11 a33 − a13 M1 = , Δ

(19)

2 a11 a22 − a12 , Δ

(20)

C1 =

M2 =

a a − a11 a23 N = 12 13 , Δ

frequency-dependent. This fluid transport coefficient κ is commonly believed to be a function of the matrix permeability, fluid viscosity, fracture spacing, and surface area flow geometry.34 In case of treating phase 2 as a porous continuum, Pride and Berryman5,12,13 suggest that it can be defined as

(16)

κ (ω ) = κ 0 1 − i

(31)

where κ0 = κ (ω = 0) is the static transport coefficient and ωc is the crossover frequency, or the relaxation frequency, which is determined by the times it takes a pressure wave to diffuse across a characteristic length scale (over which the fluid-pressure gradient still exists in phase 1). The explicit formula for κ (ω) and ωc can be found in Pride and Berryman.5,12,13 With the aid of the constitutive Eqs. (1)–(15) and the fluid mass conservation Eqs. (28) and (29), the dynamic equilibrium Eqs. (23)– (25) can be reformulated in terms of displacements u , w1 and w2 as

(21)

where

a11 − a12 − a13 a23 . Δ = det A = − a12 a22 − a13 a23 a33

ω , ωc

4

[(b11 + 3 μ) + b12 χ1 (C1 − χ2 ) + b13 χ1 (C2 − χ2 )]∇∇⋅u − μ∇ × ∇ × u

(22)

+ χ1 [b12 (M1 − χ3 ) + b13 (N − χ3 )]∇∇⋅w1 + χ1 [b12 (N − χ3 ) + b13 (M2 − χ3 )]∇ ∇⋅w2 + ω 2 (ρu + ρf w1 + ρf w2) = 0 ,

2.2. Dynamic equilibrium equations

(32)

Assuming an e−iωt time dependence, the dynamic equilibrium equations for double-porosity materials can be written in the form9,12,13:

∇⋅σ = −ω 2 (ρu + ρf w1 + ρf w2),

∇pf 1 −

ω 2 (ρf

u + m22 w1) = 0 ,

∇pf 2 − ω 2 (ρf u + m33 w2) = 0,

χ1 (C1 − χ2 )∇∇⋅u + χ1 (M1 − χ3 )∇∇⋅w1 + χ1 (N − χ3 )∇∇⋅w2 + ω 2 [ρf u + m22 w1] = 0 ,

χ1 (C2 − χ2 )∇∇⋅u + χ1 (N − χ3 )∇∇⋅w1 + χ1 (M2 − χ3 )∇∇⋅w2

(24)

+ ω 2 [ρf u + m33 w2] = 0 ,

m22

τ2 ρf ρ33 i η i η m33 = + = + , (v 2 ϕ 2 ) 2 ω κ2 v 2 ϕ2 ω k2

(34)

in which

(25)

χ1 =

where σ is the average stress tensor action on the bulk materials, u is the average displacement of the solid phase, w1 ≡ (1 − v2 ) ϕ1 (U1 − u) and w1 ≡ v2 ϕ2 (U2 − u) (with U1 and U2 being the average displacements of the fluid in phase 1 and 2) are the average flow of fluid in the two phases relative to the solid, ρf and ρ are the densities of the pore fluid and the bulk materials. In addition, m22 and m33 are defined as

τ1 ρf ρ22 i η i η = + = + , [(1 − v2 ) ϕ1]2 ω κ1 (1 − v2 ) ϕ1 ω k1

(33)

(23)

1 , 1 − κ (M1 + M2 − 2N )/ iω

(35)

χ2 = κ (C1 (M2 − N ) + C2 (M1 − N ))/ iω,

(36)

χ3 = κ (M1 M2 − N2 )/ iω,

(37)

are defined to show the influence of fluid transport between the two phases on the equilibrium equations. In other words, if such a mesoscopic flow is neglected, χ1 = 1, χ2 = χ3 = 0 . In Eqs. (32)−(34), bij are material constants which can be determined in terms of aij 33:

(26)

b11 =

(27)

where ρ22 = τ1 (1 − v2 ) ϕ1 ρf and ρ33 = τ2 v2 ϕ2 ρf are the induced masses resulting from inertial drag caused by the relative acceleration of the solid frame and the pore fluid,9τ1 and τ2 are structure factors of phase 1 and 2, ϕ1 and ϕ2 are the porosities of phase 1 and 2, η is the fluid viscosity, and k1 and k2 are the permeabilities of phase 1 and 2. It should be mentioned that, following Berryman and Wang,9 we have omitted the permeability coupling term k12 in Eqs. (24) and (25).

1 , a11

b12 = −

a12 , a11

b13 = −

a13 . a11

(38)

3. Dispersion equations and their solution 3.1. Wave potentials We now consider a Helmholtz resolution of each of the displacement fields in terms of

2.3. Fluid mass conservation equations

u = ∇φs + ∇ × ψs with

∇⋅ ψs = 0;

(39)

The fluid mass conservation equations capable of modeling fluid transfer between the two phases are given in the form5,12,13

w1 = ∇φf 1 + ∇ × ψf 1

−iωζ1 + ∇⋅q1 = iωζint ,

(28)

w2 = ∇φf 2 + ∇ × ψf 2

−iωζ2 + ∇⋅q2 = −iωζint ,

(29)

according to which, scalar potentials φs , φf 1 and φf 2 are associated with the dilatational motions in the double-porosity materials while vector potentials ψs , ψf 1 and ψf 2 are associated with the rotational motions in the media. Introducing (39)–(41) into (32)–(34), two uncoupled (dilationrotation) sets of equations given in matrix form are obtained

where q1 ≡ ẇ1 and q2 ≡ ẇ 2 are the filtration velocities of phase 1 and 2, −iωζint , the average rate at which fluid volume is being transferred from phase 1 into phase 2, is defined as

−iωζint = κ (pf 1 − pf 2 ).

(30)

[Kp ∇2 + ω 2 M] φ = 0 ,

In which κ is a proportionality coefficient that is generally 106

with

with

∇⋅ ψf 1 = 0;

∇⋅ ψf 2 = 0,

(40) (41)

(42)

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111

P. Zheng et al.

[Ks ∇2 + ω 2 M] χ = 0 ,

(43)

Table 1 Material properties for Berea sandstone.

where

⎡ k11 k12 k13 ⎤ ⎢ ⎥ Kp = ⎢ k12 k22 k23 ⎥ , ⎢⎣ k13 k23 k33 ⎥⎦

⎡ μ 0 0⎤ Ks = ⎢ 0 0 0 ⎥ , ⎢ ⎥ ⎣ 0 0 0⎦

φT = ( φs φf 1 φf 2 ),

ρf ⎤ ⎡ ρ ρf ⎢ ⎥ 0 ⎥, ρ m 22 M=⎢ f ⎢⎣ ρf 0 m33 ⎥⎦

χT = ( ψs ψf 1 ψf 2 ),

(44) (45)

with

⎛ 4 ⎞ k11 = ⎜b11 + μ⎟ + b12 χ1 (C1 − χ2 ) + b13 χ1 (C2 − χ2 ), ⎝ 3 ⎠

(46)

k12 = χ1 [b12 (M1 − χ3 ) + b13 (N − χ3 )] = χ1 (C1 − χ2 ),

(47)

k13 = χ1 [b12 (N − χ3 ) + b13 (M2 − χ3 )] = χ1 (C2 − χ2 ),

(48)

k22 = χ1 (M1 − χ3 ),

(49)

k23 = χ1 (N − χ3 ), k33 = χ1 (M2 − χ3 ).

Parameter

Berea sandstone

K , GPa μ , GPa Ks , GPa ρs , kg/m3 K1, GPa (Ks )1, GPa α1 ϕ1 B1 τ1 τ2 v2 k1, m2

6.0a 0.875b 39.0a 2650.0 10.0a 39.0a 0.744 0.178a 0.505 2.17 1.0a 0.0178a

a

b12 N + b13 M2 = C2,

(53)

χ3 (b12 + b13) = χ2 ,

0.001a

From Berryman and Wang 9 From Lewallen and Wang 34

First compressional wave 1/Q

b

(54)

in deriving (47) and (48). It is seen that Kp , Ks , and M are symmetric. Thus, by formulating equilibrium equations in terms of displacements u , w1, w2 , we can obtain two second-order symmetric dynamic systems (the first one is the P-system and the second one the S-system), each of them consisting of three coupled wave equations. 3.2. Dispersion equation for compressional waves and its solution

10

-1

10

-2

10

-3

10

-4

10

-5

10

-6

10

-7

10

-8

Consider time-harmonic plane waves of the form

φ (x) = φ A (x)exp(−i kp⋅x),

(55)

10

Mesoscopic flow present Mesoscopic flow absent

-3

10

-2

10

-1

First compressional wave velocity (m/s)

(56)

where kp ≡ kp . The determinant (56) can be expressed as a cubic polynomial D (x ) with respect to x = (kp / ω)2 , that is

D (x ) = Ax 3 + Bx 2 + Cx + F = 0,

(57)

where the complex coefficients A, B , C , F are defined as 2 A = k11 (k22 k33 − k 23 ) + k12 (k13 k23 − k12 k33) + k13 (k12 k23 − k13 k22 ),

(58)

B = −k11 (k22 m33 + k33 m22 ) + k22 (k13 ρf − k33 ρ) + k23 (k23 ρ − k12 ρf ) + k12 (k12 m33 − k23 ρf + 2k33 ρf ) + k13 ρf (k22 − k23) + k13 (k13 m22 − k23 ρf ), C = ρ (m33 k22 + m22 k33) + m22 (m33 k11 − 2ρf k13) (59)

−2m33 ρf k12 + ρf2 (2k23 − k22 − k33),

(60)

F = −ρm22 m33 + ρf2 (m22 + m33).

(61)

10

0

10

1

10

2

10

3

10

4

10

5

10

5

Frequency (Hz)

where kp is the wave vector, or the propagation vector. Introducing into (42), we can find the dispersion relation for compressible waves to be

det(kp2 Kp − ω 2 M) = 0,

b

1000.0a

η , Pa s

In passing, use has been made of the identities

b

b 10−12)

(6.0 × 10−2) 2.3a

Kf , GPa ρf , kg/m3

(51)

(52)

(1.0 ×

κ , (Gpa s)−1

(50)

b12 M1 + b13 N = C1,

(1.0 × 10−16)

k 2 , m2

2600 2580 2560 2540 2520 2500

Mesoscopic flow present Mesoscopic flow absent

2480 10

-3

10

-2

10

-1

10

0

10

1

10

2

10

3

10

4

Frequency (Hz) Fig. 1. P1-wave speed (a) and inverse quality factor (b) as a function of frequency for fractured Berea sandstone.

From Eq. (57), three complex roots xi = (kpi / ω)2 (i = 1, 2, 3) can be obtained which indicates three distinct types of compressional waves, namely the so-called P1-, P2-, P3-waves, exist in the double porosity medium as previously reported by Berryman and Wang.9 Since Eq. (57) is a cubic equation with complex coefficients, there is no general

arithmetic or algebraic method for finding the roots of such equations except in the case of a cubic equation with real coefficients which roots can be obtained by using Cardan's formula. Therefore, we will use 107

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111

10

1

10

0

10

-1

10

-2

Mesoscopic flow present Mesoscopic flow absent

Third compressional wave 1/Q

10

2

Mesoscopic flow present Mesoscopic flow absent

4.0 3.5 3.0 2.5 2.0

10

10

2

10

1

10

0

10

-1

10

-2

-3

10

-2

10

-1

0

1

2

10 10 10 Frequency (Hz)

10

3

10

4

10

5

10

Mesoscopic flow present Mesoscopic flow absent

10

-3

10

-2

10

-1

0

1

2

10 10 10 Frequency (Hz)

10

3

10

4

10

5

-2

10

-1

10

0

10

1

10

2

10

3

10

4

10

3

10

5

4

10

k13 = k23 = k33 = 0,

5

0

Vpi (ω) =

-3

10

-2

10

-1

10

0

10

1

10

2

10

5

1 ≡ si

1 , xi

i = 1, 2, 3,

(67)

where si (i = 1, 2, 3), reciprocal of the velocity, is the complex slowness. Note that the choice of square root of xi should meet the requirement that the real and imaginary parts of the root are both positive (see, e.g., Berryman and Wang9). In addition, the inverse of the quality factor Q , as the measure of attenuation loss, is defined as

(62)

1 Im[s (ω)] (ω ) = 2 , Q Re[s (ω)]

(68)

with s (ω) being the frequency dependent slowness of an inhomogeneous wave. Using the linear transformation

(63)

where 2 B′ = −k11 k22 m33 + k22 k13 ρf + k12 m33,

10

Fig. 3. P3-wave speed (a) and inverse quality factor (b) as a function of frequency for fractured Berea sandstone.

we thus find two possible starting values to be

C′2 − 4B′F ) , 2B′

Mesoscopic flow present Mesoscopic flow absent

Frequency (Hz)

Newton's iteration algorithm to solve this equation. Since Newton's iteration process is operated in the complex plane (which would be relatively more difficult than that is implemented in the real domain), it is important to choose good starting values for at least two of the three roots of (57). By considering ∇⋅w2 = 0 , which corresponds to the case for which phase 2 is entirely embedded within phase 1, one can choose two possible starting values. In this case, noting

(−C′ ∓

10

15

10

Fig. 2. P2-wave speed (a) and inverse quality factor (b) as a function of frequency for fractured Berea sandstone.

(0) = x1,2

-3

Frequency (Hz)

Third compressional wave velocity (m/s)

Second compressional wave velocity (m/s)

Second compressional wave 1/Q

P. Zheng et al.

φ A = Λp φ Ap ,

(64)

(69)

where

C′ = k22 (ρm33 − ρf2 ) + m22 (m33 k11 − ρf k13) − 2ρf m33 k12.

⎛ 1 1 1 ⎞ Λp = ⎜ Λ21 Λ22 Λ23 ⎟ , ⎜ ⎟ ⎝ Λ31 Λ32 Λ33 ⎠

(65)

Once the final values of the two roots x1 and x2 are determined, the third roots can be obtained directly as (66)

with

Thus, the compressional wave velocities can be determined from each value of x as

Λ2i =

x3 = −(B / A + x1 + x2 ).

108

⎞ ⎛ T (φ Ap ) = ⎜ φpA1 φpA2 φpA3 ⎟ , ⎠ ⎝

k23 xi (k13 xi − ρf ) − (k12 xi − ρf )(k33 xi − m33) 2 2 (k22 xi − m22 )(k33 xi − m33) − k 23 xi

,

(70)

(71)

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111

P. Zheng et al.

0.005 Mesoscopic flow present Mesoscopic flow absent

0.004 Shear wave 1/Q

0.003 0.002 0.001 0.000 -0.001 -3 10

10

-2

10

-1

10

0

10

1

10

2

3

10

3

10

10

4

10

4

10

5

Frequency (Hz)

Shear wave velocity (m/s)

608.5

Mesoscopic flow present Mesoscopic flow absent

608.0 607.5 607.0 606.5 606.0 -3 10

10

-2

10

-1

10

0

10

1

10

2

10

5

Frequency (Hz) Fig. 4. S-wave speed (a) and inverse quality factor (b) as a function of frequency for fractured Berea sandstone.

Λ3i =

k23 xi (k12 xi − ρf ) − (k13 xi − ρf )(k22 xi − m22 ) 2 2 (k22 xi − m22 )(k33 xi − m33) − k 23 xi

, (i = 1, 2, 3),

Fig. 5. Variations of P1-wave speed (a) and inverse quality factor (b) as a function of frequency for fractured Berea sandstone with different fluid transport coefficients.

with sˆ being the complex slowness. The choice of square root of ξ should obey the same rule as for the choice of square root of xi . Using the relation

(72)

the coupled P-system (42) can be decoupled into three scalar Helmholtz equations

∇2 φpiA

kpi2 φpiA

+

= 0,

(i = 1, 2, 3),

χA = Λs ψ A s,

(73)

(78)

where

which are obviously the wave equations for P1-, P2-, and P3-waves, respectively.

ΛTs = (1 Λf 1 Λf 2 ),

(79)

with the amplitudes radios Λf 1, Λf 2 given by 3.3. Dispersion equation for shear wave and its solution

Λf 1 = −m22 / ρf ,

χ (x) =

∇2 ψ sA + ks2 ψ sA = 0,

(74)



ω 2 M)

= 0,

(75)

4. Examples

where ks ≡ ks . From the determinant (75), we have 2 μm22 m33 ξ − m11 m22 m33 + m12 (m22 + m33) = 0,

Vs (ω) =

In this section, we perform some numerical calculations to show the effect of mesoscopic flow on the attenuation and velocity dispersion of body waves (P- and S-waves) in double-porosity materials. In our subsequent examples, phase 2 will be treated as throughgoing joints embedded within the host phase 1 and the double-porosity material is assumed to be made up of a certain rock, i.e. the Berea sandstone. The parameters used in the calculations are taken from Table 1. Since the present model is used to account for the attenuation in the seismic

(76)

/ ω )2 .

Thus, the phase velocity of the shear wave is found to

1 1 ≡ , sˆ ξ

(77)

where ξ = (ks be

(81)

which corresponds to the wave equation for S-wave.

Substitution into (43) leads to

det(ks2 Ks

(80)

the coupled S-system (43) can be reduced to

In the same manner let

χA (x)exp(−i ks⋅x).

Λf 2 = −m33 / ρf ,

109

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111

P. Zheng et al.

First compressional wave velocity (m/s)

First compressional wave 1/Q

10

-1

10

-2

10

-3

10

-4

10

ambient water oil hot water

-3

10

-2

10

-1

0

1

10 10 10 Frequency (Hz)

2

10

3

10

4

2620 2600 2580

ambient water oil hot water

2560 2540 2520 2500 2480 10

-3

10

-2

10

-1

0

1

10 10 10 Frequency (Hz)

2

10

3

10

4

Fig. 7. Variations of P1-wave speed (a) and inverse quality factor (b) as a function of frequency for fractured Berea sandstone with different matrix porosity levels.

Fig. 6. Variations of P1-wave speed (a) and inverse quality factor (b) as a function of frequency for fractured Berea sandstone with different pore-fluid viscosities.

Furthermore, from Fig. 2, we observe that P2 wave is diffusive at low frequencies, but becomes propagatory at about 10 kHz. This fact is consistent with findings previously reported by Berryman and Wang.9 Fig. 4 shows the attenuation and dispersion of shear wave (S-wave) for fractured Berea sandstone. At first sight, it is unexpected to find that the effect of mesoscopic flow on S-wave is completely unobservable. But by re-checking the dispersion equation for shear wave (76), we understand immediately that the presence of mesoscopic flow, of course, has no influence on the attenuation curve since the incorporation of mesoscopic flow has no influence on the matrix elements of M which is defined in (44). The dependence of attenuation and dispersion of P1-wave on the values of fluid transport coefficient for fractured Berea sandstone is plotted in Fig. 5. High κ values means high rates of fluid flow between phases and, consequently, shortens the time required for pressure gradients between the matrix and fracture to establish equilibrium in a wave cycle. It is indeed that the relaxation peaks induced by mesoscopic flow shifts with the values of κ , as observed from Fig. 5. As expected, the peak to the left in the attenuation curve occurs at higher frequencies with higher κ values. It is also noted that the levels of attenuation, Q−1, are insensitive to the κ values. Fig. 6 displays the effect of pore-fluid viscosity on the attenuation and dispersion of P1-wave. In the example, the pore-fluid viscosities are respectively η = 10−3 Pa s (ambient water), 2 × 10−4 Pa s (hot water), and 5 × 10−3 Pa s (oil). From Fig. 6, we observe that the effect of

band, we assume that the proportionality coefficient κ is frequency independent and take a constant value. Therefore, we let κ take its static value, namely κ (ω) = κ0 , in our examples. Figs. 1, 2 and 3 show the attenuation and dispersion of P1-, P2-, P3waves, respectively, for a fractured Berea sandstone. To show the effects of fluid flow between two phases on wave attenuation, we compare attenuation curves computed with such mesoscopic flow present to those with the mesoscopic flow absent. It is clearly observed that there are two relaxation peaks in the attenuation curve of P1-wave (solid line), whereas in the case of no fluid flow between phases (which precisely corresponds to the model studied by Berryman and Wang9) there are only one relaxation peak known as the Biot loss in the attenuation curve (dash line). The peak to the left, namely the dominating one, centered in the region 0.01–10 Hz is due to the socalled mesoscopic flow which has been previously reported by Pride and Wang5,12,13 by treating phase 2 as a porous continuum. This effect can also be clearly observed from the attenuation curves of P2- and P3waves. In Fig. 1B, we note that the predicted level of intrinsic loss induced by mesoscopic flow is roughly in the range 10−2 < Q−1 < 10−1 which is about at least one-order-of magnitude higher than that induced by Biot loss (which is centered in the region 1–10 kHz). This fact can be further confirmed from the attenuation curve of P2-wave. Since the attenuation level of P3 wave is almost Q−1 > 2 observed from Fig. 3B, we thus conclude that the P3 wave is essentially diffusive rather than propagatory at all frequencies in the range considered. 110

International Journal of Rock Mechanics & Mining Sciences 91 (2017) 104–111

P. Zheng et al.

Low-frequency range. J Acoust Soc Am. 1956;28:168–178. 2 Biot MA. Theory of propagation of elastic waves in a fluid-saturated porous solid. II. Higher frequency range. J Acoust Soc Am. 1956;28:179–191. 3 White JE. Computed seismic speeds and attenuation in rocks with partial gas saturation. Geophysics. 1975;40:224–232. 4 Dvorkin J, Mavko G, Nur A. Squirt flow in fully saturated rocks. Geophysics.. 1975;60:97–107. 5 Pride SR, Berryman JG, Harris JM. Seismic attenuation due to wave-induced flow. J Geophys Res. 2004;109:B01201. 6 White JE, Mikhaylova NG, Lyakhovitskiy FM. Low-frequency seismic waves in fluidsaturated layered rocks. Izv Acad Sci USSR Phys Solid Earth. 1975;11:654–659. 7 Norris AN. Low-frequency dispersion and attenuation in partially saturated rocks. J Acoust Soc Am. 1993;94:359–370. 8 Gelinsky S, Shapiro SA. Poroelastic backus averaging for anisotropic layered fluidand gas-saturated sediments. Geophysics. 1997;62:1867–1878. 9 Berryman JG, Wang HF. Elastic wave propagation and attenuation in a doubleporosity dual-permeability medium. Int J Rock Mech Min Sci. 2000;37:63–78. 10 Johnson DL. Theory of frequency dependent acoustics in patchy saturated porous media. J Acoust Soc Am. 2001;110:682–694. 11 Gurevich B. Elastic properties of saturated porous rocks with aligned fractures. J appl Geophys. 2003;30:203–218. 12 Pride SR, Berryman JG. Linear dynamics of double-porosity and dual-permeability materials. I. Governing equations and acoustic attenuation. Phys Rev E. 2003;68:036603. 13 Pride SR, Berryman JG. Linear dynamics of double-porosity and dual-permeability materials. II. Fluid transport equations. Phys Rev E. 2003;68:036604. 14 Brajanovski M, Gurevich B, Schoenberg M. A model for P-wave attenuation and dispersion in a porous medium permeated by aligned fractures. Geophys J Int. 2005;163:372–384. 15 Berryman JG. Seismic waves in rocks with fluids and fractures. Geophys J Int. 2007;171:954–974. 16 Gurevich B, Brajanovski M, Galvin R, Müller TM. Toms-Stewart J. P-wave dispersion and attenuation in fractured and porous reservoirs – Poroelasticity approach. Geophys Prospect. 2009;57:225–237. 17 Ba J, Carcione JM, Nie JX. Biot-Rayleigh theory of wave propagation in doubleporosity media. J Geophys Res. 2011;116:B06202. 18 Müller TM, Gurevich B, Lebedev M. Seismic wave attenuation and dispersion resulting from wave-induced flow in rocks – A review. Geophysics.. 2010;75 75A147–75A164. 19 Warren JE, Root PJ. The behavior of naturally fractured reservoirs. Soc Pet Eng J. 1963;3:245–255. 20 WilsonRK, AifantisEC.On the theory of consolidation with double porosity. Int J Eng Sci. 1982; 20: 1009–1035. (Correction, Int J Eng Sci. 1983; 21: 279–281.) 21 Wilson RK, Aifantis EC. A double porosity model for acoustic wave propagation in fractured-porous rock. Int J Eng Sci. 1984;22:1209–1217. 22 Khaled MY, Beskos DE, Aifantis EC. On the theory of consolidation with double porosity. III. A finite element formulation. Int J Numer Anal Methods Geomech. 1984;8:101–123. 23 Beskos DE, Aifantis EC. On the theory of consolidation with double porosity II. Int J Eng Sci. 1986;24:1697–1716. 24 Elsworth D, Bai M. Flow-deformation response of dual-porosity media. ASCE J Geotech Eng. 1992;118:107–124. 25 Auriault JL, Boutin C. Deformable media with double porosity – Quasi-statics. I: Coupling effects. Transp Porous Med. 1993;10:153–169. 26 Bai M, Elsworth D, Roegiers J-C. Modelling of naturally fractured reservoirs using deformation dependent flow mechanism. Int J Rock Mech Min Sci Geomech Abstr. 1993;30:1185–1191. 27 Tuncay K, Corapciouglu MY. Wave propagation in fractured porous media. Transp Porous Med. 1996;23:237–258. 28 Corapciouglu MY, Tuncay K. Wave propagation in fractured porous media saturated by two immiscible fluidsThimus J-F, Abousleiman Y, Cheng AH-D, Coussy O, Detournay E, eds. Poromechanics: a tribute to Maurice A Biot, Rotterdam: Balkema; 1998, pp.197–203. 29 Olny X, Boutin C. Acoustic wave propagation in double porosity media. J Acoust Soc Am. 2003;114:73–89. 30 Sams MS, Neep JP, Worthington MH, King MS. The measurement of velocity dispersion and frequency-dependent intrinsic attenuation in sedimentary rocks. Geophysics. 1997;62:1456–1464. 31 Berryman JG, Wang HF. The elastic coefficients of double-porosity models for fluid transport in jointed rock. J Geophys Res. 1995;100:24611–24627. 32 Berryman JG, Pride SR. Models for computing geomechanical constants of doubleporosity materials from the constituents’ properties. J Geophys Res. 2002. http:// dx.doi.org/10.1029/2000JB000108. 33 Zheng P, Gao Z, Ding B.The elastic coefficients of double-porosity materials: a revisit. Transp Porous Med. 2016; 111: 555–571. (Erratum, Transp Porous Med. 2016;112: 783–784.) 34 Lewallen KT, Wang HF. Consolidation of a double-porosity medium. Int J Solids Struct. 1998;35:4845–4867.

viscosity is to shift the relaxation peaks induced by Biot loss remark2 / D ,12 ably. Since the peaks induced by Biot loss occurs roughly at f = vp1 with vp1 being the P1-wave velocity and D is the pore-pressure diffusivity, and since the diffusion coefficient D is inversely proportional to viscosity, the Biot-loss peaks shift rightward with increasing viscosity. On the other hand, the peaks induced by mesoscopic flow seem rather insensitive to viscosities changes. Fig. 7 shows the attenuation and dispersion of P1-wave for an increase in matrix porosity assuming no other parametric changes. The degree of matrix porosity depends largely on the consolidation degree of sandstones, namely, a poorly-consolidated sandstone is, in general, with a high degree of matrix porosity while a well-consolidated one is, in most cases, with a low degree of matrix porosity. It is observed form Fig. 7 that the effect of matrix porosity is to change the attenuation and dispersion levels of P1-wave considerably. Since larger porosity allows for lager fluid flow from fractures into matrix pores and vice versa, the increasing matrix porosity causes an increase in the magnitude of attenuation and dispersion. For Biot loss, the attenuation level has just the opposite dependence on the matrix porosity than does the mesoscopic flow. 5. Discussion and conclusions We have presented an analytical method to exactly solve the field equations for double-porosity materials. The proposed approach is based on the reformulation of the field equations in terms of solid displacement u and the relative fluid-solid displacements w1, w2 . With this new form of wave equations, we have obtained two second-order symmetric dynamic systems, namely, the P-system for compressible waves and the S-system for shear wave. To understand the loss mechanism caused by mesoscopic flow, we have compared the attenuation cures with the mesoscopic flow present to those with the mesoscopic flow absent. Numerical calculations indicate that the amount of attenuation caused by mesoscopic flow is at least oneorder-of magnitude higher than that induced by Biot loss for fast compressional wave, and thus confirms the previous findings reported by Pride and Berryman using an approximation scheme. Since the theoretical model for the frequency-dependent proportionality coefficient κ (ω) is not available for the case of phase 2 being throughgoing joints, we assume that κ (ω) is frequency independent and take its static value. The impact of such static assumption on the final results needs to be clarified. We leave this analysis to future work. The present work is limited by the assumption of an isotropic double-porosity medium. In contrast to the present isotropic model to describe fractures, Gurevich et al.16 presented two anisotropic models to account for the effect of fractures. In the first model, the factures are modeled as planes of weakness and in the second model fractures are treated as thin penny-shaped voids of finite radius. The effect of fracture-induced anisotropy is expressed through anisotropic Gassmann's equations.11 An extension of the present study to anisotropic media, e.g. a transversely isotropic double-porosity material, will be also left to the future. Acknowledgements The work of PZ and BD was supported by the National Natural Science Foundation of China under Grant Nos. 11402150 and 51478435. The work of XS was supported by Shanghai Sailing Program under Grant No. 16YF 1408000 and Natural Science Foundation of Shanghai under Grant No. 16ZF 1423600. References 1 Biot MA. Theory of propagation of elastic waves in a fluid-saturated porous solid. I.

111