M O D E L I N G E L A S T I C W A V E V E L O C I T I E S I N P O R O U S M E D I A : F R E Q U E N C Y - D E P E N D E N T E F F E C T S O F H E T E R O G E N E I T Y A T T H E P O R E - A N D P A T C H - S C A L E By Richard Taylor B. Sc. (Hons. Physics) University of British Columbia A T H E S I S S U B M I T T E D I N P A R T I A L F U L F I L L M E N T O F T H E R E Q U I R E M E N T S F O R T H E D E G R E E O F M A S T E R O F S C I E N C E in T H E F A C U L T Y O F G R A D U A T E S T U D I E S D E P A R T M E N T O F E A R T H A N D O C E A N S C I E N C E S and I N S T I T U T E O F A P P L I E D M A T H E M A T I C S We accept this thesis as conforming to the required standard T H E U N I V E R S I T Y O F B R I T I S H C O L U M B I A August 1999 © Richard Taylor, 1999 In presenting this thesis in partial fulfillment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of Earth & Ocean Sciences The University of British Columbia Vancouver, Canada Date: Abstract The prediction of relationships between elastic wave velocities in a porous medium and the properties of the fluid and solid constituents therein is a longstanding problem in geophysical exploration. Previous authors have shown that such relationships depend on both the wave frequency and degree of heterogeneity present. This frequency depen-dence arises via the state of relaxation of the pore fluids. At sufficiently low frequencies, deformation-induced flow leads to equilibration of the fluid pressure, and the pore fluids are said to be in a "relaxed" state. At sufficiently high frequencies, there is insufficient time for equilibration to occur, and the fluids are said to be "unrelaxed". Current models of elastic wave velocities in porous media are, for the most part, confined to either the relaxed limit (e.g., the poroelastic Biot-Gassmann theory) or the unrelaxed limit (e.g., inclusion-based effective medium theory). In this thesis we incor-porate an explicit description of the relaxation mechanism into inclusion-based effective medium theory, so as to extend the theory toward the relaxed limit. Analysis of the mechanisms of relaxation leads to a description of the effective elastic behavior of the porous medium in terms of effective complex elastic moduli of the medium's constituents. Previous authors have identified two distinct scales of fluid distribution heterogene-ity: the pore scale and the patch scale. Accordingly, we treat these scales separately, describing relaxation in terms of Poiseuille flow at the pore scale, and Darcy's law at the patch scale. The results of our analyses are effective medium theories that provide a consistent approach to the prediction of elastic wave velocities, as well as attenuation due to the relaxation mechanism, over a broad range of frequencies and length scales of heterogeneity. In particular, our model is applicable to the regime where the pore fluids are in a state of relaxation intermediate between the completely relaxed and completely unrelaxed end members. ii Table of Contents Abstract ii List of Tables v List of Figures vi Acknowledgement ix 1 Introduction 1 2 Heterogeneity and Relaxation at the Pore Scale 8 2.1 Introduction 8 2.2 Introduction to Inclusion-Based Effective Medium Theory 13 2.3 Defining Effective Elastic Moduli for Inclusions 17 2.4 A Simple Pore-to-Pore Relaxation Model 20 2.4.1 Model Development 21 2.4.2 Modeling Results and Discussion 31 2.4.3 Incorporating Pore Compliance Effects 39 2.4.4 Limitations of the Pore-to-Pore Relaxation Model 42 2.5 Incorporating Long-Range Fluid Pressure Communication 43 2.5.1 Model Development 43 2.5.2 Modeling Results and Discussion 45 2.6 An Explicit Pore Network Model 47 2.6.1 Model Development 47 2.6.2 Modeling Results and Discussion 50 2.7 Discussion and Comparison of the Pore-Scale Relaxation Models 54 iii 3 Heterogeneity and Relaxation at the Patch Scale 57 3.1 Introduction 57 3.2 Equations of Pressure Diffusion in a Porous Medium 62 3.2.1 Pressure Diffusion Equation for Imposed Strain 64 .3.2.2 Pressure Diffusion Equation for Imposed Stress 65 3.3 Effective Bulk Modulus of a Patch 67 3.3.1 Solution of the Pressure Diffusion Equation for a Patch 67 3.3.2 Definition of the Effective Bulk Modulus of a Patch 74 3.3.3 Accounting for Partial Saturation 77 3.4 Discussion 80 3.4.1 Frequency Response of a Fluid-Saturated Patch 80 3.4.2 Effect of Permeability Contrast Across a Patch Boundary 81 3.4.3 Effect of Patch Size 84 3.4.4 Multiple Scales of Heterogeneity 87 3.5 An Application: Modeling Saturation Hysteresis 90 3.6 Summary 99 4 Conclusion 101 References 109 Appendices 112 A Probability Distribution of Distances Between Neighboring Pores 112 List of Notation 117 iv List of Tables 2.1 Physical properties of the constituents of a model sandstone 32 2.2 Parameters used to quantify the pore space and hydraulic properties of a model sandstone 32 3.1 Physical parameters for Spirit River sandstone saturated wi th water and air 85 v List of Figures 1.1 Computed P-Wave and S-wave velocities vs. water saturation for a model sandstone, illustrating the end members corresponding to the poroelas-tic model of Gassmann (1951) and the inclusion-based model of Berry-man (1980a; 1980b) 4 2.1 A sample of a porous solid as conceptualized in the inclusion-based effective medium theory approach 13 2.2 Incorporating fluid pressure communication into an inclusion-based model. 18 2.3 Schematic illustration of the replacement of hydraulically connected inclu-sions with equivalent inclusions having effective bulk moduli 20 2.4 Schematic illustration of an idealized pore-to-pore relaxation model: two fluid-filled inclusions connected by a narrow duct 21 2.5 Graphs of the real and imaginary parts of the effective bulk modulus K* vs. ui/u)0 for the pore-to-pore relaxation model 24 2.6 Phase angle 6 (by which pressure leads dilatation) vs. u>/u0 , for the pore-to-pore relaxation model 24 2.7 Probability distribution for the distance r from the center of a given water-filled inclusion to the center of its nearest gas-filled neighbor, for various levels of water saturation Sw 27 2.8 Shear-induced flow between two plates 28 2.9 Computed effective elastic moduli vs. water saturation at various frequen-cies for a model sandstone, using the pore-to-pore relaxation model. . . . 32 2.10 Computed P-wave velocity vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model 34 2.11 Computed P-wave Q'1 vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model 34 2.12 Computed S-Wave velocity vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model 37 vi 2.13 Computed S-wave Q 1 vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model 37 2.14 Computed P-wave velocity vs. frequency at various saturations for a model sandstone, using the pore-to-pore relaxation model 38 2.15 Computed Q'1 vs. frequency at various saturations for P-waves in a model sandstone, using the pore-to-pore relaxation model 38 2.16 Computed P-wave velocity vs. water saturation at various frequencies for a model sandstone, using the modified pore-to-pore relaxation model. . . 46 2.17 Computed P-wave Q"1 vs. water saturation at various frequencies for a model sandstone, using the modified pore-to-pore relaxation model. . . . 46 2.18 Schematic illustration of the explicit pore network relaxation model: every inclusion is connected to a number of other inclusions by cylindrical ducts through which fluid can flow 48 2.19 Computed P-wave velocity vs. water saturation at various frequencies for a model sandstone, using the explicit pore network relaxation model. . . 52 2.20 Computed P-wave Q~l vs. water saturation at various frequencies for a model sandstone, using the explicit pore network relaxation model. . . . 52 2.21 Computed P-wave velocity vs. frequency at 80% water saturation for a model sandstone, using the explicit pore network relaxation model. . . . 53 3.1 Schematic illustration of a fluid-saturated porous medium, of which we consider a particular small volume V 62 3.2 Schematic illustration of a model porous medium with a spherical patch of radius embedded in an infinite background medium 68 3.3 Graphs of the real part of the solution Pf(r)/po at various phases during a wave cycle, for UT = 103 72 3.4 Graphs of the real part of the solution Pf(r)/po at various phases during a wave cycle, for LOT = 10 _ 1 73 3.5 Graph of the real part of K*/Kd vs. <J0T/(2-K) for the same values of the various parameters as used in Figures 3.3-3.4, with oti = 0.5 81 3.6 Real part of K* / Kd vs. U)T/(2IT), for a water-saturated patch in a typical sandstone that is otherwise saturated with air, and for various values of the permeability contrast K\jK2 83 vii 3.7 Graphs of the real and imaginary parts of K*/Kd vs. patch radius at wave frequencies of 100 Hz and 50 kHz, for a spherical water-saturated patch in a Spirit River sandstone that is otherwise saturated with air 85 3.8 P-wave velocity vs. wave frequency, computed using the patch-scale relax-ation model for a particular low-permeability sandstone having 80% water saturation, with half of the water in patches with radius 0.2 mm and the other half with radius 1 cm 89 3.9 P-wave attenuation vs. wave frequency for the same medium considered in Figure 3.8 89 3.10 Elastic P-wave velocity vs. water saturation in Spirit River sandstone (data given in Table 3.1) at 100 Hz, computed using the patch-scale relaxation model, illustrating the dependence on the length scale of saturation het-erogeneity 95 3.11 Elastic P-wave velocity vs. water saturation in Spirit River sandstone (data given in Table 3.1) at 10 kHz, computed using the patch-scale relaxation model, illustrating the dependence on the length scale of saturation het-erogeneity 95 3.12 Computed saturation hysteresis of elastic P-wave velocity in Spirit River sandstone, at 100 Hz and 10 kHz 97 3.13 Computed saturation hysteresis of P-wave attenuation Q~1 for Spirit River sandstone, at 100 Hz and 10 kHz 97 4.1 Schematic illustration of the ranges of wave frequency and length scale of heterogeneity for which various effective medium theories are valid, for a particular typical low-permeability rock 106 4.2 Schematic illustration of the ranges of wave frequency and length scale of heterogeneity for which the pore- and patch-scale relaxation models developed in this thesis are valid, for the same medium considered in Figure 4.1 107 A . l Geometrical construction used in deriving the probability distribution for the nearest-neighbor distance in a spatially random distribution of inclusions. 113 A.2 Graph of the probability distribution for the nearest-neighbor distance in a spatially random distribution of inclusions 116 viii A c k n o w l e d g e m e n t Foremost, I extend my gratitude to my supervisor, Rosemary Knight, for making this research possible, and for her indispensable insight, helpful criticism, flexibility, and confidence in my abilities. Thanks also to the other members of the Rock Physics group at UBC, who have contributed their energies and insight to my research. I am particularly indebted to Christian Scullard, for much fruitful conversation and for his enthusiasm in providing an independent derivation of the results of Appendix A. Many of my friends have had a part in bringing this thesis to completion, all of whom deserve my gratitude. Among these, I especially thank my wife, Amy, for her continuing support and encouragement. Much credit is also due to my teachers, too numerous to name, whose contributions are present in every aspect of this thesis. This research was supported by an NSERC Industrially-Oriented Research Grant to Rosemary Knight with funding from Imperial Oil, Petro-Canada and Western Atlas. I have also been supported in part by a University Graduate Fellowship from the University of British Columbia. ix Chapter 1 Introduction Modeling the dependence of elastic wave velocities in a fluid-saturated porous medium on the properties of the medium's constituents (and particularly those of its saturat-ing fluids) is a long-standing problem in geophysics. This is a critical issue in seismic exploration, where such models are required in order that measurements of the elastic wave velocities can be used to indirectly measure other properties of interest in situ. For example, it is a common practice to use velocity measurements to infer the relative pro-portions of liquid and gas in porous rock, by assuming a particular relationship between the level of liquid saturation and the elastic wave velocities. An approximation which frequently plays a role in the development of such models, and one which we use throughout this thesis, is that the elastic behavior of the porous medium is adequately described by effective medium theory. In the effective medium the-ory approximation, the length scale of heterogeneity of the porous medium is assumed to be much smaller than the elastic wavelength. Scattering of elastic waves by hetero-geneities can then be neglected, and the elastic response of the porous composite can be approximated as that of an equivalent homogeneous solid. The problem of modeling elastic wave velocities is then reduced to that of estimating the elastic moduli of this effective medium, from which the elastic wave velocities can be computed using known formulae. The chief difficulty that arises when fluid-saturated media are considered, is in esti-mating the contribution of the various pore fluids to the elastic response of the composite. 1 Chapter 1. Introduction 2 The problem involves coupled fluid dynamics and elasticity, with non-trivial (and usu-ally poorly quantified) pore space geometry. To render this problem tractable, previous authors have typically employed one of two simplifying assumptions. These are: (1) no flow of the pore fluids is induced by the elastic deformation of the medium; and (2) on a time scale much shorter than the wave period, the pore fluids flow so as to equilibrate the fluid pressure throughout the pore space. These assumptions effectively define two "end members" of the range of elastic behavior possible in a given fluid-saturated porous medium. In case (1) we say that the pore fluids are completely "unrelaxed" and con-tribute maximally to the elastic response of the medium. This is in contrast with case (2) where the pore fluids are said to be completely "relaxed", and contribute minimally to the elastic response of the medium (i.e., the effective elastic moduli of the medium, and hence the elastic wave velocities, are minimized). At one end member we have the assumption that the pore fluids are in a completely unrelaxed state. That is, there is insufficient time during a wave period for significant flow to occur to communicate the fluid pressure throughout the pore space. This as-sumption therefore requires that the wave period be sufficiently short (i.e., that the frequency be sufficiently high) to restrict flow of the pore fluids. This is the assumption used in inclusion-based effective medium theory models like those of Berryman (1980a; 1980b), Kuster and Toksoz (1974a), Budiansky and O'Connell (1976), and le Ravalec and Gueguen (1996). In these models the pore space is considered to be composed of hydraulically isolated fluid-filled inclusions in a solid porous matrix. The success of these models has been demonstrated (Coyner, 1984; Murphy, 1985; Wang et ai, 1991) in mod-eling laboratory data obtained at ultrasonic frequencies, where the assumption of no fluid pressure communication is generally thought to be valid. At the other end member we have the assumption that the pore fluids are in a com-pletely relaxed state. That is, during a wave period there is sufficient time for any Chapter 1. Introduction 3 non-uniform fluid pressure to induce flow so as to equilibrate the pore fluid. In a medium with non-zero hydraulic permeability this equilibration/relaxation occurs on a finite time scale, so that this assumption will be valid if the wave period is sufficiently long (i.e., if the frequency is sufficiently low). The relaxation time scale depends on the length scale of fluid distribution heterogeneity, such that the relaxation occurs more quickly if the spatial distribution of pore fluids is more homogeneous. The state of relaxation of the pore fluids therefore depends not only on the wave frequency but also on the length scale of heterogeneity of the fluid distribution. That is, the assumption of complete relaxation can be valid at higher frequencies if the pore fluids are distributed more homogeneously. The assumption of complete relaxation is used in poroelasticity models like those of Gassmann (1951) and Biot (1956a; 1956b). The Gassmann poroelastic theory in particular is frequently used in seismic exploration, where laboratory data (Murphy, 1985) have demonstrated its validity for a number of porous materials. Gassmann's formula (Mavko et al., 1996, p. 185), gives the effective bulk modulus, Keff, of a porous medium saturated with a fluid having bulk modulus Kf. Here 0 is the porosity, Kd is the dry frame bulk modulus (i.e., the bulk modulus of the porous matrix in the absence of pore fluids), and KS is the bulk modulus of the mineral solid of which the matrix is composed. The parameter a is the "poroelastic parameter" defined as a = 1 — KD/KS. In the Gassmann approach, the effective shear modulus, pefi, of the fluid-saturated medium is assumed to be equal to the dry frame shear modulus, i.e., «ef f = « d + cj) a — 0 ' (1.1) Kf Ks Peff — Pd- (1.2) In the case of multiple pore fluids having bulk moduli «/(2),---i Kf(n) and levels of Chapter 1. Introduction 4 2.2 0.88 w 0.84 E ^ » 0.82 0.86 1.6 sJJiapn 0.8 1.5 0.78 0 0.5 S 0 0.5 s w Figure 1.1: Computed P-Wave and S-wave velocities vs. water saturation for a model sandstone, illustrating the end members corresponding to the poroelastic model of Gassmann (1951) (—) and the inclusion-based model of Berryman (1980a; 1980b) (—). partial saturation Si, S2, ... , Sn (with Yli=i^i = -0' the- assumption of fluid pressure equilibration allows one to define the "effective fluid" bulk modulus, K*f, defined by (Domenico, 1977) which is then used in the place of Kf in equation (1.1). Consider, for example, the relationship between the P- and S-wave velocities (vp and vs, respectively) and the level of water saturation Sw (i.e., the volume fraction of the pore space that is filled with water), in a typical water/air-saturated sandstone. We use the data for the model sandstone considered in Endres and Knight (1989), except that we assume that all pores have spherical geometry. Figure 1.1 presents graphs of vp and vs vs. Sw for the two end members corresponding to the Gassmann (1951) and Berryman (1980a; 1980b) effective medium theories. It is clear that the two end members are considerably different in form and magnitude, with the velocities computed with the (1.3) Chapter 1. Introduction 5 Berryman model (corresponding to the unrelaxed state of the pore fluids) being higher than those of the Gassmann model (corresponding to the relaxed state), for most of the range of water saturation. Hence, in applications, the choice of assumption for the state of the pore fluids (either relaxed or unrelaxed) is critical. Moreover, the difference between the two models is so great that we expect there to exist geophysically relevant situations in which neither assumption is valid. Since, at present, there exists no model that is valid at both end members, elastic wave velocities cannot be consistently modeled across the full range of possible behavior. The present inconsistency between the poroelastic and inclusion-based models has been addressed by previous authors and has led to attempts to reconcile them (Endres and Knight, 1997). An issue of particular importance is that there exists no theory which can be used to model both ultrasonic wave velocities measured in the laboratory and field seismic velocities, so that it is difficult to apply laboratory findings to seismic exploration. To this end, Endres and Knight (1997) have extended the inclusion-based approach to the limit of complete fluid pressure relaxation, but unfortunately their model is restricted to the completely relaxed and completely unrelaxed limits. The purpose of this thesis is to attempt to continuously bridge the gap between the poroelastic and inclusion-based models, and develop a model that is both consistent with these end members and able to estimate elastic behavior in the intermediate regime. We accomplish this by incorporating an explicit description of the physical mechanism of fluid pressure relaxation/communication into inclusion-based effective medium theory. Our focus is on modeling the frequency-dependent effects of heterogeneity, of either the distribution of pore fluids or the properties of the medium—the key feature being the inhomogeneous fluid pressure that is induced during elastic deformation. Laboratory data (Cadoret et al, 1995) suggests that heterogeneity at a variety of length scales is inherent in many fluid-saturated porous media, and can depend on the nature of the Chapter 1. Introduction 6 drainage/imbibition process. Our results can be used to model elastic wave velocities and attenuation over a broad range of wave frequencies and length scales of heterogeneity. Furthermore (and importantly) our results are consistent with the previous low- and high-frequency end members (i.e., the poroelastic and inclusion-based models). In Chapter 2 we focus on heterogeneity at the pore scale. Here fluid pressure relax-ation/communication is governed by flow between individual pores, which depends on pore-scale properties of the porous medium. We develop a pore-scale flow model, based on Poiseuille flow, that allows us to predict relationships like those in Figure 1.1 in the in-termediate regime between the end members we have identified. Furthermore, our model allows us to predict the viscous losses (and hence the wave attenuation) associated with the pore-scale relaxation process. In Chapter 3 we focus on heterogeneity at the "patch" or macroscopic scale, where the pore fluids saturate the medium in "patches" of a size on the order of many pore diameters. Here fluid pressure relaxation/communication is governed by pressured diffu-sion between patches at the macroscopic scale, which depends on macroscopic properties of the porous medium. We develop a patch-scale flow model, based on Darcy's law for flow in porous media, that allows us to predict wave velocities in the intermediate regime between the end members in the case where the heterogeneity of the fluid distribution is at the patch scale only. We account for the elastic effects of patch-scale heterogeneity by using an inclusion-based approach similar to that of le Ravalec and Gueguen (1996). Again, our model also allows us to predict the viscous losses and wave attenuation asso-ciated with the patch-scale relaxation process. Our approach at both scales of heterogeneity is to introduce complex effective elastic moduli to describe the elastic/anelastic behavior of particular fluid-saturated regions (individual pores in the case of pore-scale heterogeneity, and "patches" in the case of patch-scale heterogeneity). An inclusion-based formulation (Berryman's (1980a; 1980b) Chapter 1. Introduction 7 e f f ec t i ve m e d i u m t h e o r y ) i s t h e n u s e d t o c o m p u t e t h e e f f e c t i v e e l a s t i c m o d u l i o f t h e c o m p o s i t e . U n d e r l y i n g o u r a p p r o a c h i s t h e a s s u m p t i o n t h a t t h e h y d r a u l i c p r o b l e m o f f l u i d p r e s s u r e r e l a x a t i o n c a n b e d e c o u p l e d f r o m t h e e l a s t i c p r o b l e m . T h a t i s , t h e h y d r a u l i c p r o b l e m i s s o l v e d s e p a r a t e l y i n o r d e r t o d e f i n e e f f e c t i v e m o d u l i f o r e a c h p o r e o r " p a t c h " ; t h e s e e f f e c t i v e m o d u l i a r e t h e n c o m b i n e d i n a s e p a r a t e e f f e c t i v e m e d i u m c o m p u t a t i o n . Chapter 2 Heterogeneity and Relaxation at the Pore Scale 2 . 1 Introduction An important issue in seismic exploration is the determination of relationships between elastic wave velocities and other parameters of interest such as the elastic properties of the various fluid and mineral solid constituents present, the relative abundances of the various pore fluids, and the spatial distribution of these fluids. Usually it is these parameters that one wishes to measure, whereas in practice it is possible to obtain direct measurements of the elastic wave velocities only. The assumed relationship between elastic wave velocities and the parameters of interest then plays a critical role in estimating these parameters from velocity measurements. In the present work we are specifically interested in the relationship between the level of liquid saturation and the elastic P-wave velocity. In particular, we focus on modeling the effects of spatial heterogeneity. Heterogeneity on a variety of length scales, from pore-scale to reservoir-scale, is found in most naturally occurring geologic systems. There is substantial experimental (Knight and Nolen-Hoeksema, 1990; Cadoret, 1993) and theoretical (Endres and Knight, 1989; le Ravalec et al, 1996) evidence that such heterogeneity can have a significant effect on elastic wave velocities. In the theoretical treatment by Endres and Knight (1989) of the relationship between P-wave velocity and water saturation, two "end members" are identified—corresponding to the extreme cases of homogeneous and heterogeneous distributions of fluids within the pore space—each 8 Chapter 2. Heterogeneity and Relaxation at the Pore Scale 9 exhibiting a radically different dependence of P-wave velocity on saturation. To avoid confusion, it should be made clear that the term "heterogeneity" as used in this work generally implies heterogeneity at length scales much smaller than the elastic wavelength. The theoretical treatment in this work and many others in fact relies on the "effective medium" assumption that the length scale of heterogeneity is much less than the wavelength, so that a porous medium appears to be homogeneous at the wavelength scale. Wave scattering can then be neglected, and elastic wave velocities are computed as they would be for a truly homogeneous sample having the same elastic moduli as the (inhomogeneous) porous medium. Thus some authors refer to the "homogeneous effective medium", a term which denotes the length scale of heterogeneity/homogeneity relative to the wavelength, but does not necessarily imply homogeneity at sub-wavelength scales. The concept of spatial heterogeneity as it applies to modeling elastic wave velocities in porous media is a relative term that depends on the wave frequency. The wave period determines the time scale of the elastic deformation, and it is the elastic response of a porous medium in that time window of deformation that defines whether it should be treated as a "homogeneous" or "heterogeneous" system. Consider, for example, a porous medium in which for some reason (such as inhomogeneity of the pore fluid distribution or of the pore space geometry) the fluid pressure induced by an elastic deformation is initially not uniform throughout the pore space. The scale of heterogeneity could range from pore-scale to reservoir-scale and might, realistically, be present at all scales. Under deformation by a passing elastic wave, the induced fluid pressure, and therefore the contribution of the fluid to the elastic stress in the medium, will be spatially heterogeneous. Given sufficient time to respond, the fluid can flow so as to equilibrate the distribution of fluid pressure, and the pore fluid eventually becomes "completely relaxed". In terms of elastic behavior, the pores fluids then act collectively as a single effective fluid, and the "heterogeneous" system can be considered to be effectively homogeneous relative to the time scale of Chapter 2. Heterogeneity and Relaxation at the Pore Scale 10 deformation. The response time, or equilibration time, is the key to the definition of heterogeneity and wi l l , in general, increase wi th increasing length scale of the heterogeneity that ex-ists in the porous medium. For a given wave frequency, a heterogeneous distribution of fluid pressure at one scale might have sufficient time to behave as a relaxed or effectively homogeneous system, while heterogeneity at a larger scale wi l l behave as "unrelaxed". Similarly, when the wave frequency is varied, the behavior of a system can change be-tween "homogeneous" and "heterogeneous". For example, a distribution of pore fluids that is homogeneous at the scale of meters, but is heterogeneous at the scale of microme-ters, exhibits effects of heterogeneity when the wave frequency becomes sufficiently high that the micrometer-scale fluid pressure heterogeneity becomes unrelaxed. This link to wave frequency has led to the use of the term "low frequency" to refer to the relaxed or effectively homogeneous case, and "high frequency" to refer to the unrelaxed or het-erogeneous case. Clearly for any given wave frequency there is a crit ical length scale above which the system can be described as heterogeneous (unrelaxed) and below which the system can be described as homogeneous (relaxed). Since most porous media are inherently heterogeneous at some sufficiently small length scale, when they are described as being spatially homo/heterogeneous it must be wi th reference to a particular length scale and frequency of deformation. A substantial amount of theoretical work has been done in predictive modeling of elastic wave velocities in porous media. One class of models, the poroelasticity theories (Gassmann, 1951; Biot , 1956a; Biot , 1956b), has successfully treated the case where the distribution of induced fluid pressure is homogeneous. This is the case either when the elastic response of the pore fluids is truly homogeneous, or the wave frequency is low enough that there is complete communication of fluid pressures throughout the pore space. Thus these models are often applied in field seismic measurements, where wave Chapter 2. Heterogeneity and Relaxation at the Pore Scale 11 frequencies can be sufficiently low that the assumption of homogeneity is valid (Mur-phy, 1985). However, the assumptions required by poroelasticity theories make these models inapplicable at higher frequencies, where incomplete fluid pressure relaxation has significant effects on elastic behavior. A second class of models, the inclusion-based effective medium theories (IBEMT's) (Kuster and Toksoz, 1974a; Berryman, 1980a; Berryman, 1980b; le Ravalec et ai, 1996), have been used to model elastic wave velocities at higher frequencies and for heteroge-neous fluid distributions (e.g., Endres and Knight (1989)). In the derivation of these models, the pore space is treated as being composed of hydraulically isolated fluid-filled inclusions. These models are therefore applicable for describing laboratory measurements of wave velocities at ultrasonic frequencies (Murphy, 1985), when pore-scale communi-cation of fluid pressure can be considered to be negligible, but are inapplicable at lower (e.g., seismic) frequencies where hydraulic connectivity may act to communicate fluid pressure at the pore scale. The inconsistency between the poroelastic and the IBEMT approaches raises the difficulty that neither model can be used to directly compare laboratory measurements made at ultrasonic frequencies with measurements made at field seismic surveys. We also do not have a model with which to predict behavior at seismic scales when the elastic response of the pore fluids is heterogeneous. A number of authors (Dvorkin and Nur, 1993; Mavko and Jizba, 1991) have attempted to extend the poroelastic approach to higher frequencies. However, neither of these models permits explicit treatment of fluid distribution heterogeneity, and neither approach is consistent with both the poroelastic and inclusion-based formulations. Endres and Knight (1997) have incorporated fluid pressure equilibration into the inclusion-based approach, but their model is restricted to the completely relaxed and unrelaxed limits, and therefore cannot be used to compute velocity estimates at intermediate frequencies. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 12 The critical physical mechanism which distinguishes the poroelastic and inclusion-based approaches is that of fluid pressure communication, or relaxation. As discussed above, compression due to a passing elastic wave in a region of heterogeneity (e.g., of the pore fluid distribution, or of the pore space geometry, or both) induces a correspondingly heterogeneous distribution of fluid pressure. At sufficiently low frequencies, the fluid pressure becomes "relaxed", and the poroelasticity approach (e.g., Gassmann's (1951) formula) is applicable. At higher frequencies, however, there is insufficient time for enough flow to occur to relax the induced fluid pressures. This unrelaxed fluid pressure heterogeneity causes a rock sample to be stiffer and therefore to exhibit higher elastic wave velocities, since the higher induced fluid pressures contribute to a greater overall stress in the sample. In the limit where the frequency is sufficiently high that no flow develops and fluid pressures are completely unrelaxed, the assumption of the inclusion-based approach is valid, and an IBEMT is applicable. From the above considerations, it is clear that for a model to be consistent with both poroelastic and inclusion-based models, it must treat the mechanisms of fluid pressure communication and relaxation. It is also clear that the effects of a heterogeneous elas-tic response of the pore fluids must be treated as a frequency-dependent phenomenon. However, in most theoretical work so far (e.g., Endres and Knight (1989) and le Ravalec et al. (1996)) the effects of heterogeneity have been treated only in the unrelaxed limit. In this chapter we attempt to model the effects of fluid pressure heterogeneity and re-laxation at the pore scale only, focusing on the particular relationship between the elastic P-wave velocity and the level of saturation of the wetting fluid (typically water). We in-corporate fluid pressure communication into the inclusion-based modeling framework by allowing Poiseuille flow between inclusions to occur in response to induced fluid pressure gradients. Having assumed this particular relaxation mechanism, we derive expressions for the frequency-dependent effective dynamic elastic moduli for the inclusions. We then Chapter 2. Heterogeneity and Relaxation at the Pore Scale 13 fluid-filled pore solid matrix O o o o o o o Figure 2.1: A sample of a porous solid as conceptualized in the inclusion-based effective medium theory approach. use these effective moduli in an IBEMT computation of the elastic moduli and wave velocities of model porous media. Our results are consistent both with previous unre-laxed (i.e., high frequency) inclusion-based results and with the low-frequency Gassmann theory. 2.2 Introduction to Inclusion-Based Effective Medium Theory In the inclusion-based approach to modeling elastic wave velocities in porous media, the porous medium under consideration is conceptualized as a solid matrix in which fluid-filled inclusions with idealized geometry (usually spherical or ellipsoidal) are embedded. This idealized porous medium is represented schematically in Figure 2.1. Assuming that the size of inclusions is much less than the elastic wavelength, we can neglect scattering and treat the medium as an effectively homogeneous material for which effective elastic moduli can be defined (this is the effective medium assumption). At shorter wavelengths, the medium does not appear homogeneous to a passing elastic wave, and velocities be-come dependent on the details of scattering by inclusions (e.g., fast-path effects). A number of theories have been developed for estimating the elastic moduli of the Chapter 2. Heterogeneity and Relaxation at the Pore Scale 14 type of porous medium represented in Figure 2.1. Notable among these are the perturba-tion results of Eshelby (1957) and Walpole (1971), the Kuster-Toksoz (1974a) model, the self-consistent theories developed by Wu (1966), O'Connell and Budiansky (1974) and Budiansky and O'Connell (1976), the fully self-consistent theory due to Berryman (1980a; 1980b), and the differential self-consistent (DSC) theories developed by Henyey and Pom-phrey (1982) and le Ravalec and Gueguen (1996). A recurring problem with a number of these models is that they are valid only for low inclusion density (i.e., low porosity). The Berryman model and the DSC theories have attempted to counter this difficulty, and are applicable at more realistic levels of porosity. Once effective elastic moduli for the medium have been defined, they are used to compute the elastic wave velocities as one would for a truly homogeneous medium having the same elastic properties. For an isotropic medium the elastic P- and S-wave velocities, vp and vs, respectively, are given by the formulas where K and p are the effective bulk and shear moduli of the medium, and p is the average density. Berryman (1980a; 1980b) uses the concept of impedance matching to derive estimates of the effective bulk and shear moduli of the type of medium depicted in Figure 2.1. The main result of these papers is the set of equations (32)-(33) in Berryman (1980b); we summarize these equations below. Indexing the n constituents of the medium by i = 1.. .n (where each constituent is defined by a particular elastic modulus «j, shear (2.1) (2.2) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 15 modulus Hi, volume fraction Cj , and inclusion geometry), Berryman gives n Y/Ci{^i-Ke{f)P*i = 0 (2.3) i=l n Y/Ci(pi-^)Q*i = 0 (2.4) 1=1 -as estimates of the effective bulk and shear moduli Keff and peg, respectively, of the composite. One of the indices (conventionally i — 1) is chosen to represent the solid matrix phase. The terms P*1 and Q*1 are functions of the elastic properties and geometry of each inclusion type. For spherical inclusions, these terms are given by (Berryman, 1980b, Table 1) P* = + 4 ^ e f f (2.5) where F = (^) (9^eff+ 8^eff\ ( 2 ? ) V 6 J V «ef f + 2//eff J ' 1 ' ' For more general ellipsoidal inclusions, these terms become more complicated; we refer the reader to Berryman (1980b). Since it is not our purpose here to specifically investigate the effects of inclusion geometry (such investigations have been carried out by other others, e.g., Endres and Knight (1989)), we restrict our considerations to the simpler case of spherical inclusions and use equations (2.5)-(2.7) throughout this work. Equations (2.3) and (2.4) implicitly define the effective elastic moduli Keg and peg of the porous composite. They generally cannot be solved analytically for K e f f and pes in terms of the other parameters, and an iterative numerical method must be used. Berry-man (1980b) suggests a simple iterative method for solving these equations numerically. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 16 However, we have found this method to be unstable in a number of instances, and have achieved more reliable results using the Newton-Raphson method (Press et ai, 1992). We note that equations (2.3)-(2.4) were derived assuming that the porous medium is isotropic, so that the effective medium description requires only two independent elastic moduli. This requires a statistically isotropic distribution (and orientation in the case of non-spherical geometry) of inclusions throughout the composite. In principle it should be possible to derive analogs of these equations for the case of anisotropic porous me-dia. However, such is not our purpose here, and throughout this work we consider only isotropic media. While a number of other inclusion-based models exist, Berryman's IBEMT offers many advantages over others. It has been shown (Berryman, 1980b) that where rigorous bounds on effective moduli (e.g., the Hashin-Shtrikman (1963) bounds) are known, the solution of equations (2.3) and (2.4) always satisfy these bounds. Other theories, such as the Kuster-Toksoz model satisfy these bounds only for low porosity. We therefore have reason to believe that the Berryman model generally provides more accurate estimates of effective moduli. Furthermore, it is critical for the work presented here that Berryman uses a dynamic formulation, so that equations (2.3) and (2.4) apply in the case where the elastic moduli of the matrix and inclusions are complex numbers. Static models for estimating effective elastic moduli, such as the interaction-energy approach (O'Connell and Budiansky, 1974), do not allow for complex moduli. A number of inclusion-based models, including the Kuster-Toksoz and differential self-consistent theories, are ambiguous in that one of the constituents must be chosen to serve as the matrix—this phase subsequently plays a role different from that of the other phases in the computation of effective moduli. For levels of porosity on the order of 50%, the choice of the matrix phase is in fact arbitrary, and any one of the constituents might equally well be chosen to serve as the matrix. The Berryman model avoids this ambiguity; Chapter 2. Heterogeneity and Relaxation at the Pore Scale 17 each constituent in the porous medium plays the same role in the computations, and consequently equations ( 2 .3 ) and ( 2 .4 ) are symmetric with respect to interchange of indices. For the reasons discussed above, Berryman's inclusion-based effective medium theory presented here is our preferred method for computing estimates of elastic wave velocities. In the following sections we incorporate hydraulic connectivity into the inclusion-based modeling framework, to extend the Berryman model to low and intermediate frequencies. 2.3 Denning Effective Elastic M o d u l i for Inclusions As discussed above, the primary limitation of inclusion-based effective medium theories (IBEMT's) is that each inclusion is hydraulically isolated from all others. This inability to account for hydraulic interaction of inclusions (i.e., "communication", or "relaxation" of pore fluid pressures) restricts the validity of these models to high frequencies, and makes IBEMT's incompatible with poroelasticity theories used to estimate low-frequency velocities. In order to develop a more realistic effective medium theory applicable to a broader range of frequencies, it is necessary to account for the mechanisms of pore pressure communication. While the models of Endres and Knight (1997) and Xu (1998) have extended the inclusion-based approach to the low frequency limit by allowing for complete fluid pressure relaxation, these models are limited to the completely relaxed and unrelaxed limits. It is necessary to include a more complete description of the relaxation mechanism to predict elastic behavior over a broader range of frequencies. Within the inclusion-based modeling framework, the incorporation of pore pressure communication can be achieved by making hydraulic connections between inclusions in the porous medium represented in Figure 2 . 1 . We conceptualize each hydraulic connec-tion as a duct through which fluid can flow in response to pressure gradients between Chapter 2. Heterogeneity and Relaxation at the Pore Scale 18 fluid-filled pore solid matrix duct (a) (b) Figure 2.2: Incorporating fluid pressure communication into an inclusion-based model, (a) A porous medium sample with each inclusion hydraulically isolated; (b) A sample with inclusions connected by ducts through which fluid can flow. inclusions. The diagrams in Figure 2.2 illustrate schematically this addition of hydraulic connectivity to the inclusion-based model. At sufficiently high frequencies, where there is insufficient time for flow between inclu-sions to develop, the sample depicted in Figure 2.2(b) will exhibit the same elastic behav-ior (i.e., effective elastic moduli and velocities) as the sample in (a), where the inclusions are isolated. At lower frequencies, pressure gradients between hydraulically connected inclusions saturated with different fluids will induce flow, leading to equilibration of fluid pressures at sufficiently low frequencies. This "relaxation" of pore pressures will cause inclusions to be more compliant; consequently, sample (b) will be more compliant—and hence exhibit lower elastic wave velocities—at lower-frequencies. We take the novel approach of treating this frequency dependent compliance of inclu-sions by defining an effective complex bulk modulus for each inclusion. Assuming a time dependence of the form elu}t for the pressure p and the dilation 9 of a given inclusion, we can write a relation of the form p = -K*d, (2.8) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 19 where the coefficient K* defines the complex effective bulk modulus of the inclusion. It should be emphasized that K* as defined here is not an inherent property of the inclusion itself, but rather must be defined with reference to the mechanism of hydraulic interaction with neighboring inclusions. The dependence of K* on frequency will be determined by the particular mechanism of pressure relaxation. In Section 2.4 and following we turn to the problem of explicitly computing K*(UJ) for particular relaxation mechanisms. However, regardless of the partic-ular mechanisms involved, the required limiting behavior of K* for high and low frequency is known. In the high-frequency limit, where there is insufficient time for flow between inclusions to develop, K* for a given inclusion will just be equal to the bulk modulus of the pore fluid in that inclusion. That is, the inclusions are effectively hydraulically isolated. In the low-frequency limit, we expect fluid pressures to equilibrate throughout the pore space. In this equilibrated state the pore fluids behave elastically like a single effective fluid of bulk modulus K*, with where n is the number of distinct pore fluids and Si denotes the fraction of the total fluid volume occupied by the zth fluid (with X^Li & = -0-Having defined an effective bulk modulus K* for each inclusion, the critical conceptual step required to incorporate fluid pressure communication into an inclusion-based model is to replace the network of interconnected inclusions represented in Figure 2.2(b), with a set of elastically equivalent, isolated inclusions with individually defined bulk moduli K* . This step is represented schematically in Figure 2.3. The essential point is that the two samples depicted in Figure 2.3 are elastically equivalent, since the elastic behavior of the individual inclusions in sample (a) is identical to that of those in (b). However, whereas in sample (a) the effect of fluid pressure communication is represented explicitly (2.9) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 20 inclusions have complex bulk moduli K* (a) (b) Figure 2.3: Schematic illustration of the replacement of hydraulically connected inclusions with equivalent inclusions having effective bulk moduli. by the network of ducts connecting the inclusions, in (b) the inclusions are isolated and the effect of fluid flow is encapsulated in the constitutive elastic relations (i.e., equation (2.8)) for the inclusions. This has the advantage that once the effective moduli K* of the inclusions have can be computed, one can use an inclusion-based effective medium theory, applied to the conceptual sample depicted in Figure 2.3, to estimate its effective elastic moduli, and hence the elastic wave velocities in the composite.. In the following sections we develop a number of models of the pore-to-pore relaxation mechanism and derive explicit expressions for the effective inclusion moduli. We then use these models, with Berryman's effective medium theory, to investigate the saturation-and frequency-dependence of elastic wave velocities in a model porous medium. 2.4 A Simple Pore-to-Pore Relaxation Model As an explicit example of a pore scale model of fluid pressure communication, we consider the simplest case of just two fluid-filled inclusions connected by a narrow duct (e.g., a crack between solid grains). Fluid can flow through the duct in response to a pressure Chapter 2. Heterogeneity and Relaxation at the Pore Scale 21 volume V bulk modulus K. L Figure 2.4: Schematic illustration of an idealized pore-to-pore relaxation model: two fluid-filled inclusions connected by a narrow duct. gradient, thus communicating fluid pressure between the inclusions. The model is rep-resented schematically in Figure 2.4. We derive expressions for the effective bulk and shear moduli of the inclusions in terms of the forcing frequency and the properties of the pore fluids, and use these results to estimate elastic wave velocities in a model porous medium, allowing for fluid pressure communication between inclusions. 2.4.1 M o d e l Development Effective Bulk M o d u l i of Inclusions For simplicity, we assume for now that both inclusions have the same volume V, and that each inclusion contains a single fluid. Labeling the inclusions with the indices 1, 2, we take the bulk modulus of the fluid in each inclusion to be Ki and K 2 , respectively. Assuming Poiseuille flow through the duct, the volumetric flow rate q from inclusion 1 to inclusion 2 is given by where the coefficient 7 is a function of the geometry of the duct and the properties of the fluid in the duct, and pi, p2 are the respective fluid pressures in the inclusions. For a cylindrical duct of radius r and length L, containing a Newtonian fluid with viscosity q = 7(Pi - P2) (2.10) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 22 77, 7 is given by (White, 1986, p. 305) 4 irr ~<=WL- (2'N) We assume that the duct radius is sufficiently small that there is enough time for the flow to become completely developed within a wave cycle. Bio t (1956b) derives the more general frequency-dependent version of equation (2.11) in the case where this assumption does not hold. Differential equations for the evolution of the pressures in the inclusions in terms of q and the volumetric dilatations 9\ and 92 of the inclusions, are given by ^ = -Kl(*h + l) ( 2 .12) dt \dt V J v ; dp2 (d82 q\ . . -dT = -K2{-oW-v)> (2-13) and follow immediately from the constitutive elastic relation for the fluids. We assume that the dilatations are equal, so 9\ — 92 = 9 (i.e., the inclusions are deformed under isos-train conditions—this assumption wi l l be justified later by the analysis in Section 2.4.3). Assuming a time dependence of the form eiut for the quantities pi, p2 and 9, equations (2.12)-(2.13) become iujPl = - K l (iu9 + (2.14) - K 2 ( i u 9 - ^ , (2.15) which, upon substitution of q from equation (2.10) can be solved to yield an expression of the form P i = —K*9, (2.16) where K* is given by Chapter 2. Heterogeneity and Relaxation at the Pore Scale 23 and (2.18) U0 = — ( « ! + K 2 ) . (2.19) Equation (2.17) gives the expression for the effective bulk modulus, K * , of inclusion 1. The limiting behavior of K* for high and low frequencies is as expected from the dis-cussion of Section 2.3. In the high-frequency limit, expression (2.17) reduces to «* = /«!. That is, the addition of hydraulic connectivity has no effect on the elastic behavior of inclusion 1, and it behaves as if it were hydraulically isolated from the other inclusion. This is consistent with the usual high-frequency assumption used in inclusion-based mod-eling. In the low-frequency limit, expression (2.17) reduces to K* = «, which is identical to equation (2.9). We investigate graphically the details of the frequency dependence of K* as defined in equation (2.17). Figure 2.5 presents graphs of the real and imaginary parts of ~ vs. ^ for the particular case where we take ^ = 0.1 (i.e., inclusion 2 is filled with a fluid that is much more compressible than that in inclusion 1). The behavior exhibited in Figure 2.5 is typical for a viscous relaxation mechanism (Bourbie et al, 1987, pp. 117-123), with a distinct transition from low-frequency to high-frequency behavior (evident in the graph of the real part of «*) occurring about the characteristic frequency u0. There is also a characteristic peak in viscous dissipation (represented by the imaginary part of K*) near u!Q, corresponding to a maximum near 45° of the phase angle by which the pressure in the pore leads the dilatation. Figure 2.6 shows a graph of this phase angle vs. tu/u>0 for the particular case considered above. It is clear from this figure that the pressure and dilatation are in phase only in the low- and high-frequency limits; otherwise the dilatation lags behind the pressure by a phase angle between 0° and 45°. A critical parameter in the effective pore fluid bulk modulus defined by equation Chapter 2. Heterogeneity and Relaxation at the Pore Scale 24 Figure 2.6: Phase angle 8 (by which pressure leads dilatation) vs. U/UQ , for the pore-to-pore relaxation model. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 25 (2.17) is the inter-pore distance L which enters through the "conductivity coefficient" 7. As adjacent inclusions become further apart and L increases, the rate of flow between them decreases, thus limiting the degree to which fluid pressure is communicated between them. In the particular case described above, at a given frequency an increase in L would cause an increase in the effective bulk modulus of inclusion 1. Clearly, before equation (2.17) can be of practical use in defining effective inclusion moduli for an effective medium computation, we require an estimate of the parameter L for our inclusion-based medium. If inclusions are placed in a regular cubic array with number density n (inclusions per unit volume), the edge length of each unit cell, or equivalently the distance between the centers of inclusions, is given by n - 1 / 3 . For inclusions of equal volume V, the number density n is conveniently specified in terms of the porosity 0 by the relation n = 0/V. While L n - 1 / / 3 provides an estimate of the average pore-to-pore distance, in a more realistic porous medium with inclusions distributed randomly throughout the matrix, L would have a different value for each inclusion. It is therefore more realistic to treat L, and by extension the effective bulk moduli of inclusions, as random variables. If the probability distribution of L is known, the distribution of effective bulk moduli K* of inclusions follows immediately from the relationship between L and K* given by equations (2.17) and (2.11). The distribution of moduli K* is incorporated naturally into Berryman's IBEMT via the terms Q in equations (2.3) and (2.4). In Appendix A we derive the probability distribution for the distance r from the center of a given inclusion to the center of its nearest neighbor. This distribution is given by the probability density function P(r) defined by P(r) = ±Tmr2e-^nr\ (2.20) such that P(r)dr specifies the probability that the center-to-center distance is between r and r + dr. The distribution of duct lengths L is obtained from equation (2.20) by Chapter 2. Heterogeneity and Relaxation at the Pore Scale 26 subtracting from r the length corresponding to the sizes of the inclusions. With equations (2.17) and (2.20), along with estimates of the other physical parame-ters of the porous medium, we can now compute effective bulk moduli for the inclusions and then apply an inclusion-based effective medium theory to compute effective elastic properties for our model porous medium. We now proceed to modeling the saturation-and frequency-dependence of elastic wave velocities for a model sandstone saturated with water and air. We assume a saturation distribution such that each inclusion is filled with either water or air, and vary saturation by increasing the number of water-filled inclu-sions. We assume for now that for a given water-filled inclusion, the relaxation mechanism is dominated by fluid pressure communication with the nearest gas-filled inclusion. At low levels of saturation, gas-filled pores will be abundantly available. Most water-filled pores will be hydraulically connected to a gas-filled pore that is near by, so that the relaxation occurs relatively quickly. As the level of saturation is increased, gas-filled inclusions become less abundant and are typically further away from a given water-filled inclusion, and the relaxation process slows down. At a given frequency, this will result in less relaxed (i.e., elastically stiller) inclusions and hence an increase in their effective bulk moduli. This in turn causes an increase in the elastic stiffness of the composite. The change in abundance of gas-filled inclusions with changing saturation, and its effect on the state of relaxation of inclusions, is accounted for heuristically through the number density n9 of gas-filled inclusions. If the water saturation is Sw, then the number density of gas-filled pores is simply ng = (1 — Sw)n, where n is the number density of inclusions, as before. The number density ng is then used in equation (2.20) to compute distances to gas-filled inclusions. The effect of saturation on the distance from a given liquid-filled inclusion to its nearest gas-filled neighbor is illustrated in Figure 2.7 with a graph of the probability Chapter 2. Heterogeneity and Relaxation at the Pore Scale 27 2 0.5 h 1.5h 0 1h 0 0.2 0.4 0.6 0.8 1.2 1.4 1.6 1.8 2 rn 1/3 Figure 2.7: Probability distribution for the distance r from the center of a given wa-ter-filled inclusion to the center of its nearest gas-filled neighbor, for various levels of water saturation Sw. distribution (2.20) for various values of n9 corresponding to different levels of saturation. This figure clearly illustrates the expected tendency of water-filled inclusions to be further away from their nearest gas-filled neighbor as saturation is increased. The distribution also becomes more spread out with increasing saturation, while it is more narrowly peaked for lower levels of saturation. This suggests that accounting for the random distribution of L (rather than using one estimate for every pair of inclusions) becomes more important at higher levels of saturation. Effective Shear M o d u l i of Inclusions A fluid cannot, by definition, support a static shear stress, and fluids therefore have static shear moduli of zero. However, a viscous fluid can support a dynamic shear stress, and for sinusoidal forcing in time this dynamic shear stress gives rise to a non-zero effective complex shear modulus. To estimate the effective shear modulus of a pore fluid, we consider an idealized model consisting of two infinite parallel plates at y — 0 and y = L, between which is a fluid of viscosity n and density p (see Figure 2.8). The plate at y = L is held stationary while the plate at y = 0 is forced sinusoidally to that its displacement Chapter 2. Heterogeneity and Relaxation at the Pore Scale 28 y = L y = 0 Figure 2 .8 : Shear-induced flow between two plates. x(t) perpendicular to the y-axis is given by x(i) = x0eiu}t. The velocity of the fluid in the x-direction is given by u(y,t). The equation of motion for the fluid in the ^-direction is du d2u (2 .21 ) Assuming a velocity of the form u(y,t) — uQ(y)elujt and defining the kinematic viscosity v = n/p, equation ( 2 . 2 1 ) yields the following equation for the velocity profile u0(y), d2u0 u dy* = I - U Q . V (2 .22) The boundary conditions which u0(y) must satisfy are the no-slip conditions UQ(L) = 0 and ito(O) — U = iuxQ. The general solution of equation ( 2 . 2 2 ) is My) = c^-v^y + c2eV^y, (2 .23 ) for arbitrary constants c i , c2. Imposing the boundary conditions and solving for these constants yields the velocity profile sinh (v^R(l - D) where the dimensionless parameter R is given by (2 .24) R=J-L. v (2 .25 ) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 29 The required force r(t) per unit area on the plate at y — 0 is given by du(y,t)\ r(t) = - 7 7 -dy y=o (2.26) r0e where . . . . ViR T 0 = IUIT] fXo\ y/XHtanh We now imagine that the fluid is replaced by a homogeneous solid having a shear modulus such that when a sinusoidal stress of amplitude r 0 is applied to the boundary at y = 0, the same 'displacement amplitude x0 is induced. The shear modulus of this effective material defines the effective complex shear modulus of the fluid. The shear strain e(£) in this effective material is given by x{t) e(t) = L ' (2.28) Defining the effective dynamic bulk modulus p* by the relation r = p*e, (2.29) it follows from equation (2.27) that [i* is given by T0L V = x0 V~iR = luin-(2.30) tanh \/iR By equation (2.29) we mean that under the specified sinusoidal shear stress, the displacement of the plate at y = 0 is the same as if the fluid were replaced by a solid material having shear modulus p* given by equation (2.30). That is, the fluid could be Chapter 2. Heterogeneity and Relaxation at the Pore Scale 3 0 replaced by this effective medium without altering the elastic behavior observed at the material boundaries. The dimensionless parameter R is a measure of the degree to which the momentum transferred from the plate at y — 0 to the fluid diffuses toward the other plate. In the limit —> 0, the momentum diffusion is complete, and the velocity profile between the plates becomes linear. The velocity shear becomes constant across the fluid layer, and p* —> iun. This is the expression used by Kuster and Toksoz (1974b) and Berryman (1980a; 1980b) for the dynamic shear modulus of pore fluids. For R 3> 1, the velocity shear is localized in a boundary layer at the plate at y = 0, and for increasing R the magnitude of p* grows without bound. We see that there are two physical processes that contribute to "stiffening" of the fluid layer with increasing frequency: (i) increasing u> increases the overall magnitude of the velocity shear, and (ii) increasing u> increases R, which concentrates the velocity shear near the plate at y = 0, giving rise to a greater shear stress at the plate. The expression (2.30) is obviously specific to the geometry of the idealized model being considered. A fluid in the pore space of a particular porous medium is subject to different (and much more complex) geometrical constraints, and the expression for the effective shear modulus will be different. Nevertheless, we regard expression (2.30) as a reasonably accurate estimate for the effective shear modulus of a pore fluid, where the parameter L is interpreted as representing a characteristic length scale of the pore space (e.g., the pore diameter). While the preceding development serves to demonstrate how frequency-dependent effective shear moduli can be defined for pore fluids, in practice we have found that the effects on the elastic behavior of model porous composites are negligible. This can be attributed to the fact that the magnitude of the effective shear modulus as defined in equation (2.30) is always very small relative to the shear modulus of the solid matrix, Chapter 2. Heterogeneity and Relaxation at the Pore Scale 31 as well as to the fact that the contribution of viscous shear to dissipation is dominated by losses due to the fluid pressure relaxation mechanism. Thus, in the modeling results that follow, in order to focus on the effects of fluid pressure relaxation we will entirely neglect the effects of non-zero effective shear moduli of fluids discussed here. 2.4.2 Model ing Results and Discussion Using the model developed in Section 2.4.1, we compute effective elastic moduli and wave velocities over a range of frequencies and levels of saturation for a "model sandstone". Our model sandstone is composed of water- and air-filled spherical inclusions in a solid quartz matrix, and is similar to that considered by Endres and Knight (1989), except that where they have used a distribution of inclusion shapes, for illustrative simplicity we consider only spherical inclusions. Pore-to-pore fluid pressure relaxation is permitted between inclusions filled with different fluids, and equation (2.17) is used to compute the resulting effective bulk moduli of the inclusions. Berryman's (1980a; 1980b) self-consistent effective medium theory is then used to compute estimates of the effective elastic moduli of the composite. Data describing the physical properties of the constituent materials are given in Table 2.1; other parameters used in the model are given in Table 2.2. The results of the self-consistent estimates of effective moduli for the model sandstone are presented in Figure 2.9 as graphs of the real parts of the effective bulk and shear moduli vs. water saturation, at various frequencies. A different curve is plotted for each of eight logarithmically sampled frequencies between 40 Hz and 150 kHz. As expected, at every level of saturation except Sw = 0 and Sw = 1, the high-frequency moduli of the composite are greater than the low-frequency moduli, owing to the greater restriction of pore pressure relaxation. The effective moduli are frequency-independent at 0 and 100% saturation because the pore fluid is homogeneous (either all-air or all-water) and Chapter 2. Heterogeneity and Relaxation at the Pore Scale 32 Matrix Brine Air bulk modulus, K [Pa] shear modulus, p [Pa] density, p [kg/m3] 2.14 x 101 0 3.05 x 109 2.42 x 103 2.10 x 109 0 1.08 x 103 1.38 x 106 0 1.17 x IO1 Table 2.1: Physical properties of the constituents of a model sandstone. (From Endres & Knight (1989)). Parameter Numerical Value porosity <f) 0.31 inclusion radius r 100 inclusion volume V = |7rr 3 0.0042 mm 3 inclusion number density n = <p/V 74 m m - 3 water viscosity 77 10~3 Pa-s duct radius a 2.0 pm Table 2.2: Parameters used to quantify the pore space and hydraulic properties of a model sandstone. Figure 2.9: Computed effective elastic moduli vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 3 3 therefore no pressure gradients arise which would cause flow and frequency-dependent behavior. The imaginary components of the effective elastic moduli are associated with wave attenuation due to viscous losses. While the effective medium computation provides estimates of these imaginary components, we do not graph them here. Instead, we prefer to quantify attenuation in terms of the "inverse quality factor" which has the physical interpretation of being the fraction of wave energy lost to viscous dissipation during a wave period. The elastic P-wave and S-wave velocities and quality factors (vp, vs, Q'1 and Qj1, respectively) for the model sandstone were computed using the equations (after Gueguen and Palciauskas ( 1 9 9 4 , p. 1 7 2 ) ) where «eff and pes are the computed complex-valued effective bulk and shear moduli of the sample, and p is the average density of the sample. The resulting graph of vp vs. Sw is presented in Figure 2 . 1 0 . Again, a different curve is plotted for each of the frequencies sampled in Figure 2 . 9 . The Gassmann curve was computed by first using the inclusion-based model to compute the dry-frame moduli—whereby the poroelastic parameter a can be computed—and then applying equations (1.1)—(1.2) to compute the effective elastic moduli of the fluid-saturated sample. As expected, the high-frequency limit curve plotted in Figure 2 . 1 0 coincides with the result that would have been obtained if the calculation had been performed with-out taking fluid pressure relaxation into effect. As the frequency is decreased, allowing (2 .31 ) (2 .32 ) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 3 4 Figure 2.10: Computed P-wave velocity vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model. For comparison, the results obtained using Gassmann's equation are plotted as the dashed ( ) curve. 0.3 0.25 0.2 0.15 0.1 0.05 0 . . . — ! 1 1 r 1.4 kHz 4.6 kHz^/ / \ . • 15kHz\ 1 1 L 0 0.2 0.4 0.6 S 0.8 w Figure 2.11: Computed P-wave Q 1 vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 35 pore-to-pore relaxation of fluid pressure to occur, the composite becomes more compli-ant and the P-wave velocity correspondingly decreases. Note that in the low-frequency limit, where fluid pressures are equilibrated (at least between pairs of inclusions), the computed vp vs. Sw curve is in close agreement with the low-frequency curve computed using Gassmann's equation. The graph of Q~1 vs. Sw resulting from the same simulation is presented in Figure 2.11. For clarity, curves are plotted for only three different frequencies intermediate between the high- and low-frequency limits. The qualitative features of each curve are similar, with Q~1 going to zero at saturations of zero and one (where there is no flow induced by saturation heterogeneity, and therefore no associated viscous losses). Each curve exhibits a loss peak at a relatively high saturation, with this peak occurring at a lower saturation with increasing frequency. Presumably this peak occurs at that saturation at which the average distance over which a given water-filled inclusion must relax is optimal so as to maximize the viscous loss over a wave cycle. As demonstrated in Figure 2.5, this occurs when the time scale of pressure relaxation is of the same order as the time scale of elastic deformation. It is perhaps significant that the Q~l vs. Sw curves shown in Figure 2.11 are identical in form to the laboratory measurements shown in Figure 6 of Murphy et al. (1986), and are of the same order of magnitude. Figure 2.12 presents graphs of the S-wave velocity vs vs. Sw at various frequencies, for the same simulation. Again, the low-frequency curve nearly coincides with the Gassmann prediction, which is almost indistinguishable in this figure. It is interesting to note that contrary to the assumption usually made (e.g., equation (1.2)), the S-wave velocity at a given saturation does not depend solely on the dry-frame shear modulus. This assumption does appear to be valid at low frequencies (as can be seen by the close agreement with the Gassmann prediction), but the model sandstone appears to have a higher S-wave velocity—resulting from a higher effective shear modulus—at higher frequencies. We Chapter 2. Heterogeneity and Relaxation at the Pore Scale 36 attribute this behavior to the fact that in the inclusion-based formulation, the elastic field induced in the sample by a pure-shear deformation at infinity has dilatational as well as pure-shear components in the vicinity of an inclusion. Thus the behavior of the sample under shear deformation depends to some degree on the bulk moduli—and not just the shear moduli—of the pore fluids. We can see from Figure 2.12 that this dependence is slight, with less than 5% variation in vs with frequency. Figure 2.13 shows the associated graph of Qj1 vs. Sw. As we have neglected the effects of non-zero shear moduli of pore fluids at high frequencies, the S-wave attenuation exhibited here must be due solely to the slight dependence, discussed in the previous paragraph, of the shear behavior on the bulk moduli of the pore fluids. We expect, then, that the S-wave attenuation for our model will be of the same form as the P-wave attenuation, but of much smaller amplitude. This is just the behavior exhibited in Figure 2.13. The transition from low-frequency (i.e., Biot-Gassmann) to high-frequency (i.e., con-ventional inclusion-based modeling) behavior is more clearly illustrated in the frequency domain. We present in Figure 2.14 the same data plotted in Figure 2.10, but as P-wave velocity vs. frequency. At each level of saturation, the vp vs. frequency curve approaches an asymptotic limit at low frequencies, closely corresponding with to the low-frequency limit given by Gassmann's equation. Each curve also approaches an asymptotic limit at high frequencies, corresponding to the high-frequency estimate obtained under the as-sumption of zero pore pressure communication (i.e., the assumption made in conventional inclusion-based effective modeling). Figure 2.15 shows similar graphs of Q~1 vs. frequency, at various levels of saturation. The behavior exhibited in this figure is typical of viscous relaxation phenomena. There is a broad attenuation peak centered about a characteristic frequency corresponding to the characteristic frequency of the pore-to-pore relaxation mechanism, with attenuation Chapter 2. Heterogeneity and Relaxation at the Pore Scale 37 0.87 r s w Figure 2.12: Computed S-Wave velocity vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model. For comparison, the results obtained using Gassmann's equation are plotted as the dashed ( ) curve. Figure 2.13: Computed S-wave Q 1 vs. water saturation at various frequencies for a model sandstone, using the pore-to-pore relaxation model. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 38 2.2 frequency [Hz] Figure 2.14: Computed P-wave velocity vs. frequency at various saturations for a model sandstone, using the pore-to-pore relaxation model. Figure 2.15: Computed Q 1 vs. frequency at various saturations for P-waves in a model sandstone, using the pore-to-pore relaxation model. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 39 tailing off to zero in the low- and high-frequency limits. Because the time-scale of pressure relaxation increases with saturation, due to the fact that the pressures in individual water-filled inclusions must relax over greater distances to their nearest gas-filled inclusion (as discussed in Section 2.4.1), the attenuation peak shifts significantly toward lower frequency with increasing saturation. 2.4.3 Incorporating Pore Compliance Effects While we have developed in Section 2.4.1 a model for the relaxation of heterogeneous fluid pressures arising from a heterogeneous distribution of pore fluids, we have neglected the possibility of heterogeneous fluid pressures arising from spatial variation of pore compliance. We have assumed so far that under deformation of the porous composite, the volumetric dilatation 8 is the same in both inclusions. In reality, the dilatation should be greater in inclusions that have more compliant geometry. In particular, high aspect ratio ellipsoidal inclusions (e.g., thin cracks) are inherently more compliant than spherical inclusions (le Ravalec et ai, 1996)—we expect that the dilatation, and therefore the induced fluid pressure, will be greater in thinner, more compliant inclusions. To derive a more general form of equation (2.17) for the case of inclusions of different geometry, we use the following result derived by Endres and Knight (1997), Kb V Kb J which relates the incremental dilatation d6mc of an inclusion to the applied incremental dilatation d6A at infinity, and the incremental change in fluid pressure pj. The terms Kb, Pb are the bulk and shear moduli of the background medium in which the inclu-sion is embedded (i.e., they are equal to Ke^ and peff of the composite). The function P(Kb, Pb, Kinc, pinc, a) is identical to the term P*1 introduced in Section 2.2. Here a is the aspect ratio of the ellipsoidal inclusion. Berryman (1980b) gives the function P for Chapter 2. Heterogeneity and Relaxation at the Pore Scale 40 a variety of other inclusion geometries. Applying equation (2.33) to inclusions 1 and 2 in Figure 2.4 yields the equations d91 + ^=(d9A + ^l)pl Kb \ Kb J d92 + (^=(d9A + d ^ ) p 2 , Kb \ Kb J (2.34) where we have used the shorthand notation (2.35) Pi = P(Kb, pb, 0, 0, CHi) P2 = P(Kb,pb,0,0,a2). We assume that both inclusions are subject to the same externally applied dilatation 9A, since both are embedded in the same sample. As in Section 2.4.1 (i.e., equations (2.12)-(2.13)), we write dPl = —K\ (de1 + l i p \ ; P 2 ) d t dP2 = - K 2 ( d92 + 1 < J ) \ ~ V l ) d t (2.36) Note that for the sake of generality we have relaxed the previous assumption that Vi — V2. Assuming a time dependence of the form eluJt for the quantities 9A, 9\, 92, P i and p2, equations (2.34) and (2.36) can be solved to yield an expression K = Pi '81 K + l — Ki (2.37) for the effective complex bulk modulus of inclusion 1. Here k and CJ 0 are given by the expressions + P2V2 K p ( H + ^ ) + ^ ( p p ) 7 KXK2 Kb~ K2 + P2K2 V 2 \ (P2 Kb 1 — + — + v2 — «1 «2/ \Pl (2.38) (2.39) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 41 The functional form of K* as defined in equation (2.37) is identical to that of the previous model (i.e., equation (2.17)), so that the discussion of Section 2.4.1 regarding the frequency dependence of K* applies here as well. Indeed, the only differences between equation (2.37) and (2.17) are in the definitions of the characteristic relaxation frequency UQ and the low-frequency effective bulk modulus k. If we impose on the new result the assumptions of the previous section that V\ — V2 = V and Pi = P2 (i.e., both inclusions have the same size and geometry), we find that equations (2.38)-(2.39) reduce to Equation (2.40) is identical to equation (2.18). Equation (2.41) provides the correction to equation (2.19) in the case where (9i ^ 62, and becomes identical to this equation in the limit that K2 <C Equation (2.37) therefore reduces to the result of the previous model (equation (2.17)) in this limit. Hence, the assumption 6i = 92 employed in the previous section is in fact valid in the limit 5 1 « 1, which is generally satisfied in the case of fluid-filled inclusions embedded in a solid matrix, where typically ^ < 10 _ 1 . Equation (2.37) is much more general than equation (2.17), and is more accurate as it obviates the need for assuming that the dilatations of the inclusions are equal. Most importantly, this definition of K* applies to the case where the inclusions have compliances. In the particular case where Ki = K2 (i.e., when the fluids in the inclusions are the same), equation (2.17) yields K* = KI independent of the wave frequency, because no flow between the inclusions is induced. However, in this case equation (2.37) yields K* = Ki only if Pi = P2 (i.e., if both inclusions are equally compliant because they have the same geometry). Otherwise, k is not equal to Ki, and frequency-dependent relaxation of fluid pressures does occur. (2.40) u0 = 1 7 l + f(Pi-l)V («1 + K2). (2.41) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 42 We have derived the result (2.37) in the interest of generality, and to illustrate that relaxation phenomena arising from a distribution of pore geometries can be accounted for with our model. However, in the following sections, we revert to the simpler (though less accurate) result of Section 2.4.1, so that we can more easily consider how the model might be adapted to include the effects of longer-range fluid pressure communication. 2.4.4 Limitations of the Pore-to-Pore Relaxation Mode l The simple model we have developed here incorporates pore pressure relaxation into an inclusion-based model in a consistent way so as to predict elastic wave velocities that reproduce previous high-frequency results, as well as agree closely with the expected ve-locities in the low-frequency limit. While the agreement of the model with the expected behavior in the high- and low-frequency limits is encouraging, the model is probably too simple to be applicable at high saturations and intermediate frequencies. The char-acteristic frequency of relaxation co0, and the relaxation behavior near this frequency, are computed using an heuristic pore-to-pore relaxation model that incorporates fluid pressure communication only to a limited extent. The principal failing of the model is that each inclusion is assumed to be hydrauli-cally connected to only one other inclusion, whereas in a more realistic porous medium we expect inclusions to be multiply connected. Consequently, our present relaxation model underestimates the permeability of the medium and the degree of fluid pressure relaxation that occurs at a given frequency. This limitation of the model causes part of the discrepancy between the computed and the model's low-frequency limit and the Gassmann prediction. In Biot-Gassmann theory (i.e., in the low-frequency limit), fluid pressures are equilibrated throughout the pore space and among all pore fluids. However, in the present model, each inclusion is hydraulically connected to only one other inclusion and is not "aware" of the fluids present in the rest of the pore space. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 43 The ability for the fluid pressure in a given inclusion to relax is also expected to depend on how the fluids are arranged in the pore space. If water-filled inclusions are distributed heterogeneously so that patches of water-filled inclusions exist, inclusions in the patches will have less access to gas-filled inclusions and will therefore be less relaxed than if water-filled inclusions were distributed more homogeneously. Cadoret et al. (1995) and Knight and Nolen-Hoeksema (1990) provide experimental evidence for this stiffening effect of heterogeneous fluid distribution. The simple pore-to-pore relaxation model does not provide a means of incorporating the effects of such global properties as patchy saturation into the relaxation dynamics. Furthermore, in a multiply-connected pore space, as saturation is increased the num-ber of water-filled inclusions relaxing to a single gas-filled inclusion increases, thus limiting the degree to which each water-filled inclusion can relax. This will cause a stiffening effect at high saturations that is not seen in the simple pore-to-pore model, which overestimates the degree of pore pressure relaxation at high saturations by assuming that liquid- and gas-filled inclusions are present in a one-to-one ratio. The low-frequency curve in Fig-ure 2.10 consequently falls below the curve predicted by Gassmann's equation. In the following sections we introduce variations on the simple pore-to-pore relaxation model presented here, in order to overcome some of these shortcomings. 2.5 Incorporating Long-Range Flu id Pressure Communication 2.5.1 M o d e l Development In order for the pore-to-pore relaxation model to agree with the low-frequency results given by Gassmann's equation, each inclusion must have information about the entire pore space so that in the low-frequency limit the effective bulk modulus of the fluid in every inclusion is given by equation (1.3). We achieve this in a somewhat ad hoc manner Chapter 2. Heterogeneity and Relaxation at the Pore Scale 44 as follows. We retain the assumption that, to a first order approximation, the relaxation behav-ior of a given inclusion is dominated by pore pressure communication with its nearest neighboring inclusion. However, we now compute the effect of this relaxation as if the inclusion were connected not just to its nearest neighbor, but to the whole pore space via its nearest neighbor. Conceptually, we replace inclusion 2 in Figure 2.4 with a pore of essentially infinite volume (representing the whole pore space) filled with an effective fluid having bulk modulus «^ f f given by equation (1.3). This has the effect of giving inclusion 1 global information about the distribution of pore fluids, without changing the assumption that the relaxation dynamics are dominated by Poiseuille flow to its nearest neighbor. At low levels of saturation, where the pore space is mostly gas-filled, this change to the model will have little effect, as the effective bulk modulus KJ1 differs negligibly from the bulk modulus of the gas phase. Only at high saturations will the presence of the liquid phase increase K*f. This increase in decreases the induced pressure gradient between pores 1 and 2, limiting the amount of flow and therefore the degree of relaxation of fluid pressure. This mimics the effect of reduced relaxation, discussed above, due to multiple water-filled inclusions relaxing to a single gas-filled inclusion at high saturations. We construct an expression for the effective bulk modulus K* of the fluid in pore 1 by again solving equations (2.10), (2.12) and (2.13), but with pore 2 now filled with a fluid having bulk modulus « / as defined in equation (1.3), and letting the volume of pore 2 go to infinity. Once again, assuming sinusoidal time dependence, solving these equations yields an expression of the form (2.16), where K* is given by KI + i—K? 1 + -*-U>0 (2.42) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 45 with 7«i (2.43) UJ0 = v : Equation (2.42) is of the same form as equation (2.17), with the only difference being that here re = re^ff. Graphs of the real and imaginary parts of re* vs. u therefore have the same form as those in Figure 2.5. The effect of inclusion 2 having infinite volume is that the low-frequency limit of re* as defined in equation (2.42) is now rey-ff, irrespective of the properties of the fluid in pore 1. As this will be true for every inclusion in the porous medium, this change to the model has the global effect of making re^ff the effective bulk modulus for every inclusion in the low-frequency limit, which is the assumption used in the low-frequency Biot-Gassmann theory. We therefore expect this model to be in agreement with the predictions of Gassmann's equation. 2.5.2 Model ing Results and Discussion As before, we have computed effective bulk moduli for inclusions using equation (2.42), used Berryman's IBEMT to compute effective elastic moduli for the composite, and then used equation (2.31) to compute the resulting elastic P-wave velocities. The details of the dependence of the effective moduli on frequency and saturation will not be presented for this model. Rather, we move directly to the consideration of the saturation- and frequency-dependence of the P-wave velocity. In Figure 2.16 we present graphs of vp vs. Sw at various frequencies. Once again, frequencies are sampled at logarithmically spaced intervals to illustrate the spectrum of behavior from the low-frequency limit to the high-frequency limit. Figure 2.16 is qualitatively similar to Figure 2.10, except that the low-frequency curve now agrees much more closely with the expected curve obtained using Gassmann's equation. The pore-to-pore relaxation model now accurately reproduces the Chapter 2. Heterogeneity and Relaxation at the Pore Scale 46 2.2 r •| 51 1 1 i i i 0 0.2 0.4 0.6 0.8 1 S w Figure 2.16: Computed P-wave velocity vs. water saturation at various frequencies for a model sandstone, using the modified pore-to-pore relaxation model. For comparison, the results obtained using Gassmann's equation are plotted as the dashed ( ) curve, which is nearly indistinguishable from the low-frequency limit of this model. Figure 2.17: Computed P-wave Q 1 vs. water saturation at various frequencies for a model sandstone, using the modified pore-to-pore relaxation model. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 47 expected results in both the high- and low-frequency limits. The associated graphs of Q~1 vs. Sw are presented in Figure 2.17. Our modification to the pore-to-pore relaxation model apparently has little effect on the attenuation characteristics, as the curves in this figure are indistinguishable from those of the previous model, given in Figure 2.11. Although this modified pore-to-pore relaxation model has been successful in reproduc-ing the expected results in the appropriate limits, the details of the predicted frequency dependence of elastic waves velocities remains questionable. In particular, the charac-teristic frequency LOQ about which the transition from low- to high-frequency behavior occurs, has probably been over-estimated. This results from the under-estimation of the permeability of the sample, as discussed in Section 2.4.4. In the following section we de-velop a more general relaxation model which accounts for multiple hydraulic connections between inclusions, in order to allow for more realistic sample permeabilities. 2.6 A n Explicit Pore Network Mode l 2.6.1 M o d e l Development In order to model more accurately the frequency- and saturation-dependence of elastic wave velocities, incorporating the increased permeability of a multiply-connected pore-space and accounting for the details of pore fluid distribution, we develop here a relax-ation model for an network of inclusions. Consider a pore space composed of inclusions connected as illustrated in Figure 2.4, but allowing a given inclusion to be hydraulically connected to any number of its neighboring inclusions, forming a hydraulically connected network. We illustrate this schematically in Figure 2.18, which shows a particular inclu-sion (labeled 1) connected to a number of other inclusions (labeled by the indices 2, 3 and 4) by thin ducts through which fluid can flow in response to pressure gradients. The diagram indicates cylindrical ducts of radii r i 2 , r i 3 and r i 4 , and lengths L 1 2 , L 1 3 and L 1 4 . Chapter 2. Heterogeneity and Relaxation at the Pore Scale 48 Figure 2.18: Schematic illustration of the explicit pore network relaxation model: every inclusion is connected to a number of other inclusions by cylindrical ducts through which fluid can flow. However, in principle, ducts of any geometry might be considered. As in the Section 2.4.1, we assume Poiseuille flow through the ducts in response to pressure gradients, such that the volumetric flow rate ^ from inclusion i to inclusion j is given by where the coefficient 7^ is a function of the geometry of the duct, and for the case of a cylindrical duct is given by equation (2.11). The differential equation analogous to equations (2.12) and (2.13), describing the evolution of the fluid pressure pi in inclusion i, is where m, n, . . . are the indices of the inclusions to which inclusion i is hydraulically connected. Again, we have assumed that the volumetric dilatations, 9, imposed on each Qn = Iii iVi ~ Pj) (2.44) (2.45) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 4 9 i n c l u s i o n b y d e f o r m a t i o n o f t h e p o r o u s m e d i u m , a r e e q u a l . A s s u m i n g a t i m e d e p e n d e n c e o f t h e form e l w i f o r t h e q u a n t i t i e s Pi, pj, . . . a n d 8, a n d s u b s t i t u t i n g for q^ f r o m e q u a t i o n ( 2 . 4 4 ) , w e c a n r e w r i t e e q u a t i o n ( 2 . 4 5 ) as lUJPi -iuKi8 - (qim{pi - pm) + qin(Pi - p « ) H ) • (2 .46 ) A n e q u a t i o n o f t h i s f o r m i s c o n s t r u c t e d f o r e v e r y i n c l u s i o n i n t h e p o r o u s m e d i u m s a m p l e ( p o t e n t i a l l y y i e l d i n g a v e r y l a r g e se t o f e q u a t i o n s ) . I n d e x i n g a l l i n c l u s i o n s i n t h e s a m p l e b y t h e n u m b e r s 1 t h r o u g h AT, w e d e f i n e f o r t h e z t h i n c l u s i o n t h e q u a n t i t y 7* b y li l<j<N,j& (2 .47 ) w h e r e w e t a k e 7^ t o b e z e r o i f i n c l u s i o n s i a n d j a r e n o t d i r e c t l y h y d r a u l i c a l l y c o n n e c t e d . D e f i n i n g t h e m a t r i x G b y G = 7 i - 7 i 2 - 7 i 3 • • • - 7 1 J V - 7 2 1 7 2 - 7 2 3 • • • - 7 2 7 V - 7 3 1 - 7 3 2 7 3 • • • - 7 3 ; v (2 .48 ) -7JVI - 7 A T 2 - 7 i V 3 • • • 7iv ( n o t e t h a t G i s a s y m m e t r i c m a t r i x , s i n c e 7^ = 7jj), a n d v e c t o r s o f p o r e p r e s s u r e s a n d p o r e f l u i d b u l k m o d u l i b y Pi P2 a n d K = «i «2 (2 .49 ) Chapter 2. Heterogeneity and Relaxation at the Pore Scale 50 the set of equations (2.46) for all inclusions can then be written in matrix form as iu>p = —iunO — « i / V i 0 0 K2/V2 Gp. (2.50) We can then solve equation (2.50) for the vector of pore pressures p to yield an expression of the form p = -K*9, (2.51) which defines the vector of effective pore fluid moduli, K*, for every inclusion in the sample. This vector is given by / L , A/ , n 1 \ " i / V i 0 0 K2/V2 G (2.52) where I is the identity matrix. In the case of only two inclusions connected by a single duct, equation (2.52) reduces to the previous result, equation (2.17). Equation (2.52) merely generalizes equation (2.17) to the case of a multiply-connected pore space. Fur-thermore, it is immediately apparent that in the high-frequency limit equation (2.52) reduces to the expected results K* = K (i.e., the inclusions behave as is they were hy-draulically disconnected). 2.6.2 Model ing Results and Discussion As before, we compute effective bulk moduli for inclusions using equation (2.52), and use Berryman's IBEMT to compute elastic wave velocities for the composite. We present in Figure 2.19 graphs of vp vs. Sw at logarithmic frequency intervals, for a model sandstone. The pore network was constructed by placing 250 spherical inclusions randomly within Chapter 2. Heterogeneity and Relaxation at the Pore Scale 51 a solid sphere with the appropriate volume to create a sample with the desired porosity, and placing a narrow duct between each inclusion and its three nearest neighboring inclusions. The parameters describing the relevant dimensions and physical properties of the model sandstone are again as given in Tables 2.1 and 2.2. Saturation was varied by starting with all inclusions gas-filled, and filling random gas-filled inclusions with water one at a time. Because the number of inclusions in our conceptual sample is finite, the effects of increases in saturation are discrete, whereas our previous models have depended con-tinuously on saturation. This causes the curves plotted in Figure 2.19 to be jagged at intermediate frequencies. These curves can be made more smooth by increasing the number of inclusions in the sample. However, because of the matrix inversion required in equation (2.52), this rapidly increases the time it takes to compute the effective pore fluid moduli K*. For the illustrative purposes of this particular realization, we found 250 inclusions to be a satisfactory compromise between accuracy and computation time. As expected, the high-frequency vp vs. Sw curve shown in Figure 2.19 reproduces the curve that would have been computed if a conventional inclusion-based effective medium theory without pore pressure communication had been used. We also note that we obtain very close agreement with the Gassmann prediction in the low-frequency limit. However, the main improvement in this model over previous models is that the interconnectedness of the pore space is treated explicitly. Therefore, we expect that the details of the frequency dependence of vp will be more realistic than that predicted with the previous models (due to the more realistic connectivity of the pore space). Furthermore, this model is also suitable for investigating the effects of macroscopic fluid distribution, as explicit control is available over which inclusions are filled with which fluids. The associated graphs of Q~1 vs. Sw are presented in Figure 2.20. Again, these curves are necessarily less smooth than those computed with the previous models. However, the Chapter 2. Heterogeneity and Relaxation at the Pore Scale 52 2.2 r 1.5' 1 1 ' 1 1 0 0.2 0.4 0.6 0.8 1 S w Figure 2.19: Computed P-wave velocity vs. water saturation at various frequencies for a model sandstone, using the explicit pore network relaxation model. For comparison, the results obtained using Gassmann's equation are plotted as the dotted (• • •) curve, which is nearly indistinguishable from the low-frequency limit of this model. Figure 2.20: Computed P-wave Q 1 vs. water saturation at various frequencies for a model sandstone, using the explicit pore network relaxation model. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 53 frequency [Hz] Figure 2.21: Computed P-wave velocity vs. frequency at 80% water saturation for a model sandstone, using the explicit pore network relaxation model. qualitative characteristics are identical to those of the curves presented in Figures 2.11 and 2.17. We illustrate in Figure 2.21 the details of the dependence of P-wave velocity on frequency for this model sandstone. This figure shows a graph of vp vs. frequency for the particular level of saturation Sw — 0.8. Figure 2.21 exhibits qualitative similarities with Figure 2.14, with the expected asymptotic limits for both high and low frequencies, and with a transition for low- to high-frequency occurring about a characteristic relaxation frequency. However, this figure differs substantially from Figure 2.14 in that the transition occurs over a much wider range of frequencies (over more than two decades of frequency rather than just one), and is less symmetric about the central relaxation frequency. Chapter 2. Heterogeneity and Relaxation at the Pore Scale 54 2.7 Discussion and Comparison of the Pore-Scale Relaxation Models In this chapter we have presented a means of incorporating fluid pressure communication (i.e., relaxation) into inclusion-based effective medium theory. We have developed a number of explicit models of the relaxation process, and the results have been shown to be consistent with both the previous high-frequency inclusion-based velocity estimates and the Biot-Gassmann low-frequency estimates. Of the three relaxation models presented here, we consider the explicit pore network model developed in Section 2.6 to be the most accurate and versatile. The previous two models can in fact be viewed as ad hoc attempts to reproduce the combination of local and global flow effects that are accurately represented in the network model. For example, the effects of the random distribution of duct lengths, which were derived analytically for the first two models, are implicitly incorporated in the construction of the inclusion network and need not be incorporated explicitly. Furthermore, the effects of long-range pore pressure communication (i.e., beyond just the nearest neighboring inclusion) and fluid distribution heterogeneity are also automatically accounted for by the network model, while they are approached in an ad hoc way in the model presented in Section 2.4, and not at all in the model of Section 2.5. The pore network model also allows one to explicitly model effects that are beyond the reach of the other models. The explicit construction of the pore network allows one to specify the distribution of fluids in the pore space, the sizes and shapes of inclusions, and the number and radii of ducts (hence the permeability of the medium). While the pore network model provides the means to investigate the effects of many more parameters, as well as correlation structures among these parameters, it has the associated disadvantage that all of these parameters must be supplied explicitly in any particular realization of the model. Unfortunately, too little information about the pore Chapter 2. Heterogeneity and Relaxation at the Pore Scale 55 space of a given porous medium is known in order to construct such a detailed model. At best, one must rely on estimate of the statistical properties of the medium. Another disadvantage of the complexity introduced by the inclusion network is its increased com-putational requirements. The time required to solve the system of equations for the effective moduli of inclusions for a realistic network is perhaps out of proportion with the limited realism achievable with the model. A related issue is that the interpretation of the analytical results is less transparent. While the previous two models provide simple equations for the frequency dependence of effective inclusion moduli, the more general equations provided by the network model are much more complicated and therefore less amenable to intuitive understanding. For these reasons, the modified pore-to-pore relax-ation model may be preferred for a number of applications, as it reproduces the same low-and high-frequency limiting behavior, without the added complexity of the full network model. The most substantial improvement of the network model over the other relaxation models presented here is that it allows for the distribution of pore fluids to be specified explicitly. For reasons already discussed, a heterogeneous distribution of pore fluids is expected to limit the degree of pore pressure relaxation and thereby increase the elastic stiffness of a porous composite. Because the network model allows the inclusion network to be specified explicitly, fluids can be distributed in the network with any desired arrangement and correlation structure (i.e., scale of heterogeneity). The computation of effective inclusion moduli and subsequent application of an inclusion-based effective medium theory allows one to investigate the effect of any partic-ular arrangement on elastic wave velocities and attenuation. However, the computational complexity of the network model makes it unwieldy for many applications, especially if one wishes to investigate the effect of patch sizes that require a network of many thou-sands of inclusions. In the following chapter we take a different approach to modeling the Chapter 2. Heterogeneity and Relaxation at the Pore Scale 56 r e l a x a t i o n o f f l u i d p r e s s u r e s w h e n t h e l e n g t h s c a l e o f s a t u r a t i o n h e t e r o g e n e i t y is o n t h e o r d e r o f m a n y p o r e d i a m e t e r s . A t t h i s s c a l e w e c a n u s e t o a d v a n t a g e a c o n t i n u u m d e s c r i p -t i o n o f t h e p o r o u s m e d i u m a n d i t s s a t u r a t i n g f l u i d s , e x p r e s s e d i n t e r m s o f m a c r o s c o p i c p a r a m e t e r s s u c h as p o r o s i t y a n d p e r m e a b i l i t y . T h i s a v o i d s t h e c o m p l e x i t y o f s o l v i n g fo r t h e fluid p r e s s u r e e v o l u t i o n i n m a n y i n d i v i d u a l i n c l u s i o n s , a n d t h e a n a l y t i c r e s u l t s a r e m o r e t r a n s p a r e n t l y i n t e r p r e t a b l e . Chapter 3 Heterogeneity and Relaxation at the Patch Scale 3.1 Introduction There is substantial evidence from laboratory experiments (Knight and Nolen-Hoeksema, 1990; Cadoret, 1993) and numerical modeling (White, 1975; Endres and Knight, 1989; Gist, 1994; le Ravalec et al., 1996) that elastic wave velocities in fluid-saturated porous media can depend strongly on the spatial distribution of pore fluids. In particular, it has been observed that a porous medium in which the fluids are distributed heterogeneously is typically less compliant—and therefore exhibits higher elastic wave velocities—than the same medium with a more homogeneous distribution of pore fluids. As discussed in greater detail in the introduction to Chapter 2, this effect of heterogeneity depends on the wave frequency. The elastic reinforcement provided by a heterogeneous distribution of pore fluids is due to the fact that the fluid pressure is initially heterogeneous when the medium is deformed by a passing elastic wave. If the wave frequency is low enough, there is sufficient time during a wave cycle for flow to occur, whereby the fluid pressure is equilibrated, and the reinforcing effect of the fluid is minimized. The effect of fluid distribution heterogeneity is therefore most pronounced at high frequencies, where there is insufficient time for fluid pressure relaxation (i.e., equilibration) to occur, and be-comes negligible at sufficiently low frequencies where the induced fluid pressures become completely relaxed. In the current literature (see, for example, Gist (1994)), two distinct scales—the pore 57 Chapter 3. Heterogeneity and Relaxation at the Patch Scale 58 scale and the "patch" scale—of spatial heterogeneity of the fluid distribution are identi-fied. With each of these scales is associated a characteristic time scale of fluid pressure relaxation. Pore scale heterogeneity, considered in detail in Chapter 2, is associated with changes in fluid properties or pore geometry over a length scale less than a few pore diam-eters. The time scale of fluid pressure relaxation at this scale depends on the pore-scale properties of the porous medium, and the properties of the pore fluids. In the present chapter we consider heterogeneity at the "patch" scale, which is as-sociated with changes in fluid and/or matrix properties over a length scale many times greater than a typical pore diameter. (We restrict ourselves, however, to the "effective medium" approximation, in which the length scale of heterogeneity is assumed to be much less than the elastic wavelength, so that scattering from heterogeneities can be ne-glected.) There is significant experimental and theoretical evidence that heterogeneity at this scale can occur in natural geologic systems. Cadoret et al. (1995) provide laboratory evidence that centimeter-scale saturation heterogeneity can arise in limestones during drainage by drying, and that the scale of saturation heterogeneity can depend strongly on the nature of the imbibition/drainage process. Knight et al. (1998) also suggest how patch-scale saturation heterogeneity might arise at much larger (e.g., reservoir) scales when the underlying lithology is heterogeneous. Consider, for example, a liquid/gas-saturated medium in which fully liquid-saturated macroscopic regions (i.e., liquid-saturated "patches") exist adjacent to fully gas-saturated regions. Under deformation by a passing elastic wave, the fluid pressure in the liquid-saturated patches (which is initially higher than in the gas-saturated regions) reinforces the porous matrix and contributes to the elastic stiffness of the sample. However, given a sufficiently long wave period, flow from liquid-saturated patches to the gas-saturated regions will be induced so as to equilibrate the fluid pressure between macroscopic regions, thus limiting the extent to which the liquid reinforces the matrix. The time scale of this Chapter 3. Heterogeneity and Relaxation at the Patch Scale 59 fluid pressure relaxation can be estimated using the analysis of Dvorkin et al. (1994), who derive a diffusion equation governing the spatial distribution of fluid pressure. In the case where the pore fluid is much more compressible than the porous matrix, the fluid pressure diffusivity D is given by D^'if, (3.1) where Kf is the bulk modulus of the pore fluid, K is the intrinsic permeability of the medium, 4> is the porosity, and n is the viscosity of the fluid. For pressure diffusion out of a liquid-saturated patch with radius R, dimensional analysis yields r>2 r « ^ (3-2) as an estimate of the relaxation time scale r. If the time scale of deformation (e.g., the wave period) is much less than r, there is insufficient time for significant relaxation to occur, and the liquid is effectively trapped within the patch. If the time scale of deformation is much greater than r, there is sufficient time for complete relaxation to occur, and the fluid becomes effectively homogeneous. In this limit, Gassmann's (1951) formula (equation (1.1)) can be used to compute the effective elastic moduli—and hence the elastic wave velocities—of the composite. However, for many realistic porous media, and for wave frequencies relevant to geo-physical exploration, the time scale of patch-scale relaxation can be of the same order as the wave period. Consider, for example, a water/gas-saturated Fontainbleau sandstone with Kf = 2.25 GPa, K = 670 mD, </> = 0.136 and n = I O - 3 Pa-s (material data from Lucet (1989)). For a patch radius of 5 cm, equation (3.2) yields r w 2 . 3 x 10 - 4 s, which corresponds to a frequency 1/r of about 4 kHz (i.e., within the sonic range typically relevant in field seismic surveys). Assuming that a time much greater than r is required for complete relaxation to occur, we infer that at seismic frequencies Gassmann's formula Chapter 3. Heterogeneity and Relaxation at the Patch Scale 60 is not applicable for the sandstone considered here if the saturation heterogeneity is of a length scale of 5 cm or greater. Clearly, a more general model is required, in order to account for the effects of heterogeneity at this scale. Le Ravalec et al. (1996) have proposed such a model, using an inclusion-based ap-proach to estimate elastic wave velocities in fluid-saturated porous media with patch-scale heterogeneity, in the unrelaxed (i.e., high frequency) limit. They consider the porous medium to be composed of distinct macroscopic regions having different levels of partial fluid saturation. A differential self-consistent inclusion-based formulation is used to individually estimate effective elastic moduli for each region. The same inclusion-based approach is then used, at the sample scale, to estimate effective elastic moduli for the composite, which is composed of these individual regions in appropriate volu-metric proportions. While their model incorporates effects of both pore- and patch-scale heterogeneity, it does not allow for fluid pressure relaxation, and therefore requires the as-sumption that the wave frequency is sufficiently high (or the length scale of heterogeneity is sufficiently great) that relaxation can be neglected. The two approaches mentioned above (i.e., the poroelastic approach of Gassmann (1951) and the inclusion-based approach of Le Ravalec et al. (1996)) to modeling the relationship between partial fluid saturation and elastic wave velocities, treat opposite extremes of the frequency spectrum, and hence are mutually inconsistent. Each is limited to either the completely relaxed or completely unrelaxed (i.e., high- or low-frequency) regime, and neither provides estimates of elastic and anelastic (e.g., wave attenuation attributed to the patch-scale flow mechanism) behavior in the intermediate regime. We desire a model that is sufficiently general to treat both the high- and low-frequency limits in a consistent way, as well as treat the intermediate ranges of frequency and length scale of heterogeneity. One such model is the "gas pocket" model of White (1975), in which an idealized Chapter 3. Heterogeneity and Relaxation at the Patch Scale 61 fluid distribution is assumed, consisting of a cubic array of spherical gas-saturated pock-ets embedded in an otherwise liquid-saturated medium. Effective elastic moduli for this idealized medium are estimated, taking into account the mechanism of flow induced be-tween the liquid- and gas-saturated regions. Within the limits of effective medium theory (i.e., in which the size of the gas pockets is assumed to be much smaller than the elastic wavelength), this model treats the full frequency range, and reduces to the expected re-sults in the low- and high-frequency limits. However, the White model assumes that the properties of the porous matrix are spatially uniform, and therefore cannot treat the case where the saturation heterogeneity arises because of lithologic heterogeneity as suggested in Knight et al. (1998). Furthermore, the White model allows for only one length scale of heterogeneity, whereas one might realistically wish to consider the case where liquid-or gas-saturated patches of a range sizes exist. The new model we present here contains elements of both the White (1975) and le Ravalec et al. (1996) models. As in le Ravalec et al., we use an inclusion-based approach to compute effective elastic moduli of the composite, which is considered to be composed of distinct regions (i.e., "patches") having different levels of partial saturation. However, rather than use an inclusion-based approach to compute the effective elastic moduli of the individual patches (which immediately restricts the model to the unrelaxed limit), we follow a similar approach to that of White (1975). Using the equations of poroelasticity and Darcy's law, we compute the frequency response of a spherical patch, allowing for fluid flow between regions, and derive expressions for the effective complex bulk moduli of the patch. We then use Berryman's (1980a; 1980b) inclusion-based formulation (rather than the differential self-consistent approach used by le Ravalec et al.) to estimate effec-tive elastic moduli of the composite composed of liquid- and gas-saturated patches. The advantage of our formulation is its ability to model the effects of patch-scale heterogene-ity over a wide range of frequencies and length scales of heterogeneity while remaining Chapter 3. Heterogeneity and Relaxation at the Patch Scale 62 fluid-saturated porous medium sample small volume V solid fluid Figure 3.1: Schematic illustration of a fluid-saturated porous medium, of which we con-sider a particular small volume V. consistent with the expected low- and high-frequency limiting behavior. 3.2 Equations of Pressure Diffusion in a Porous Medium In this section we derive partial differential equations describing the diffusion of fluid pressure in a porous medium. Our results are used in the following section to analyze the relaxation of fluid pressure in a "patch", and to derive an expression for the complex effective bulk modulus of a patch. We consider a sample of fluid-saturated porous medium, and focus on a particular small region within the sample, having volume V and saturated with a single pore fluid (see Figure 3.1). We assume that the porosity within this region is homogeneous, in that when a uniform pressure is applied at the region's boundary the induced fluid pressure is spatially uniform. In particular, this implies either that the pore-scale fluid pressure heterogeneity induced by spatial variation of pore space geometry is negligible, or that the flow generated by this heterogeneity leads to pore-scale fluid pressure equilibrium on a time scale that is much shorter than the time scale of deformation. This is just the assumption used by Gassmann (1951) and many later authors, except that we require Chapter 3. Heterogeneity and Relaxation at the Patch Scale 63 that it hold only for the local volume V rather than for the whole sample. Given this assumption, Gueguen and Palciauskas (1994) (following Nur and Byerlee (1971)) derive the following constitutive relation of poroelasticity, p-apf = -Kd9, (3.3) which relates the induced fluid pressure pf to the pressure p applied at the boundary of the region and the volumetric dilatation 9 of the region. Here Kd is the "dry frame" bulk modulus (i.e., the bulk modulus of the porous material in the absence of pore fluids), and the "poroelastic parameter" a is defined as a = 1 - ^ , (3.4) Ks where KS is the bulk modulus of the mineral solid of which the porous matrix is composed. The simple volumetric average 6 = <j)9f + (1 - (b)9s (3.5) can be used to expressed the overall dilatation 9 of the region in terms of the individual dilatations 9f and 9S of the fluid and solid components, respectively, where 0 is the porosity. Assuming that the fraction of the region's boundary occupied by the fluid phase is identical to the porosity (this assumption is discussed in Biot (1956a)), the same average can be used to relate the total pressure at the sample surface to the pressures in the individual components. That is, P = <ti>f + (1 - 4>)Ps, (3.6) where ps is the pressure in the solid component. Assuming that the solid component is subject to a homogeneous pressure due to the pore fluid, the constitutive elastic relation for the solid component is given by ps = -Ks9s. (3.7) Chapter 3. Heterogeneity and Relaxation at the Patch Scale 64 The constitutive elastic relation for the fluid component can be written as dpf _ d9f KfQ ~dT ~ -Kf-dT ~ W ( } where Kf is the bulk modulus of the fluid. Here Q represents the volume rate of fluid flow out through the region's boundary dV, and can be written as the surface integral Q = I v-nds, (3.9) J&V where v is the volume-averaged volumetric flow rate and n represents the outward unit normal to the sample surface. We now consider the limit of equation (3.8) as the volume of the particular region under consideration goes to zero, yielding the equation describing the evolution of the fluid pressure throughout the sample. Darcy's Law (Bourbie et al, 1987, p. 31), v = - - V P / ) (3.11) V can be used to express the flow velocity v in terms of the fluid pressure gradient V p / , the intrinsic permeability K of the porous medium, and the viscosity rj of the fluid. Substituting expression (3.11) into equation (3.10), we arrive at the final form of the constitutive relation for the fluid, 3.2.1 P r e s s u r e D i f f u s i o n E q u a t i o n for I m p o s e d S t r a i n If we assume that the time evolution of the dilatation field 9 within a sample of porous medium is known, we can solve the system of equations (3.3), (3.5), (3.6), (3.7) and Chapter 3. Heterogeneity and Relaxation at the Patch Scale 65 (3.12) for the time derivatives of the unknown quantities p, ps, pf, 9S and 9f in terms of | | . This yields the equation dpf 2 a 89 — = D e V P f - - F 0 - (3.13) for the evolution of the fluid pressure throughout the sample. The constants Do and Fe are given by De = ^ (3.14) 077 F9 Kf Ks (pK,s If the dilatation in a given region of porous medium is held fixed (i.e., 89/dt = 0), equation (3.13) reduces to the simple diffusion equation 8Pf dt = DeV2Pf. (3.16) We therefore interpret Dg as the fluid pressure diffusivity under conditions of constant dilatation. De is, in fact, identical to the "rock diffusivity" K derived by Dvorkin et al. (1994), and equation (3.13) is the generalization of their diffusion equation for time-varying strain. 3.2.2 Pressure Diffusion Equation for Imposed Stress Alternatively, one can consider the pressure field p to be known throughout the sample, and solve the same system of equations for the time derivatives of the unknown quantities ps, Pf, 9, 9S and 9f in terms of This yields the equation Chapter 3. Heterogeneity and Relaxation at the Patch Scale 66 for the evolution of the fluid pressure. The constants DP and Fp are given by D>=^r (318) F p Kf Ks (j)Kd If the overall pressure in a given region of porous medium is held fixed (i.e., dp/dt — 0), equation (3.17) again reduces to a simple diffusion equation, with diffusivity DP. We therefore interpret DP as the fluid pressure diffusivity under conditions of constant overall pressure. Although the terms DP and Fp differ from DQ and FQ, in practice it is generally true that Kf <C «d < K S , so that FQ f» Fp ~ Kf, and hence D$ ~ Dp. Both diffusivities are then identical to the diffusivity D defined in equation (3.1). It is interesting to note that when V2Pf = 0 (i.e., when the fluid pressure is equili-brated throughout the medium), both of equations (3.13) and (3.17) can be reduced to an expression of the form • which defines an effective bulk modulus, K* , for the sample. In both cases K* is can be shown to be given, after simplification, by — + — -Kf Ks which is just Gassmann's formula (compare with equation (1.1)), and is the expected result in the case of complete fluid pressure relaxation. We can therefore regard equations (3.13) and (3.17) as extensions of Gassmann's formula to the case of non-equilibrated pore fluids, of which equation (3.21) is merely a limiting case. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 67 3.3 Effective Bulk Modulus of a Patch In this section we use the diffusion equation (3.17) derived in Section 3.2 to analyze the frequency response of a fluid-saturated patch undergoing pressure relaxation. We then proceed as in Chapter 2 and seek a relationship between the stress and dilatation of the patch, which we use to derive an expression for the complex effective bulk modulus of the patch. 3.3.1 Solution of the Pressure Diffusion Equation for a Patch We consider an infinite sample of porous medium, of which a spherical patch of radius R (denoted region 1) is saturated with a single fluid having bulk modulus «/(i) . The rest of the sample (denoted region 2) is saturated with a fluid having bulk modulus «/( 2 ) (see Figure 3.2). While the material properties (e.g., porosity, permeability, etc.) are assumed to be spatially uniform within each region, we allow for the possibility that these properties might change across the patch boundary. This is consistent with the expectation that heterogeneity of the fluid distribution can be caused by—and coincides with—heterogeneity of the underlying lithology (Knight et al., 1998). We consider the response of the sample to a uniform sinusoidal pressure p(t) = p0elujt applied at the sample boundary (i.e., at infinity). We assume that the pressure field p is spatially uniform within the sample. While in general this is only approximately true, it is true for the equilibrium elastic field within a spherical patch embedded in a spherical sample under uniform compression. Non-uniform pressures are introduced when the geometrical symmetry of the problem is broken (e.g., when the patch is ellipsoidal rather than spherical). The alternative formulation, which is to impose a spatially uniform dilatation field, is certainly untenable since the dilatation is expected to be discontinuous across the patch boundary. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 68 applied pressure p(t) > region 2 (infinite porous background medium) region 1 (spherical "patch") A A A F i g u r e 3 .2 : S c h e m a t i c i l l u s t r a t i o n o f a m o d e l p o r o u s m e d i u m w i t h a s p h e r i c a l p a t c h o f r a d i u s R ( d e n o t e d as r e g i o n 1) e m b e d d e d i n a n i n f i n i t e b a c k g r o u n d m e d i u m ( d e n o t e d as r e g i o n 2 ) . Chapter 3. Heterogeneity and Relaxation at the Patch Scale 69 Since we regard the pressure field p as prescribed, we use equation (3.17) to describe the evolution of the fluid pressure in the sample, i.e., . Dp(1)V2Pf + r < R (inside the patch), £ P / = J <Md(i) at , 3 2 2 v Dv(2<S72pf + - — — r > R (outside the patch) ^ 02«d(2) at where r is the radial coordinate measured from the center of the patch, and the indices 1 and 2 denote material properties in the regions inside and outside the patch, respectively. The terms D p(i), Dp@), Fp^ and Fp(2) are defined for their respective regions as in equations (3.18)-(3.19). Defining the dimensionless variables f = i , f = T - (3.23) where the time scale r is given by R2 and defining the dimensionless parameters (3.24) _ Qi^p(i) r 0=2^2 n % faFpi?) Kxn2 01 « d ( l ) 02«d(2) £>p(2) 01 ^p(l) • ^ 2 ?7l we can rewrite equation (3.22) in the simpler form (3.26) Exploiting the assumed spherical symmetry, we can write the dimensionless Laplacian operator V 2 as f 2 <9r V df Chapter 3. Heterogeneity and Relaxation at the Patch Scale 70 Assuming a solution of the form Pf(f,i) = pf(r)elbJt, equation (3.26) becomes ' ^ + ! * £ + i i r l P o f < i . \ n \^FT + -^r + ^ r 2 P „ f > i , ^ \ ar 2 r dr J where the dimensionless angular frequency a) is given by u = UT. (3.28) The general solution of equation (3.27) is given by {T i + + %e-Jw f < 1, r (3.29) r r for arbitrary constants c i , C2, C3 and C4. The boundary conditions that Pf(f) must satisfy are dp/ dr = 0 (3.30) lr=0 (ie., no flow through r = 0, as required by symmetry), which implies that c 2 = — C i , and lim |p/(r)| < 00 (3.31) f—+00 (i.e., the solution must be bounded as r —> 00), which implies that c 3 = 0. At f = 1 we must also have continuity of the fluid pressure, i.e., pf(r)=pf(l+), (3.32) and of the fluid flux defined by equation (3.11), i.e., Kx dpf I _ K2dPf % dr f = i - rj2 dr (3.33) f = i + Chapter 3. Heterogeneity and Relaxation at the Patch Scale 71 Applying the conditions (3.32)-(3.33) and solving for the remaining coefficients Ci and c 4 yields the solution ( ( r 2 - Ti)( l + \/iCjQ) sinhVw^r „ 11 H j== j=— , j== r < 1, £viu> cosh vza) + (1 — £ + y/iuSl) sinh \J%CJ r " ^ f ( r i - r 2 ) ( V ^ cosh y/i&- sinh v 7 ^) e v^f2 ( i-r) 1 2 H / = = T = , 7= ^ r > 1. £vza) cosh Vza) + (1 — £ -f- vza)f2) sinh v i £ r (3.34) Figure 3.3 presents graphs of the real part of the solution (3.34) for the fluid pressure inside and outside the patch, at different phases during a wave cycle. The particular values Ti = 0.2, T 2 = 0.02, £ = 0.1 and Q = 0.5 have been chosen for the various dimensionless parameters of the model, corresponding to a patch saturated with a fluid that is ten times less compressible than the fluid saturating the surrounding medium.' In this case we expect the greater fluid pressure induced within the patch to cause diffusion of pressure out of the patch. The wave frequency is chosen so that ur = 103 (i.e., the wave period is about 103 times shorter than the time scale of fluid pressure relaxation in the patch). These graphs illustrate the kind of behavior expected from the discussion in Sec-tion 3.1. Since the time scale of deformation is too short for significant relaxation to occur during a wave cycle, the greater fluid pressure induced inside the patch remains al-most entirely localized within the patch throughout the wave cycle. In other words, there is little fluid pressure communication between the patch and the surrounding medium. Other than a small deviation from a uniform distribution of fluid pressure near the bound-ary of the patch (i.e., near r = 1), the patch and the surrounding medium behave as if they are hydraulically isolated from each other. The patch is almost completely "unre-laxed" , and the induced fluid pressure contributes maximally to reinforcing the porous matrix. Figure 3.4 presents similar graphs, but for the case UJT — 10 - 1 , corresponding to Chapter 3. Heterogeneity and Relaxation at the Patch Scale 72 Figure 3.3: Graphs of the real part of the solution Pf(r)/Po at various phases during a wave cycle, for o»r = 103. Note that r has been scaled so that the patch corresponds to the domain [0,1]. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 73 0.02 -0.02 0.02 -0.02 0.02 € o CL -0.02 0.02 -0.02 cor= 2TI/4 -0.02 Figure 3.4: Graphs of the real part of the solution P/(r)/p 0 at various phases during a wave cycle, for CUT — 10 - 1 . Note that r has been scaled so that the patch corresponds to the domain [0,1]. Note also that the vertical scale is different than that of Figure 3.3. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 74 a wave period that is much longer than the time scale of fluid pressure relaxation in the patch. In this case there is sufficient time during a wave period for much more fluid pressure communication between the patch and the surrounding medium to occur, again as expected from the discussion of Section 3.1. The degree to which fluid pressure is communicated between the two regions is clear from the fact that the fluid pressure inside the patch differs only slightly from that outside the patch. The more compressible fluid saturating the meclium surrounding the patch acts as a low-pressure reservoir into which the fluid pressure induced inside the patch can diffuse. Consequently the fluid pressure inside the patch is much less throughout the wave cycle than in the case LOT = 1 0 3 . The decreased fluid pressure within the patch reduces the degree to which the fluid reinforces the porous matrix, and we therefore expect the sample as a whole to be more compliant in this case. 3.3.2 Definition of the Effective Bulk Modulus of a Patch Using the solution (3.34) for the induced fluid pressure field, along with the poroelastic relation given by equation (3.3), we can write the dilatation field inside the patch as 9(f,i) = e{f)eiC>\ where 0(f) Po - aiPf(f) Po 1 - a i T i + £\fiCb~ cosh \fiCj + (1 — £ + ViuQ) sinh y/iCJ r (1?2 — T i ) ( l + y/ubti) sinh y/i&f )] 3.35) Chapter 3. Heterogeneity and Relaxation at the Patch Scale 75 The total dilatation of the patch, 90, which is just the average of 8(f) over the patch, is then given by 60 = 3 f f26(f) df Jo 10 Po i - a i r\ + -3 (r 2 — Ti)(l + \/i£bQ)(\/iili cosh \fiCj — sinh \fiCS) VJJ £y/iCJ cosh y/iCj + (! — £ + y/iuCt) sinh y/itj (3.36) By relating the total dilatation of the patch to the applied pressure, we can define a complex effective bulk modulus, K*, for the patch, i.e., Po = H*9Q. (3.37) In other words, under sinusoidal deformation with frequency u>, the patch behaves iden-tically to a homogeneous elastic solid with bulk modulus K * . It follows immediately from equation (3.36) that K* is given by n - i l - « i Ti + 3 (T2 — Ti)(l + \JiCjQ)(\/id) cosh — sinh \fiCS) ^Vidj cosh y/iod + (! — £ + y/iuiD) sinh Vioj~ (3.38) After some simplification, and rewritten in terms of the original wave frequency LU, equa-tion (3.38) can be written K « d ( l ) i - a i i r i + zMtr-a-ri)] (3.39) where 3 \ / (1 + V'iu>T£l)(y/'iu>T — tanhy/iur) ILOTJ \ £>yJuJT' + (1 - f + yWrft) tanh V^wr (3.40) Chapter 3. Heterogeneity and Relaxation at the Patch Scale 76 Limiting Forms of the Effective Bulk Modulus We note that in the limit u —> oo (i.e., when there is insufficient time for communication of fluid pressure between the patch and the rest of the sample to occur), f(u) —> 0 and equation (3.39) reduces, after simplification, to K* = KD(L) + 0 I A ^ _ ^ (UJ - C O ) . (3.41) «s(l) This is just Gassmann's formula (compare with equation (1.1)) for the bulk modulus of the patch when the boundary of the patch is taken to be impermeable, and is the expected high frequency limit of the model. The high-frequency limit of our model is therefore identical to the patch-scale heterogeneity model proposed by Knight et al. (1998), where Gassmann's formula is used to estimate the effective bulk modulus of individual patches, under the assumption that the wave frequency is sufficiently high that no fluid pressure communication between adjacent patches occurs. Equation (3.39) can be regarded as the extension of their model to lower frequencies for which fluid pressure communication does occur. We note, however, that our model reduces to Gassmann's formula in the high frequency limit only because we have assumed that the porosity is homogeneous, so that even at high frequencies the induced fluid pressure is spatially uniform within the patch. The presence of pore-scale heterogeneity would act to further stiffen the patch at sufficiently high frequencies. This limitation also applies to the model of Knight et al. (1998). In the limit uo —> 0 (i.e., in the static limit), f(ui) —> 1 and equation (3.39) reduces, after simplification, to ^-^ + ^r^EA {u^0) (3,42) Chapter 3. Heterogeneity and Relaxation at the Patch Scale 77 in the particular case where K S(X) = KS(2), «d(i) = «d(2) and 0i = 02- This is just Gassmann's formula for the bulk modulus of the patch when fluid 1 is replaced by fluid 2. Again, this is the expected low frequency limit of the model. The result is independent of the bulk modulus of the fluid saturating the patch because we have assumed that the region saturated with fluid 2 is infinite in extent. In the case of total fluid pressure communication throughout the sample, the induced fluid pressure (and therefore the contribution of the fluid to the elastic stiffness of the composite) is determined entirely by the bulk modulus of the fluid in region 2. 3.3.3 Accounting for Partial Saturation As noted above, the low- and high-frequency limits of equation (3.39) are in agreement with the known Gassmann predictions for our patch model. However, at low frequencies the model actually underestimates the effective bulk modulus of a liquid-saturated patch in a real finite porous medium sample. This error arises because in our model construction we assume that the background medium surrounding the patch is infinite in extent. Consequently, the overall level of saturation of fluid 1 is effectively zero, and the effective fluid bulk modulus in the case of complete fluid pressure communication is just ft/(2)> independent of In order to use our results to estimate effective patch moduli for partial saturation in realistic finite porous materials, we must refine the model in order to take the effect of finite saturation into account. There are two immediately apparent means of accomplishing this. One is to take the approach of White (1975), and embed the patch at the center of a finite spherical sample, the boundary of which is closed, and the radius R2 of which is determined by the desired level of saturation of fluid 1 . That is, as the overall level of saturation of fluid 1 in the sample approaches 1, R2 approaches R, while as the overall level of saturation approaches 0, R2 approaches infinity. In this case the only change to the model derivation would be Chapter 3. Heterogeneity and Relaxation at the Patch Scale 78 the replacement of the boundary condition (3.31) with the condition 0 (3.43) dPf df \f=R2/R (i.e., the boundary at r = R2 is impermeable). This change serves to simulate the effects of finite saturation, and in the low frequency limit (i.e., in the limit of complete fluid pressure communication) will yield an effective fluid bulk modulus Kef given by 4f = — + — > (3-44) where Si is the overall level of saturation of fluid 1 in the sample. The effective bulk modulus of the patch in the low-frequency limit will then be given by the Gassmann prediction = ""W + <pi t - f r < 3 - 4 5 ) Kf in the particular case where = «d(2)) «s(i) = « s ( 2 ) and 0i = 4>2. This is the expected low-frequency result for partial saturation, and differs from equation (3.42) only in the replacement of Kj( 2 ) with K*f. Unfortunately, when the new boundary condition (3.43) is implemented and the rest of the derivation is carried out, the modified expression for K* becomes impractically complicated. The heuristic model we have proposed hardly justifies the use of such a complex result, and we prefer the following alternative method of incorporating finite saturation effects, which does not increase the complexity of the resulting expression for K*. As described above, we expect the low-frequency limit of our model to yield the modified result (3.45). However, finite saturation effects will not alter the high-frequency limit of the model, where no fluid pressure communication occurs and the elastic behavior Chapter 3. Heterogeneity and Relaxation at the Patch Scale 79 of the patch is independent of the components of region 2. Furthermore, we do not expect partial saturation effects to significantly alter the form of the frequency response, but only to alter the low-frequency limit determined by the parameter r2. We simulate this behavior by considering the patch to be embedded in an infinite background saturated with a fluid having bulk modulus K*f (as given by equation (3.44)), rather than «/ ( 2 ) . That is, we consider the patch to be embedded in an infinite sample having an overall level of saturation S\ of fluid 1 (rather than zero as before), and treat the region outside the patch as if its fluids were distributed homogeneously. This guarantees that the low-frequency limit of the model agrees with equation (3.45). Furthermore, when S\ = 1 (i.e., when the whole sample is completely saturated with the fluid saturating the patch) and the properties of the matrix are constant throughout the sample, we have the expected result that Ti = T 2 , so that equation (3.39) is independent of o>. That is, no patch-scale flow of fluid is induced by the deformation of the sample. The model is easily generalized to the case of n distinct pore fluids having bulk moduli « / ( 2 ) ) - • • iKf{n)- I n this case the generalization of equation (3.44) is given by (Domenico, 1977) K f i=i K / M where Sj is the level of saturation of the ith fluid (defined such that Yli=i & = -0- The final form of our model for the effect complex bulk modulus K* of a spherical patch saturated with fluid 1, in a porous medium having a level of saturation S% of fluid 1, is then summarized as follows, K* = M l ) fa A7\ l - a i t T i + Z M f T a - r O ] 1 ' } IUT — tanh y/iur) f(u) = IUJTJ \ £v/iZ~Jr~+ (1 - £ + Viurri) tanh Viurr Chapter 3. Heterogeneity and Relaxation at the Patch Scale 80 where the various dimensionless parameters are defined by aiFp(i) 0 ^ ( 2 ) n Dp{1) Kxr]2 R2 ^l«d(l) 02«d(2) A>(2) if2??l Dpii) and we have further defined the quantities 1 1 1 a i 01^ 71 ^p(l) Ks(l) 0i«d(l) n _. ^ 2 F p ( 2 ) 1 _ 1 1 a 2 02??2 Fp(2) K f KS(2) 02«d(2) 3.4 Discussion 3.4.1 Frequency Response of a Fluid-Saturated Patch In Figure 3.5 we present graphs of the real and imaginary parts of K*/KJ vs. LUT/(2IT), for the same particular choice of values of the various dimensionless parameters as used in Figures 3.3-3.4, and with c\i = 0.5. The parameter CUT/(2TV) can be interpreted as a measure of the wave frequency relative to the relaxation frequency of the patch, so that when LOT I (2ir) « 1 the time scale of deformation and the time scale of fluid pressure relaxation are of the same order. The graphs exhibit the expected qualitative behavior, with the patch stiffening with increasing frequency. As CUT approaches infinity, Re{«;*} approaches the asymptotic value associated with the upper Gassmann limit given by equation (3.41). As LOT approaches zero, Re{K*} approaches an asymptotic value associ-ated with the lower Gassmann limit given by equation (3.42). The transition between the low- and high-frequency limits occurs near the characteristic relaxation frequency, where LOT/(2TT) ?n 1. The imaginary part of K* (which represents the component of the stress in the patch that is out of phase with the dilatation, and is associated with viscous losses caused by fluid flow) exhibits a peak near the relaxation frequency, and goes to zero in the low- and high-frequency limits. In these limits the pressure and dilatation within Chapter 3. Heterogeneity and Relaxation at the Patch Scale 81 1.25 0.08 y oT 1.1 cr vi°1.15 1.05 1.2 10 - 2 COT/(2TC) on/(27t) Figure 3.5: Graph of the real part of K*/'nd vs. UIT/(2TT) for the same values of the various parameters as used in Figures 3.3-3.4, with ct\ = 0.5. the patch are in phase. We note that, in contrast to the pore-scale relaxation models developed in Chapter 2, the transition from low- to high-frequency behavior is much more gradual, occurring over more than three frequency decades. Significant relaxation of fluid pressure occurs even for u>r/(27r) FU 100, and the patch becomes essentially un-relaxed only for UJT/(2IT) > 103. Note that this was a feature observed when long-range fluid pressure communication was accounted for in the pore-scale relaxation model of Section 2.6. 3.4.2 Effect of Permeability Contrast Across a Patch Boundary The customary estimate for the time scale of pore fluid relaxation in a patch, as given by equation (3.2), assumes that the properties of the porous matrix and fluids are con-stant across the patch boundary. Consequently, only the properties of the fluids and porous matrix within the patch itself appear in equations (3.1)-(3.2). However, many naturally occurring porous media are composed of heterogeneous lithology. In particu-lar, the permeability of porous materials can vary over many orders of magnitude (e.g., from the sub-/xD scale in tight gas sandstones to the scale of Darcys in unconsolidated Chapter 3. Heterogeneity and Relaxation at the Patch Scale 82 sands), and in a porous medium composed of multiple lithologic units (as considered in Knight et al. (1998)) significant spatial heterogeneity of permeability can be expected. Furthermore, it has been suggested (Knight et ai, 1998) that patch-scale saturation heterogeneity might even be caused by the underlying lithologic heterogeneity. Spatial heterogeneity of permeability clearly must have an effect on the time scale of patch-scale fluid pressure relaxation. At one extreme, one can consider a permeable patch surrounded by an impermeable background medium, in which case the patch will be unrelaxed at all wave frequencies and the time scale of relaxation given by expression (3.2) becomes meaningless. This effect of permeability heterogeneity on the time scale of patch-scale fluid pressure relaxation has received little attention in the current literature. In this section we briefly investigate this effect using the patch relaxation model developed in Section 3.3. We consider a water-saturated sandstone patch with permeability K~i, embedded in a gas-saturated background porous medium having a permeability K2 that is different than that of the patch. The remaining parameters required in our patch-scale relaxation model are assigned the particular values 771/772 = 20, Fp^/Fp^) = 400, (j>i/(p2 — 1, I \ = 0.9 and T 2 = 0.002, which are typical for a water/air-saturated sandstone. We have used equation (3.39) to compute the effective bulk modulus, «*, of the patch. The results are presented in Figure 3.6 as graphs of the real part of K* / Kd vs. U)T/(2TT), for various values of the permeability contrast K\jK2. The figure illustrates that the patch is almost completely relaxed when ur/(2n) « 1 in the case where Ki/K2 = 1 (i.e., if there is no permeability contrast). This indicates that, in this case, r is in fact a reasonable estimate of the time scale of relaxation. The same is true for the case where K\jK2 = 102. However, when the permeability contrast is increased to 104, the relaxation of the patch becomes significantly affected by the decreased permeability of the medium in which it is embedded. The relaxation frequency, Chapter 3. Heterogeneity and Relaxation at the Patch Scale 83 C0T/(27C) Figure 3.6: Real part of K* / KD vs. LUT/(2IT), for a water-saturated patch in a typical sandstone that is otherwise saturated with air, and for various values of the permeability contrast Ki/K2 (each curve is labeled with its associated value of K~i/K2). about which the transition from relaxed to unrelaxed behavior occurs, is almost one decade less than when there is no permeability contrast. In this case r underestimates the actual time scale of relaxation by a factor of about 10. When the permeability contrast is further increased to 106 this effect becomes accentuated, and r underestimates the actual time scale of relaxation by a factor of more than 100. As K~i/K2 approaches infinity, the wave frequency required for significant relaxation of the patch to occur approaches zero, at which point relaxation can no longer occur at all. On the other hand, we find that when Ki/K2 < 1 (i.e., when the permeability of the background medium is greater than that of the patch) K* is not significantly altered from the case where there is no permeability contrast, even as K\jK2 approaches zero. We conclude that, at least for the particular medium considered here, equation (3.2) provides an accurate estimate of the time scale of relaxation of a liquid-saturated patch when then the permeability contrast K\jK2 is less than 102. However, care must be Chapter 3. Heterogeneity and Relaxation at the Patch Scale 84 taken when using this expression to estimate the state of relaxation of a patch when the permeability contrast is greater than 102, where the reduced permeability of the background medium begins to limit the relaxation of the patch. 3.4.3 Effect of Patch Size In this section we investigate the effect of patch size on our model for the effective bulk modulus of a patch. We consider the particular example of a spherical water-saturated patch in a Spirit River sandstone that is otherwise saturated with air. The input data corresponding to these materials are summarized in Table 3.1. We use equation (3.39) to estimate the effective bulk modulus of the water-saturated patch for a range of patch sizes and wave frequencies. Figure 3.7 presents the resulting graphs of the real and imaginary parts of the computed values of K*/KD vs. patch size, for wave frequencies of 100 Hz and 50 kHz. This figure illustrates the expected dependence on patch size and frequency. At a given frequency Re{re*} for a water-saturated patch increases with increasing patch size, owing to the decreased degree of fluid pressure relaxation that can occur during a wave period. When the wave frequency is increased (from 100 Hz to 50 kHz in Figure 3.7), a smaller patch size is required for significant relaxation to occur within the time scale of deformation. Figure 3.7 provides an explicit example of the observation made in Chapter 2 that the very term "heterogeneity", in describing the distribution of pore fluids, is a frequency-dependent concept. For a wave frequency of 100 Hz, we find that water-saturated patches with radius less than about 0.2 mm are completely relaxed (i.e., the fluid distribution is effectively homogeneous, and Gassmann's (1951) approach gives a valid estimate of the effective elastic properties of the sample). This is in contrast with the 50 kHz case, for which a patch radius of less than about 10 pm is required for complete fluid pressure re-laxation to occur. In other words, a water/air-saturated sample of Spirit River sandstone Chapter 3. Heterogeneity and Relaxation at the Patch Scale 85 Parameter Numerical Value Porous Matrix (Spirit River Sandstone) KD [GPa] 5.6a pd [Gpa] 12.6a 0 0.0526 K [pD] . \ 1*_ Mineral Solid (Quartz) KS [GPa]. 38c ps [kg/m3] 2630c Water K w a t e r [GPa] 2.2 fato [10~3 Pa-s] 1_ K a i r [MPa] 0.8 77air [IO'3 Pa-s] 0.05 "Estimated from the velocity measurements of Knight and Nolen-Hoeksema (1990) 6From Knight and Nolen-Hoeksema (1990) cFrom le Ravalec et al. (1996). Table 3.1: Physical parameters for Spirit River sandstone saturated with water and air. R [m] R [m] Figure 3.7: Graphs of the real and imaginary parts of K*/Kd vs. patch radius at wave frequencies of 100 Hz and 50 kHz, for a spherical water-saturated patch in a Spirit River sandstone that is otherwise saturated with air. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 86 with fluid distribution heterogeneity of length scale on the order of 0.2 mm is effectively homogeneously saturated at wave frequencies less than 100 Hz, while the same sample behaves heterogeneously if the wave frequency is increased to 50 kHz or greater. At 50 kHz, the sample can be regarded as effectively homogeneous—so that the Gassmann estimate of the effective elastic moduli is valid—only if the length scale of heterogeneity is decreased to less than about 10 pm. These observations can be simply understood in terms of the scaling of the relaxation time scale r (as defined in equation (3.2)) with patch size. Both LO and R appear in equation (3.39) only through the dimensionless term LOT, in which LO and R2 play identical roles. That is, for a given set of matrix and fluid properties, increasing the wave frequency LO has the same effect on the effective bulk modulus of a patch as increasing R2 by the same relative amount. It comes as no surprise, then, that the effect of increasing the length scale of heterogeneity is the same as increasing the wave frequency, and that whether or not a sample behaves as if the fluids were distributed homogeneously depends on the wave frequency as well as the length scale of heterogeneity. Figure 3.7 also provides an explicit illustration of the observation, made by Knight et al. (1998), that if the pore fluids are arranged in sufficiently large patches, the fluids can have a stiffening effect even at low frequencies. For the particular combination of porous medium and fluids considered here, and for a wave frequency of 100 Hz, this stiffening effect occurs as the patch radius is increased above 0.1 mm, and is maximized for a patch radius greater than 1 cm. Thus, the conventional assumption that the pore fluids are completely relaxed in this frequency range—and that Gassmann's (1951) approach can be used to estimate the elastic wave velocities—is not valid even for this modest degree of saturation heterogeneity. Clearly, this could have serious consequences for the interpretation of seismic data. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 87 3.4.4 Mult iple Scales of Heterogeneity In real porous media, instead of a single length scale of heterogeneity, we expect hetero-geneity to exist at a wide range of scales, from the pore scale (on the order of microme-ters) through the patch scale (millimeters to centimeters) and extending to the scale of lithologic variation (on the order of meters or greater). While the models developed in this thesis do not provide a consistent framework for treating the simultaneous effects of multiple scales of heterogeneity (where the treatment is complicated by the fact that successively smaller scales of heterogeneity may themselves be embedded within larger scales of heterogeneity), they do provide insight into the expected qualitative behavior. Consider the following thought experiment in which a number of distinct length scales of heterogeneity are present in a fluid-saturated porous medium. From the our previous discussion it should be clear that there will be a frequency threshold above which the pore fluids can be regarded as completely unrelaxed, i.e., above which there is no sig-nificant flow of the pore fluids during a wave period. (We assume that all of the scales of heterogeneity are much smaller than the elastic wavelength at this frequency, so that scattering is negligible.) Starting at this frequency, if the wave frequency is gradually decreased, eventu-ally the wave period will be sufficiently long for fluid pressure relaxation to occur at the smallest length scale of heterogeneity (e.g., the pore scale, or the smallest patch scale). From our previous discussion it should be clear that it is the smallest scale of heterogeneity—by virtue of the correspondingly short length over which fluid pressures must be communicated—that is the first to relax. The frequency threshold at which this relaxation occurs is predicted by the pore-scale relaxation model of Chapter 2. As the wave frequency is decreased below this threshold, relaxation of fluid pressures at the pore scale weakens the elastic response of the composite, and there is an associated decrease Chapter 3. Heterogeneity and Relaxation at the Patch Scale 88 in the wave velocities. If the wave frequency is decreased further, the wave period eventually becomes suf-ficiently long that fluid pressure at the next greater length scale of heterogeneity (e.g., the patch scale) has sufficient time to relax. Associated with this relaxation is a further decrease in the elastic wave velocities. The frequency threshold below which this relax-ation occurs is predicted by the patch-scale relaxation model of Section 3.3. As the wave frequency is further decreased a cascade of relaxations occur, as the relaxation frequency associated with each successively greater length scale of heterogeneity is traversed. Our patch-scale relaxation model predicts that this relaxation frequency decreases as R~2, where R is a characteristic size of a fluid-saturated patch. When the wave frequency becomes sufficiently low that the fluid pressure is completely relaxed at the greatest length scale of heterogeneity present, Gassmann's (1951) approach to modeling the elastic properties of the composite becomes valid. To explicitly illustrate the thought experiment described above, we consider the par-ticular example of a water/air-saturated porous rock in which the water is present in distinct patches of 0.2 mm and 1 cm. We use the same material properties given in Table 3.1, and assume that the medium is 80% water-saturated, with the water being distributed in equal volumetric proportions between the larger and smaller patches. We use the patch-scale relaxation model to estimate effective elastic moduli for the patches, and use Berryman's (1980a) inclusion-based effective medium theory to compute effective elastic moduli for the composite at various wave frequencies. Figures 3.8-3.9 present the resulting graphs of P-wave velocity and attenuation, re-spectively, vs. wave frequency. Both graphs distinctly show the two successive relaxations that occur as the wave frequency is decreased from the completely unrelaxed limit to the completely relaxed limit. First the smaller patches become relaxed as the frequency is decreased through 104 Hz (where there is a distinct drop in P-wave velocity as well as a Chapter 3. Heterogeneity and Relaxation at the Patch Scale 89 frequency [Hz] Figure 3.8: P-wave velocity vs. wave frequency, computed using the patch-scale relaxation model for a particular low-permeability sandstone having 80% water saturation, with half of the water in patches with radius 0.2 mm and the other half with radius 1 cm. 0.08 0.07 h frequency [Hz] Figure 3.9: P-wave attenuation vs. wave frequency for the same medium considered in Figure 3.8. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 90 peak in attenuation), then the larger patches relax as the frequency is further decreased through about 10 Hz (where there is again a drop in P-wave velocity and a peak in attenuation). 3 . 5 An Application: Modeling Saturation Hysteresis As an application of the patch relaxation model developed in Section 3.3, we consider the effect of saturation hysteresis on elastic P-wave velocity measurements. Laboratory measurements of elastic wave velocities have demonstrated hysteresis in the relationship between P- and S-wave velocities and the level of water saturation in porous rocks (Knight and Nolen-Hoeksema, 1990; Cadoret, 1993). This velocity hysteresis has been associated with hysteresis of the scale of fluid distribution heterogeneity resulting from an imbibition or drainage process. The imbibition process tends to result in a more homogeneous distribution of the pore fluids, while the drainage process tends to result in a more heterogeneous, "patchy" distribution of fluids; this in turn causes fluid pressures to be less relaxed and the elastic wave velocities in the medium to be greater than during imbibition. In this section we use the patch-scale relaxation model developed in Section 3.3 to extend the inclusion-based approach of le Ravalec et al. (1996) to a broader range of frequencies. Because our relaxation model relies on the assumption that the porosity is homogeneous on the scale of an individual patch, we must restrict our consideration to porous media which satisfy this assumption. That is, we cannot account for the effects of additional heterogeneity at the sub-patch scale. In this sense our model is less general than that of le Ravalec et al, which accounts for pore-scale heterogeneity. We consider a liquid/gas-saturated porous medium in which the liquid and gas phases are segregated into 100% liquid-saturated and 100% gas-saturated regions (i.e., patches). Chapter 3. Heterogeneity and Relaxation at the Patch Scale 91 We associate with each region a length scale R representing the patch size or "radius", and use equation (3.47) to estimate an effective complex bulk modulus for the patch. The effective shear modulus of each patch is defined using the usual assumption that the effective shear modulus of a fluid-saturated medium is identical to the dry frame modulus of that medium, i.e., equation (1.2) (although, as discussed in Section 2.4.2, this is only an approximation). With effective elastic moduli defined for all liquid- and gas-saturated regions, we then compute effective bulk and shear moduli for the composite, using Berryman's (1980a; 1980b) inclusion-based effective medium theory, assuming that the patch sizes are all much less than the elastic wavelength. (In le Ravalec et al. this step is carried out using a differential self-consistent approach. The advantage of using Berryman's effective medium theory, beyond its consistency with all known rigorous bounds on effective moduli, is that is allows for the elastic moduli of the constituents of the composite to be complex, which is critical in the implementation of our model.) Although Berryman's effective medium theory allows for the geometry of inclusions to be quite general, at present we restrict ourselves to the case of approximately spherical patches. P-wave velocities and attenuation for the composite are then computed from the effective elastic moduli of the composite with equation (2.31). There appear to be at least the three following mechanisms by which the level of liquid saturation in a sample can be varied: 1. Varying of the level of liquid saturation in existing patches. 2. Varying the number of liquid-saturated patches of a given size. 3. Varying the size of existing liquid-saturated patches. In the model of le Ravalec et al. (1996), mechanisms 2 and 3 are combined into a single mechanism. A single parameter—the volume fraction of the sample occupied by patches Chapter 3. Heterogeneity and Relaxation at the Patch Scale 92 of high liquid saturation—is used to simultaneously vary the size and number of patches, without distinguishing between these mechanisms. This is possible because their model treats only the high-frequency limit, in which only the volume fraction occupied by patches, not the size of individual patches, appears as a parameter in the inclusion-based effective medium theory. However, in our frequency-dependent patch-scale heterogeneity model, we must distinguish between mechanisms 2 and 3, because the effective bulk modulus of a patch depends on its size. Given the above three independent means of varying the level of liquid saturation in a sample of a porous medium, we can expect a wide variety of possible saturation distributions—and therefore a wide variety of functional forms for the dependence of elastic wave velocities on saturation—arising from a given imbibition or drainage process. To simplify matters, we eliminate mechanism 1 by assuming that the sample can be divided into distinct regions of 100% liquid saturation and 100% gas saturation (i.e., liquid-saturated patches and gas-saturated patches). We make this assumption only to more clearly elucidate the effects of the length scale of saturation heterogeneity on P-wave velocities measured during imbibition/drainage experiments, and in practice it is easily relaxed. For the same reasons, we also eliminate mechanism 2 by assuming that during an imbibition-drainage process the number of liquid-saturated patches in a given sample is fixed, while only the sizes of individual patches vary. To model the variation of patch size with the level of liquid saturation in a sample, we define a characteristic length, Ro, representative of the length scale of saturation heterogeneity. We then assume that the typical size, Rw (e.g., the radius), of a given liquid-saturated patch can be expressed in terms of the level of liquid saturation SW as Rw = SWRQ. (3.48) That is, we assume that the size of a given liquid-saturated patch increases linearly with Chapter 3. Heterogeneity and Relaxation at the Patch Scale 93 saturation. The size of a patch is assumed to be zero at Sw — 0, and is assumed to approach R0 as Sw approaches 1. This relationship is arbitrarily chosen as the simplest relationship between Sw and R^. By defining the following linear estimate of liquid saturation, S w = T ^ T T T ' ( 3 - 4 9 ) where Rg is the typical size of a gas-saturated patch, we can substitute for Ryj from equation (3.48), and rearrange equation (3.49) to write Rg = (l- Sw)Ro- (3.50) That is, as Sw approaches zero the size of gas-saturated regions approaches RQ, while as Sw approaches 1 the size of gas-saturated regions approaches 0. Equations (3.48) and (3.50) together provide a physically reasonable parameterization of the variation of the distribution of pore fluids as a function of Sw. The parameterization involves the single parameter R0 which determines the length scale of saturation heterogeneity. Following the discussion of Cadoret et al. (1995), we can characterize a given imbi-bition or drainage process by the length scale of heterogeneity that is associated with it. In our model this length scale is associated with the parameter RQ. In particular, we expect imbibition to result in a more homogeneous distribution of pore fluids, and we therefore model the imbibition process as one for which the associated value of R0 is relatively small. That is, during imbibition the level of liquid saturation is increased via the increasing size of a large number of relatively small liquid-saturated patches. By contrast, we expect the drainage process to result in a more heterogeneous distribution of fluids, and we therefore model the drainage process as one for which the value of Ro is relatively large. That is, during drainage the level of liquid saturation is decreased via the reduction in size of fewer and larger liquid-saturated patches. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 94 We consider the particular example of a sample of Spirit River sandstone saturated with water and air. Again, the relevant physical parameters are given in Table 3.1. Using the inclusion-based model described above, we compute the elastic P-wave velocity in the sample as a function of Sw, for various values of R0. In Figure 3.10 we present the resulting graphs of vp vs. Sw for a wave frequency of 100 Hz, and for logarithmically spaced values of R0 ranging from 79 pm to 1.8 cm (these particular values are chosen only in order to most clearly illustrate the full range of dependence on Ro). As expected, we find that a sample having a greater length scale of saturation het-erogeneity exhibits a higher elastic P-wave velocity at all levels of liquid saturation. This is due to the increasing time required for a patch to relax as the size of the patch is increased. For the Spirit River sandstone, we find that for a length scale of heterogeneity of about 79 pm or less, there is sufficient time for complete fluid pressure relaxation to occur at a wave frequency of 100 Hz. At this scale of heterogeneity the fluid pressure is equilibrated throughout the sample, and in fact the curve corresponding to Ro = 79 pm coincides with the Gassmann prediction for this sample, as would be expected. When the scale of heterogeneity is increased to about 1.8 cm or greater, the opposite is true: there is insufficient time for significant relaxation of fluid pressures to occur. Hence the curve corresponding to RQ = 1.8 cm coincides with the result that would be predicted by the high-frequency model of le Ravalec et al. (1996) (given the additional assumption of homogeneous porosity that we require). Figure 3.11 presents similar graphs, except that a wave frequency of 10 kHz, and logarithmically spaced values of Ro between 7.9 pm and 0.18 cm, were chosen. Because the time scale of relaxation r defined in equation (3.24) scales as R2, the curves in this figure are exactly the same as those in Figure 3.10, since we have increased the wave frequency by two orders of magnitude while decreasing the patch size by one order of magnitude. Complete relaxation of fluid pressures at this wave frequency requires a scale Chapter 3. Heterogeneity and Relaxation at the Patch Scale 95 4 r 3.8 ^ 3.6 JL 3 4 3.2 3 2.8 -1 . 8 / ^ ^ j -^ ^ ^ ^ ^ ^ - ^ 0 1 2 ^ ^ ^ ^ ^ ^^jsmX-^J 0.0079 1 1 i ; i i 0 0.2 0.4 0.6 S 0.8 w Figure 3.10: Elastic P-wave velocity vs. water saturation in Spirit River sandstone (data given in Table 3.1) at 100 Hz, computed using the patch-scale relaxation model, illus-trating the dependence on the length scale of saturation heterogeneity. Each curve is labeled with its corresponding value of Ro, in centimeters. 4 3.8 3.6 J 3.4 3.2 3 2.8 • 0 . 1 8 / 0 ^ J - ^ ^ ^ 0 0 4 6 ^ ^ - ^ 1 X 0 1 2 ^^^^ ^ o s m x ^ J 0.00079 1 1 1 1 1 0 0.2 0.4 0.6 s 0.8 w Figure 3.11: Elastic P-wave velocity vs. water saturation in Spirit River sandstone (data given in Table 3.1) at 10 kHz, computed using the patch-scale relaxation model, illus-trating the dependence on the length scale of saturation heterogeneity. Each curve is labeled with its corresponding value of R0, in centimeters. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 96 of heterogeneity of 7.9 pm or less. Gist (1994) estimates for this rock sample a pore diameter of 0.3 pm, so that the concept of a "patch" is valid even at this small length scale. However, it is clear that if the wave frequency is increased to 1 MHz (i.e., to ultrasonic levels), complete fluid pressure relaxation cannot occur unless the "patches" have a radius of about 0 . 7 9 ( i . e . , of the same order as pore diameter, for which the very concept of a "patch" obviously fails). At this frequency he assumption used by le Ravalec et al. (1996) that patches are completely unrelaxed becomes valid for all patches with radius greater than 0.18 mm. Turning now to the particular application of modeling saturation hysteresis associated with patch-scale heterogeneity, we interpret the curves in Figures 3.10-3.11 corresponding to smaller values of Ro as being those that would result from processes that yield a more homogeneous distribution of fluids (e.g., imbibition). The curves corresponding to larger values of Ro are interpreted as those which would result from processes that yield a more heterogeneous distribution of fluids (e.g., drainage by evaporative drying). As an explicit example, we again consider an imbibition/drainage experiment with a water/air-saturated Spirit River sandstone. We assume that the imbibition proceeds leads to a distribution of pore fluids with a length scale of heterogeneity, Ro, of 0.1 mm. This is consistent with the scale of heterogeneity observed by Cadoret et al. (1995) for imbibition of various limestones by depressurization. We assume that the drainage process leads to a more heterogeneous distribution of fluids, with an associated length scale of heterogeneity of 1 cm. This is also consistent with the scale of heterogeneity observed by Cadoret et al. (1995) for drainage of various limestones by evaporative drying. In Figure 3.12 we present the resulting computed graph of P-wave velocity vs. Sw, for both imbibition and drainage at wave frequencies of 100 Hz and 10 kHz. The greater heterogeneity of the fluid distribution during drainage causes the model rock sample to exhibit hysteresis, with P-wave velocities being higher during drainage than during Chapter 3. Heterogeneity and Relaxation at the Patch Scale 97 4r 2 . 8 1 ' 1 1 1 ' 0 0.2 0.4 0.6 0.8 1 s w Figure 3.12: Computed saturation hysteresis of elastic P-wave velocity in Spirit River sandstone, at 100 Hz (—) and 10 kHz (—). Assumed length scales of saturation hetero-geneity are RQ = 0.1 mm for imbibition and R0 = 1 cm for drainage. 0 . 1 4 r 0 .12h ^ \ w Figure 3.13: Computed saturation hysteresis of P-wave attenuation Q~x for'Spirit River sandstone, at 100 Hz (—) and 10 kHz (—). Assumed length scales of saturation hetero-geneity are R0 = 0.1 mm for imbibition and R0 = 1 cm for drainage. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 98 imbibition. The computed hysteresis is most pronounced at 100 Hz, at which frequency the 0.1 mm scale liquid- and gas-saturated patches are completely relaxed, while the 1 cm scale patches are almost completely unrelaxed. At 10 kHz the 0.1 mm scale patches are less relaxed and the exhibited saturation hysteresis is less pronounced, though significant. The 1 cm scale patches are completely unrelaxed at this frequency, and the computed velocities during drainage are therefore slightly higher than at 100 Hz. The associated graphs of Q~1 vs. Sw, presented in Figure 3.13, also exhibit significant saturation hysteresis at both of the sampled wave frequencies. Naturally, the greatest attenuation is exhibited by the imbibition curve at 10 kHz, where the state of relaxation of the patches is intermediate between completely relaxed and completely unrelaxed, as can be seen in Figure 3.12. The curves with minimum attenuation are the one for imbibition at 100 Hz (where the patches are completely relaxed) and the one for drainage at 10 kHz (where the patches are completely unrelaxed). Interestingly, each of these curves exhibits a maximum attenuation at some intermediate saturation, corresponding to the optimal level of saturation such that the viscous losses due to flow between water-saturated and gas-saturated regions are maximized. The hysteresis curves presented in Figure 3.12 do not accurately model the particular form of the experimentally measured curves given for various limestones in Cadoret et al. (1995) and for Spirit River sandstone in Knight and Nolen-Hoeksema (1990). These data are more accurately reproduced when pore-scale effects of the imbibition/drainage process are taken into account, as in the model of le Ravalec et al. (1996). However, our computations do illustrate that saturation hysteresis due to patch-scale heterogeneity alone can be significant in realistic and geophysically relevant situations. Furthermore, the predicted behavior exhibits a frequency dependence than cannot be predicted by the high-frequency model of le Ravalec et al. The model presented here also provides esti-mates of the wave attenuation caused by viscous losses due to the patch-scale relaxation Chapter 3. Heterogeneity and Relaxation at the Patch Scale 99 mechanism, which have previously been provided only by the White (1975) model. 3.6 Summary In this chapter we have developed a frequency-dependent model for estimating elastic wave velocities in fluid-saturated porous media with "patch"-scale saturation hetero-geneity (i.e., media in which the length scale of saturation heterogeneity is on the order of many pore diameters). Our approach is similar to that taken in Chapter 2, in that we have sought to describe the frequency-dependent elastic response of a patch in terms of frequency-dependent, complex effective elastic moduli. With effective elastic moduli defined for each patch, an inclusion-based effective medium theory is used to estimate the effective moduli of the composite. The model is similar to that of le Ravalec et al. (1996), who also use an inclusion-based formulation to estimate elastic moduli of the composite from the effective elastic moduli of the individual patches. However, their model is lim-ited to the unrelaxed (i.e., high frequency or large patch size) limit,'whereas the model presented here considers the full frequency range within the limits of effective medium theory. Our approach has been to derive partial differential equations describing the diffusion of fluid pressure between patches. We construct a solution of these equations for the particular case of a spherical patch in an infinite porous medium subject to deformation that is sinusoidal in time (in this our model is similar to that of White (1975)). This solution is then used to relate the pressure p applied to a patch to its dilatation 6, so as to define an effective bulk modulus K* = —p/6 for the patch. Our results are consistent with the expected low-frequency limit as given by Gassmann's (1951) formula, and with the expected high-frequency, "unrelaxed" limit where there is no communication of fluid pressure between patches. Chapter 3. Heterogeneity and Relaxation at the Patch Scale 100 The chief limitation of the model presented here is that we consider saturation het-erogeneity at the patch scale only. In fact we require that the fluid distribution be homogeneous at the pore scale (or at least effectively homogeneous, in the sense that the time scale of relaxation of an initial pore-scale fluid pressure inhomogeneity is much shorter than the wave period). This is a serious limitation, since we expect that in many realistic porous materials heterogeneity will exist at many length scales ranging from the pore scale to the patch scale. However, under appropriate circumstances the model should provide a useful correction to the usual Gassmann estimate of effective elastic moduli (and hence wave velocities), for the case of saturation heterogeneity on a length scale intermediate between the pore-scale and the elastic wavelength. Because the effective patch moduli defined by our model are complex, they can be used to predict not only the velocities of elastic waves but also their attenuation due to viscous losses associated with flow of the pore fluids. This requires that the effective medium theory used to form the composite of all patches must allow for the constituents of the composite to have complex moduli. One such inclusion-based effective medium theory is that of Berryman (1980a; 1980b) used here. As a particular application of our model of patch-scale heterogeneity, we have consid-ered the effects of saturation hysteresis on elastic P-wave velocities. This application is similar to the saturation hysteresis model of le Ravalec et al. (1996) except that we have attempted to incorporate fluid pressure relaxation, whereas le Ravalec et al. consider only the unrelaxed limit. Our patch-scale relaxation model allows us to estimate the saturation hysteresis of both the velocity and attenuation of elastic waves for a broad range of wave frequencies. Chapter 4 Conclusion The objective of this thesis has been to develop predictive models of elastic wave veloci-ties in fluid-saturated porous media, with the particular aim of modeling the frequency-dependent effects of fluid distribution heterogeneity at various length scales. The models presented here complement previous models—which provide end-member velocity esti-mates in the limits of high and low frequency—by explicitly incorporating the dynamics of the fluid phase that determine the elastic behaviour of the porous composite at interme-diate frequencies. Our approach has been to incorporate a description of the mechanism of fluid pressure communication and relaxation into existing inclusion-based effective medium theories, thus extending their range of validity to lower wave frequencies. Two distinct length scales of saturation heterogeneity have been identified: the pore scale and the "patch" scale. We have considered each of these scales separately in Chap-ters 2 and 3, respectively. Although we choose different models for the dynamics of the pore fluid at these two scales, our approach to incorporating the effects of fluid pressure relaxation has been similar. In each case, we identify those constituents of the porous composite that exhibit frequency-dependent behavior due to deformation-induced flow of the pore fluids. In the particular cases of pore-scale and patch-scale heterogeneity, these constituents are the pores and the patches, respectively. We proceed by defining effec-tive complex elastic moduli for these constituents, and consider an equivalent medium in which the pores or patches are conceptually replaced by homogeneous materials having these elastic moduli. Berryman's (1980a; 1980b) inclusion-based effective medium theory 101 Chapter 4. Conclusion 102 is then used to estimate the effective elastic moduli (and hence the elastic wave velocities and attenuation) of the composite. The majority of our modeling efforts have been in deriving expressions for the ef-fective elastic moduli of the frequency-dependent constituents. In each case we seek a relationship between the overall stress and strain of a pore or patch. Noting that the overall behavior of a pore or patch under sinusoidal deformation is identical to that of some equivalent homogeneous solid, we use this relationship to define the elastic moduli of an equivalent medium with which a given pore or patch is then replaced. The overall stress-strain relationships used to define these effective moduli are derived by assuming a particular mechanism of fluid pressure relaxation (Poiseuille flow at the pore scale, and Darcy's law at the patch scale), and accounting for the contribution of this flow to the elastic response of a pore or patch. At both scales of heterogeneity, the frequency-dependent behavior of the composite arises because the degree to which fluid pressure communication between adjacent re-gions (i.e., between pores or between patches) occurs depends on the wave period. This frequency dependence can be described in terms of a characteristic relaxation time scale. At the pore scale, we have analyzed the communication of fluid pressure between two pores that are saturated with different fluids and are hydraulically connected by a narrow duct. For the particular case where the pores have identical volume V and the duct is cylindrical, we find that the time scale of pore-scale fluid pressure relaxation, r p o r e , is given by 8nLV ''"pore 7 r r 4 ( « i + « 2 ) ' where K\ and K2 are the bulk moduli of the two fluids, r and L are the radius and length of the duct, and n is the viscosity of the fluid in the duct. Previous models have for the most part corresponded to the "end members" in which Chapter 4. Conclusion 103 the wave period, r, satisfies either r >^ r p o r e (e.g., the Gassmann (1951) theory), or T ^ r pore (e.g., the inclusion-based formulations of Kuster and Toksoz (1974a) and Berryman (1980a; 1980b)). These end members correspond to the limits in which fluid pressures induced by a passing elastic wave are either completely relaxed or completely unrelaxed. By explicitly accounting for the mechanism of pore-scale fluid pressure relax-ation, the models developed in Chapter 2 are able to treat the full range of frequencies between these end members. In particular, one can use these models to estimate the elas-tic wave velocities and attenuation in the previously untreated regime where r ~ r p o r e . At the patch scale, we have analyzed the relaxation of fluid pressure between in a spherical "patch", e.g., a macroscopic region where the compressibility of the saturating fluid differs from that of the fluid in the surrounding medium. Our analysis yields the time scale of patch-scale fluid pressure relaxation, r p a t c h, given by <pr)R2 Tpatch ~ ^ K , where <p and K are the porosity and permeability of the porous matrix, n and Kf are the viscosity and bulk modulus of the fluid saturating the patch, and R is the radius of the patch. (This expression and its variants have been reported previously by other authors, e.g., Gist (1994) and Knight et al. (1998)). Again, previous models that allow for patch-scale heterogeneity provide only "end member" predictions in the cases where either r 3> r p a t c h , so that induced fluid pressures in patches are completely relaxed (e.g., as in the Gassmann (1951) theory), or r <?C r p a t c h, so that fluid pressures are completely unrelaxed (e.g., as in the IBEMT approach of le Ravalec et al (1996)). By explicitly accounting for the mechanism of patch-scale fluid pressure relaxation, the model developed in Chapter 3 treats the full frequency range between these end members. In particular, elastic wave velocities and attenuation can now be modeled in the potentially realistic regime where r ~ r p a t c h. Chapter 4. Conclusion 104 For both length scales of heterogeneity we have derived analytic expressions for the transition from low- to high-frequency effective elastic behavior that occurs when the wave period is decreased below the time scale of fluid pressure relaxation. In the low-frequency limit our models reduce to the expected results given by Gassmann's (1951) approach, and in the high-frequency limit reduce to the expected results given by conventional inclusion-based effective medium theory. Thus, our models can be seen as uniting the previously inconsistent poroelastic and inclusion-based approaches in a single unified model. Another significant advance of our models is that they provide estimates of the atten-uation of elastic waves due to the fluid pressure relaxation process. The effective elastic moduli that we define for pores and patches embody not only their elastic behavior, but also their attenuation characteristics via the imaginary parts of the moduli. While the physical mechanism of attenuation is viscous flow of the pore fluids, in our IBEMT ap-proach it is accounted for as intrinsic absorption by the equivalent elastic media used to conceptually replace the pores or patches. Our contributions to modeling the elastic behavior of heterogeneous, fluid-saturated porous media for a broad range of frequencies are summarized graphically in Figures 4.1-4.2. For a typical low-permeability porous rock in which a single length scale of hetero-geneity R (e.g., the typical pore size or the patch size) can be identified, Figure 4.1 illustrates the ranges of R and wave frequency, / , for which previous effective medium theories are applicable. The region of validity of the conventional pore-scale inclusion-based effective medium theory, where the fluid pressures in individual inclusions are assumed to be completely unrelaxed, is labeled "IBEMT". The other theories considered are the patch-scale IBEMT model of le Ravalec et al. (1996) and the Gassmann (1951) theory. The range of validity of each approach has been discussed in greater detail in previous chapters, and in Figure 4.1 we merely summarize our previous discussion. Chapter 4. Conclusion 105 A subtle point with respect to the patch-scale heterogeneity model of le Ravalec et al. is that the scale of heterogeneity ought strictly to be sufficiently large that the concept of a "patch" is meaningful. This requires that the scale of heterogeneity be at least many times greater than a typical pore diameter. However, such rigor is probably unnecessary in practice, and we anticipate no problems in applying their model to the case of intermediate scales of heterogeneity where the concept of a "patch" is ambiguous. Hence, in the interest of illustrative clarity, we neglect this possible limitation in the schematic pictorial representation of Figure 4.1. In contrast with Figure 4.1, Figure 4.2 illustrates the range of validity of the two new models developed in this thesis (i.e., the pore-scale and patch-scale relaxation models). By accounting for the mechanism of fluid pressure relaxation, the pore-scale relaxation model eliminates the restriction of conventional IBEMT to the high-frequency (i.e., unre-laxed) limit. Similarly, limitations on frequency—other than that of the effective medium theory assumption—have also been removed in the case of patch-scale heterogeneity. The key difference between Figures 4.1 and 4.2 is that, for a given scale of heterogeneity, the full frequency range can now be treated with a single unified model. As in Figure 4.1, there remains the subtle issue that in order to be strictly applicable, our patch-scale relaxation model requires that the scale of heterogeneity be sufficiently large for a "patch"—as well as the concept of permeability in treating the flow of the pore fluids—to be defined. However, in practice we suspect that our patch-scale relaxation model can be used to obtain meaningful results even as the "patch" size approaches the pore scale, and in Figure 4.2 we leave this possible limitation out of consideration. As noted in the introduction, the distribution of pore fluids in many porous media is inherently heterogeneous, on a length scale that can depend on the nature of the imbibi-tion/drainage process. Previous theoretical work on the effects of this heterogeneity on the elastic properties of porous media has focused on the end member corresponding to Chapter 4. Conclusion 106 wave frequency,/ [Hz] Figure 4.1: Schematic illustration of the ranges of wave frequency and length scale of heterogeneity for which various effective medium theories are valid, for a typical low-permeability rock. Here A is the elastic wavelength, and D is the fluid pressure diffusivity as defined in Chapter 3. Chapter 4. Conclusion 107 10° 102 104 106 10' wave frequency,/ [Hz] Figure 4.2: Schematic illustration of the ranges of wave frequency and length scale of heterogeneity for which the pore- and patch-scale relaxation models developed in this thesis are valid, for the same medium considered in Figure 4.1. Chapter 4. Conclusion 108 t h e h i g h f r e q u e n c y (i.e., u n r e l a x e d ) l i m i t , w h e r e t h e s e effects a r e m o s t p r o n o u n c e d . H o w -e v e r , i t i s c l e a r t h a t t h e effect o f h e t e r o g e n e i t y i s a f r e q u e n c y - d e p e n d e n t p h e n o m e n o n — i n d e e d , t h e ef fec ts o f h e t e r o g e n e i t y d i s a p p e a r a t s u f f i c i e n t l y l o w f r e q u e n c i e s , w h e r e t h e h e t e r o g e n e i t y i s " s m o o t h e d o u t " b y f l o w o f t h e p o r e f l u i d s , a n d G a s s m a n n ' s ( 1 9 5 1 ) t h e -o r y p r o v i d e s a n a c c u r a t e m o d e l o f t h e e f f e c t i v e e l a s t i c b e h a v i o r . T h e p r e v i o u s h i g h - a n d l o w - f r e q u e n c y m o d e l s a r e i n s u f f i c i e n t t o p r o v i d e a d e s c r i p t i o n o f t h e f u l l r a n g e o f e l a s t i c b e h a v i o r . O u r a n a l y s e s i n d i c a t e t h a t t h i s i s e s p e c i a l l y t r u e f o r p a t c h - s c a l e h e t e r o g e n e i t y , w h e r e t h e i n t e r m e d i a t e r e g i m e b e t w e e n t h e r e l a x e d a n d u n r e l a x e d e n d m e m b e r s o c c u p i e s m a n y d e c a d e s o f w a v e f r e q u e n c y . T h e m o d e l s d e v e l o p e d i n t h i s t h e s i s c o n t r i b u t e t o r ec -o n c i l i n g t h e s e e n d m e m b e r s , p r o v i d i n g e s t i m a t e s o f n o t o n l y t h e e l a s t i c w a v e v e l o c i t i e s , b u t a l s o a t t e n u a t i o n d u e t o t h e r e l a x a t i o n p r o c e s s , i n t h e i n t e r m e d i a t e r e g i m e b e t w e e n p r e v i o u s p o r o e l a s t i c a n d i n c l u s i o n - b a s e d m o d e l s . References Berryman, J. G. (1980a). Long-wavelength propagation in composite elastic media. I. Spherical inclusions. J. Acoust. Soc. Am., 68(6):1809-1819. Berryman, J. G. (1980b). Long-wavelength propagation in composite elastic media. II. Ellipsoidal inclusions. J. Acoust. Soc. Am., 68(6):1820-1831. Biot, M . A. (1956a). Theory of propagation of elastic waves in a fluid-saturated porous solid. I. Low-frequency range. J. Acoust. Soc. Am., 28(2):168-178. Biot, M . A. (1956b). Theory of propagation of elastic waves in a fluid-saturated porous solid. II. Higher frequency range. J. Acoust. Soc. Am., 28(2):179-191. Bourbie, T., Coussy, 0., and Zinszner, B. (1987). Acoustics of Porous Media. Gulf Publishing Co., Houston. Budiansky, B. and O'Connell, R. J. (1976). Elastic moduli of a cracked solid. Int. J. Solids & Structures, 12:81-97. Cadoret, T. (1993). Effet de la saturation eau-gaz sur les proprietes acoustiques des roches. Etude aux frequences sonores et ultrasonores. PhD thesis, Universite de Paris VII. Cadoret, T., Marion, D., and Zinszner, B. (1995). Influence of frequency and fluid distribution on elastic wave velocities in partially saturated limestones. J. Geophys. Res., 100(6) :9789-9803. Coyner, K. B. (1984). Effects of stress, pore pressure, and pore fluids on bulk strain, velocity, and permeability in rocks. PhD thesis, Massachusetts Institute of Technology. Domenico, S. N. (1977). Elastic properties of unconsolidated sand reservoirs. Geophysics, 42(7):1339-1368. Dvorkin, J . , Nolen-Hoeksema, R., and Nur, A. (1994). The squirt flow mechanism: Macroscopic'description. Geophysics, 59(3):428-438. Dvorkin, J. and Nur, A. (1993). Dynamic poroelasticity: A unified model with the squirt and the Biot mechanisms. Geophysics, 58(4):524-533. Endres, A. L. and Knight, R. J. (1989). The effect of microscopic fluid distribution on elastic wave velocites. The Log Analyst, 30(6):437-445. 109 References 110 Endres, A. L. and Knight, R. J. (1997). Incorporating pore geometry and fluid pres-sure communication into modeling the elastic behavior of porous rocks. Geophysics, 62(1):106-117. Eshelby, J. D. (1957). The determination of the elastic field of an ellipsoidal inclusion, and related problems. Proc. R. Soc. London, Ser. A, 241(l):376-396. Gassmann, F. (1951). Uber die elastizitat poroser medien. Vierteljahrsschrift der Natur-forschenden Gesellschaft in Zurich, 96:1-23. Gist, G. A. (1994). Interpreting laboratory velocity measurements in partially gas-saturated rocks. Geophysics, 59(7):1100-1109. Gueguen, Y. and Palciauskas, V. (1994). Introduction to the physics of rocks. Princeton University Press. Hashin, Z. and Shtrikman, S. (1963). Self-consistent elastic moduli of a cracked solid. J. Mech. Phys. Solids, 11:127-140. Henyey, F. S. and Pomphrey, N. (1982). Self-consistent elastic moduli of a cracked solid. Geophys. Res. Lett, 9(8):903-906. Knight, R., Dvorkin, J. , and Nur, A. (1998). Acoustic signatures of partial saturation. Geophysics, 63(1): 132-138. Knight, R. and Nolen-Hoeksema, R. (1990). A laboratory study of the dependence of elas-tic wave velocities on pore scale fluid distribution. Geophys. Res. Lett., 17(10):1529-1532. Kuster, G. T. and Toksoz, M. N. (1974a). Velocity and attenuation of seismic waves in two-phase media: Part I. Theoretical formulations. Geophysics, 39(5):587-606. Kuster, G. T. and Toksoz, M . N. (1974b). Velocity and attenuation of seismic waves in two-phase media: Part II. Experimental results. Geophysics, 39(5):607-618. le Ravalec, M . and Gueguen, Y. (1996). High- and low-frequency elastic moduli for a saturated porous/cracked rock—differential self-consistent and poroelastic theories. Geophysics, 61(4):1080-1094. le Ravalec, M . , Gueguen, Y. , and Chelidze, T. (1996). Elastic wave velocities in partially saturated rocks: saturation hysteresis. J. Geophys. Res., 101(Bl):837-844. Lucet, N. (1989). Vitesse et attenuation des ondes elastiques soniques et ultrasoniques dans les roches sous pression de confinement. PhD thesis, Universite de Paris VI. References 111 Mavko, G. and Jizba, D. (1991). Estimating grain-scale fluid effects on velocity dispersion in rocks. Geophysics, 56(12):1940-1949. Mavko, G., Mukerji, T. , and Dvorkin, J. (1996). Rock Physics Handbook. Rock Physics Laboratory, Stanford University. Murphy, W. F. (1985). Sonic and ultrasonic velocities: Theory versus experiment. Geo-phys. Res. Lett, 12(2):85-88. Murphy, W. F., Winkler, K. W., and Kleinberg, R. L. (1986). Acoustic relaxation in sedimentary rocks: Dependence on grain contacts and fluid saturation. Geophysics, 51(3):757-766. Nur, A. and Byerlee, J. D. (1971). An exact effective stress law for elastic deformation of rocks with fluids. J. Geophys. Res., 76(26):6414-6419. O'Connell, R. J. and Budiansky, B. (1974). Seismic velocities in dry and saturated cracked solids. J. Geophys. Res., 79(35):5412-5426. Press, W. FL, Teukolsky, S. A., Vetterling, W. T., and Flannery, B. P. (1992). Numerical Recipes in C, pages 362-367. Cambridge University Press, New York, second edition. Walpole, L. J. (1971). The elastic behaviour of a suspension of spherical particles. Quart. J. Mech. & App. Math., 25:153-160. Wang, Z., Hirsche, W. K., and Sedgwick, G. (1991). Seismic velocities in carbonate rock. J. Can. Petr. Tech., 30:112-122. White, F. M . (1986). Fluid Mechanics. McGraw-Hill, New York, second edition. White, J. E . (1975). Computed seismic speeds and attenuation in rocks with partial gas saturation. Geophysics, 40(2):224-232. Wu, T. T. (1966). The effect of inclusion shape on the elastic moduli of a two-phase material. Int. J. Solids & Structures, 2:1-8. Xu, S. (1998). Modelling the effect of fluid communication on velocities in anisotropic rocks. Int. J. Solids & Structures, 35(34-35):4685-4707. Appendix A Probability Distribution of Distances Between Neighboring Pores In this appendix we derive the probability density function for the distance between nearest neighbors in a spatially random distribution of inclusions. More specifically, we consider a spatially random distribution of points in a porous medium which represent the centers of inclusions, and seek the probability distribution for the distance from a given point to its nearest neighboring point, irrespective of the actual shape or size of the inclusion a given point represents. Consider a spherical porous medium sample, having volume V, in which is scattered a random distribution of N points labeled 7 T 0 . . . 7 T / v - i , and consider a particular point 7 r 0 . We wish to determine the probability density function P(r), defined so that the probability of the distance from ir0 to its nearest neighboring point being between r and r+dr is given by P(r)dr, where dr is taken to be an infinitesimal quantity. This situation is represented schematically in Figure A . l . In this figure the spherical sample is divided into three labeled sub-regions: QA, in which the distance to the point 7 r 0 is less than r; Q£, where this distance is between r and r + dr; and 0 B , where this distance is greater than r. For convenience in performing the necessary logical manipulations, we define the symbols R(Q), S(Cl), and T(Q) to represent the following statements, = "None of the points {TTI ... 7rjv_i} is in the region 0" S(fl) = "At least one of the points { n i . . . 7T/v-i} is in the region f2" T(Cl) = "All of the points {ni.. .7r^_i} are in the region 0", 112 Appendix A. Probability Distribution of Distances Between Neighboring Pores 113 Figure A . l : Geometrical construction used in deriving the probability distribution for the nearest-neighbor distance in a spatially random distribution of inclusions. and use the notation to represent the probability that statement X is true. The problem described above is now, as expressed in this abbreviated notation, to find p{R(ttA) and 5(f2e)}. In order to facilitate finding an analytic expression for this quantity, we perform the following series of logical manipulations: p{R(nA) and 5(a)} = P{R(^A)} X p{5(0£), given R(QA)} = p{T(nE+B)} x p{S(Q£), given T(Sle+B)} (A.l) = p{T(a+fl)} x (1 - p{T(ttB), given T(0 ( £ + B ) } ) . The shorthand notation 0 £ + B has been used to represent the union of the regions Q£ and tin-For a single point chosen randomly within some volume V, the probability of it being chosen in some sub-volume V0 is given by V0/V. For a collection of points labeled 7 T i . . . 7T/v, the probability p{{7Ti . . . 7TJV} in Vo} that all points are in the sub-volume Vo is Appendix A. Probability Distribution of Distances Between Neighboring Pores 114 given by p{{7Ti . . . T T i v } in Vo} = p{7Ti in Vo} x p{^2 in Vb} x • • • x p{nN in V0} (A.2) Vb V Using this result, the final expression in equation (A.l) can then be written p{R(QA) and S(fl6)} = Ve + VB V V-VA V - 2 N-l X 1 - B N-l N-l N-l ve + vB v-vA N-l (A.3) 1 - 1 V, N-l V-VA. where VA, V£ and VB are the volumes of the regions QA, a and r2s, respectively. To avoid "boundary effects" due to finite sample size, it is convenient to take the limiting form of the final expression in equation (A.3) as V and N go to infinity, in such a way that the number density of inclusions, n = N/V, is held fixed. This approximation simplifies the results, and is valid in the case where the sample dimensions are much greater than the average pore-to-pore distance. To this end, this expression for n is used to substitute V = N/n for all occurrences of V in equation (A.3) to yield p{R(ClA) and S(tte)} = 1 nVA N N-l 1 - 1 -nV N — nVA N—l (A.4) We then take the limit of expression (A.4) as N (and simultaneously V) approaches infinity to obtain lim p{R(QA) and 5 ( a ) } = e~nVA (l - e~nVs) N-*oo (A.5) Appendix A. Probability Distribution of Distances Between Neighboring Pores 115 Making use of the fact that V£ is an infinitesimal quantity, we compute a first order expansion of expression (A.5) for small V£ to obtain the following result for an infinite sample: p{R(nA) and S(fi e)} = nV£e-nVA (A.6) = n^irr2dre n 3 7 r r ' 3 . From equation (A.6) it is clear that the expression for P(r) that we are seeking is given by P(r) = 47mr 2e-i™ r 3 (A.7) It is easily verified that P(r) as defined in equation (A.7) satisfies the normalization criterion roo / P{r) dr = 1. (A.8) Jo That is, for finite inclusion density n, the probability that the nearest neighbor is some-where between 0 and oo is equal to 1 (i.e., the nearest neighbor must be somewhere). It is convenient to non-dimensionalize equation (A.7) by defining the length scale r 0 = n - 1 / 3 and the dimensionless variable r* = r/r0. Then we can write P(r)dr = P*{r*)dr* (A.9) where the probability density function P*(r*) is given by P*(r*) = 4ivr2e--3*r*. (A. 10) The fact that this expression is independent of n demonstrates that distances in this problem scale simply as rC1^. A graph of P*(r*) is presented in Figure A . l . We can see from this graph that most Appendix A. Probability Distribution of Distances Between Neighboring Pores 116 2 CL 0 0 0.5 1 r. 1.5 Figure A.2: Graph of the probability distribution for the nearest-neighbor distance in a spatially random distribution of inclusions. inclusions have their nearest neighbor at a distance less than n - 1 / 3 , but within this range the distribution is quite broad, so that the effect of this distribution on pore-scale relaxation phenomena is expected to be significant. The probability distribution tails off rapidly for r > n - 1 / 3 , and only a small fraction, equal to of inclusions have their nearest neighbor at a distance greater than n 1^3. The mean distance (r) from a given inclusion to its nearest neighbor, given by (A.ll) ( 3 ) 1/3 (r) = / rP(r) dr Ann (A.12) Jo significantly less than the naive estimate of n 1 / / 3 proposed in section 2.4.1. List of Notation Ci Volume fraction of the z t h constituent of a porous composite 15 D Diffusivity of fluid pressure 59 Dp Diffusivity of fluid pressure under constant sample pressure 66 De Diffusivity of fluid pressure under constant sample dilatation 65 / Wave frequency 104 F Parameter used in Berryman's (1980b) IBEMT 15 K Intrinsic permeability of a porous medium 59 L Length of the duct connecting two pores 21 n Number density of inclusions 25 V Pressure 18 Pf Pressure in the pore fluid 39 ps Pressure in the mineral solid component 63 P*1 A geometrical term used in Berryman's (1980b) IBEMT 15 q Volumetric flow rate between inclusions 21 QH A geometrical term used in Berryman's (1980b) IBEMT 15 r Radius of a cylindrical duct connecting two pores 21 R Radius of a spherical patch 59 Si Volume fraction of the pore space occupied by the ith pore fluid 19 Sw Volume fraction of the pore space occupied by liquid (usually water) 4 t Time 18 V Volume of a pore 21 vp Elastic P-wave velocity 14 117 List of Notation 118 vs Elastic S-wave velocity 14 a Poroelastic parameter equal to 1 — K d / K s 3 7 Poiseuille flow coefficient, equal to 7rr 4/(8nL) for a spherical duct 21 TJ Fluid viscosity 22 0 Volumetric dilatation 18 9f Volumetric dilatation of the pore fluid 63 6S Volumetric dilatation of the mineral solid component 63 K Bulk modulus of an isotropic medium 14 K* Effective complex bulk modulus of a pore or patch 19 k Parameter used in the pore-to-pore relaxation model 23 Kd Dry frame bulk modulus 3 fteff Effective bulk modulus of a porous composite 3 Kf Bulk modulus of the pore fluid 3 ttyff Bulk modulus of the "effective pore fluid" 4 Ki Bulk modulus of the z t h constituent of a porous composite 14 KS Bulk modulus of the mineral solid component of a porous matrix 3 A Elastic wavelength 106 p Shear modulus of an isotropic medium 14 p* Effective complex shear modulus of a pore 29 Pd Dry frame shear modulus 3 ^eff Effective shear modulus of a porous composite 3 Pi Shear modulus of the ith constituent of a porous medium 15 T Time scale of fluid pressure relaxation 59 v Kinematic fluid viscosity equal to n/p 28 p Density 14 List of Notation 119 0 Porosity (i.e., the volume fraction occupied by the pore space) 3 to Angular wave frequency equal to 2nf 18 UJQ Characteristic frequency of a relaxation mechanism 23
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Modeling elastic wave velocities in porous media :...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Modeling elastic wave velocities in porous media : frequency-dependent effects of heterogeneity at the… Taylor, Richard 1999
pdf
Page Metadata
Item Metadata
Title | Modeling elastic wave velocities in porous media : frequency-dependent effects of heterogeneity at the pore- and patch-scale |
Creator |
Taylor, Richard |
Date Issued | 1999 |
Description | The prediction of relationships between elastic wave velocities in a porous medium and the properties of the fluid and solid constituents therein is a longstanding problem in geophysical exploration. Previous authors have shown that such relationships depend on both the wave frequency and degree of heterogeneity present. This frequency dependence arises via the state of relaxation of the pore fluids. At sufficiently low frequencies, deformation-induced flow leads to equilibration of the fluid pressure, and the pore fluids are said to be in a "relaxed" state. At sufficiently high frequencies, there is insufficient time for equilibration to occur, and the fluids are said to be "unrelaxed". Current models of elastic wave velocities in porous media are, for the most part, confined to either the relaxed limit (e.g., the poroelastic Biot-Gassmann theory) or the unrelaxed limit (e.g., inclusion-based effective medium theory). In this thesis we incorporate an explicit description of the relaxation mechanism into inclusion-based effective medium theory, so as to extend the theory toward the relaxed limit. Analysis of the mechanisms of relaxation leads to a description of the effective elastic behavior of the porous medium in terms of effective complex elastic moduli of the medium's constituents. Previous authors have identified two distinct scales of fluid distribution heterogeneity: the pore scale and the patch scale. Accordingly, we treat these scales separately, describing relaxation in terms of Poiseuille flow at the pore scale, and Darcy's law at the patch scale. The results of our analyses are effective medium theories that provide a consistent approach to the prediction of elastic wave velocities, as well as attenuation due to the relaxation mechanism, over a broad range of frequencies and length scales of heterogeneity. In particular, our model is applicable to the regime where the pore fluids are in a state of relaxation intermediate between the completely relaxed and completely unrelaxed end members. |
Extent | 7760978 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Date Available | 2009-06-16 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0089086 |
URI | http://hdl.handle.net/2429/9356 |
Degree |
Master of Science - MSc |
Program |
Geophysics |
Affiliation |
Science, Faculty of Earth, Ocean and Atmospheric Sciences, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 1999-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_1999-0445.pdf [ 7.4MB ]
- Metadata
- JSON: 831-1.0089086.json
- JSON-LD: 831-1.0089086-ld.json
- RDF/XML (Pretty): 831-1.0089086-rdf.xml
- RDF/JSON: 831-1.0089086-rdf.json
- Turtle: 831-1.0089086-turtle.txt
- N-Triples: 831-1.0089086-rdf-ntriples.txt
- Original Record: 831-1.0089086-source.json
- Full Text
- 831-1.0089086-fulltext.txt
- Citation
- 831-1.0089086.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0089086/manifest