Oliver S. Krüger

Finite-difference modeling of wave propagation on microscale: A snapshot of the work in progress

Digital rock methodology combines modern microscopic imagingwithadvancednumericalsimulationsofthephysicalproperties of rocks. Modeling of elastic-wave propagation directly from rock microstructure is integral to this technology. We survey recent development of the rotated staggered grid RSG finite-difference FD method for pore-scale simulation of elastic wavepropagationindigitalrocksamples,includingthedynamic elastic properties of rocks saturated with a viscous fluid. Examination of the accuracy of this algorithm on models with known analytical solutions provide an additional accuracy condition for numerical modeling on the microscale. We use both the elastic and viscoelastic versions of the RSG algorithm to study gas hydratesandtosimulatepropagationofBiot’sslowwave.Weapply RSG method ology to examine the effect of gas hydrate distributions in the pore space of a rock. We compare resulting P-wave velocities with experimentally measured data, as a basis for buildinganeffective-mediummodelforrockscontaininggashydrates. We then perform numerical simulations of Biot’s slow wave in a realistic 3D digital rock model, fully saturated with a nonviscous fluid corresponding to the high-frequency limit of poroelasticity, and placed inside a bulk fluid.The model clearly demonstrates Biot’s slow curve when the interface is open between the slab and bulk fluid.We demonstrate slow wave propagation in a porous medium saturated with a viscous fluid by analyzing an idealized 2D porous medium represented alternating solid and viscous fluid layers. Comparison of simulation results withtheexactsolutionforthislayeredsystemshowsgoodagreementoverabroadfrequencyrange.

A numerical study on reflection coefficients of fractured media

Effective-medium theories can be used to predict reflection coefficients of an interface between an unfractured layer overlying a fractured half-space. In 2D and 3D computer simulations, we analyze wavefields that are emitted by an explosion line or point source and reflected from a fractured area in a digital rock model. The reflection coefficients from the simulations are compared to several predictions given by static effective-medium theories. The agreement between our numerical results and the theoretical reflection coefficients is best for differential effective-medium schemes in 2D as well as in 3D. This result agrees with previously published numerical static and numerical transmission-time experiments. By varying the wavelength to crack-length ratio, we consider the application range of different effective-medium theories. We observe good agreement with theoretical predictions even for a ratio equal to 9. An angle-dependent analysis of the reflected amplitude of our numerical results is also compared to results given by effective-medium theory in combination with exact reflection-coefficient formulas. As expected, the fluctuations of the reflection coefficients decrease for wider angles. In our 3D digital rock models, we vary the crack density and the infill of the crack (i.e., water [cold and hot] and oil [cold and hot]).

Leaky mode: A mechanism of horizontal seismic attenuation in a gas-hydrate-bearing sediment

The leaky mode is a possible attenuation mechanism of seismic waves propagating along lamination in gas-hydrate-bearing sediment layers. This horizontal propagation attenuation mechanism occurs when a high-velocity layer is embedded in a low-velocity zone. This is a typical situation for gas hydrate occurrences. To quantify this attenuation mechanism, a 2D digital rock model based on the crosswell data of the Mallik 2002 Gas Hydrate Production Research Well Program is used. For simplicity, our elastic simulations exclude attenuation mechanisms like scattering loss or intrinsic absorption. We demonstrate that the leaky mode is a significant horizontal attenuation mechanism that cannot be neglected. The effective attenuation of gas-hydrate-bearing sediments is a combination of intrinsic and scattering attenuation by small-scale heterogeneties and the leaky mode.

Reflection coefficients of fractured rocks: A numerical study

In this work we estimate the effective reflection coefficients of an interface between a cracked and an uncracked material. The study is based on computer simulations using the rotated staggered grid finite difference method. The numerically obtained reflection coefficients are compared to several theoretical predictions from static effective medium formulations. The agreement between our numerical data and the theoretical predictions is best for the differential schemes. This result is supported by previous studies dealing with transmission experiments. INTRODUCTION Theoretical effective medium descriptions of fractured media are commonly derived from static considerations (see e.g. Mavko et al., 1998, and references therein). Alternatively, numerical estimations of the effective elastic properties of such media have been successful Arns et al. (2002); Saenger and Shapiro (2002); Orlowsky et al. (2003). As the acoustic impedance is altered by the presence of cracks, interfaces between fractured and unfractured materials act as a reflecting boundary for elastic waves. For example, fault zones usually are associated with fractured zones. Recognising their characteristics in seismograms and the knowledge of their reflection coefficients are therefore of great interest in seismic exploration. These reflection coefficients, giving the ratio of the amplitudes of reflected to incident waves, describe a dynamic process. The question we pose in this paper is whether reflection coefficients, which are derived from the static elastic moduli gained from effective media theories, can be used to understand the dynamic process of wave reflection, including effects like the dependence of the reflection coefficients on the angle of incidence (i.e. AVO or AVOZ, see e.g. Castagna and Backus, 1993). As no analytical solutions to the wave equation exist for complex situations like strongly cracked solids, we rely on a numerical finite difference (FD) technique, the rotated staggered grid, to numerically simulate wave propagation and to determine reflection coefficients of the interface between cracked and uncracked areas. The rotated staggered grid as described in Saenger et al. (2000) allows one to simulate wave propagation in highly heterogeneous media without implementing explicit boundary conditions and without averaging elastic moduli. It has been proven to yield stable and realistic results for cracked media Krüger et al. (2002). Our experiments consist of simulations of two dimensional models with a plane wave source at the top of the model illuminating a cracked region at some distance from the source. We present numerical results for the reflection coefficient for normal incidence on cracked areas with different • crack densities, • inclination of the cracks, • dominant frequencies of the source wavelet and compare them with theoretical predictions. Further parameter variations, which remain to be done, may include the extensions of the cracked region, an inclination of the cracked region and the position of the source. 214 Annual WIT report 2004 THEORETICAL PREDICTION OF REFLECTION COEFFICIENTS IN CRACKED MEDIA Effective medium theories usually predict static elastic moduli. From these moduli, wave propagation velocities, acoustic impedances and hence reflection coefficients can be calculated. Here we focus on three theoretical predictions for parallel cracks, aligned along the x-axis. These are the non-interacting approximation (NIA), the differential scheme (DS) and an extension of the differential scheme (EDS). All theories are discussed in detail in Orlowsky et al. Orlowsky et al. (2003) and references therein. For the NIA, the energy that is needed to insert a single crack into an uncracked medium is simply added to the elastic potential of the medium. The DS recalculates the effective elastic moduli after the insertion of each crack, thus taking the effects of the so far inserted cracks into account. The EDS works like the NIA but multiplies the energy increment for each crack with a function which is derived from expressions for effective moduli corresponding to the DS results for a isotropic crack distribution. Predictions from these theories are given as a function of the crack density ρcd, which is defined forN cracks of radii ai distributed over an area A Bristow (1960): ρcd = 1/A N ∑ i=1 ai . (1) Young’s modulus E2 along the z-axis and shear modulus G for a crack density ρcd parallel to the x-axis are then given by: E2(NIA) = E · [1 + 2πρcd]−1, E2(DS) = E · e−2πρcd , E2(EDS) = E · [1 + 2πρcdecd ]−1, G(NIA) = G · [1 + πρcd(1− ν)]−1, G(DS) = G · e−π(1−ν)ρcd , G(EDS) = G · [1 + πρcd(1− ν)ecd ]−1, (2) where E denotes the Young’s modulus and ν the Poisson’s ratio of the background matrix. Note that E1 (along the x-axis) equalsE, which means that it is not affected by the cracks. Using these effective moduli, the stiffness matrix for cracks parallel to the x-axis can be calculated according to   < c11 > < c12 > 0 < c12 > < c22 > 0 0 0 < c44 >   =   1−ν E1 − ν(1+ν) E 0 − ν(1+ν) E 1−ν 2 0 0 0 1   −1 . (3) The corresponding effective velocities follow from: vP,eff = √ c22/ρg,eff , vS,eff = √ c44/ρg,eff . (4) The normal-incidence reflection coefficientRPP of an interface between two materials with velocities vP,1 and vP,eff , resp., and densities ρ1 and ρeff , resp., is given by: RPP = vP,effρeff − vP,1ρ1 vP,1ρg,1 + vP,effρeff . (5) MODELING PROCEDURE The two dimensional models consist of a homogeneous medium. Into the lower part of that medium dry cracks are placed, varying in number and inclination for different simulations. Different types of wavelets (Gaussian type and Ricker I) and dominant frequencies (fdom) are used for the plane wave, which propagates from the top of the model downwards (see Fig. 1Tab. 2). The dimension of the models is 1501 Annual WIT report 2004 215

Numerical Considerations of Fluid Effects on Wave Propagation

This paper is concerned with numerical considerations of fluid effects on wave propagation. The focus is on effective elastic properties (i.e. velocities) in different kinds of dry and fluid-saturated fractured media. We apply the so-called rotated staggered finite-difference grid (RSG) technique. Using this modified grid it is possible to simulate the propagation of elastic waves in a 2D or 3D medium containing cracks, pores or free surfaces without explicit boundary conditions and without averaging elastic moduli. Therefore the RSG allows an efficient and precise numerical study of effective velocities in fractured structures. This is also true for structures where theoretically it is only possible to predict upper and lower bounds. We simulate the propagation of plane P- and S-waves through three kinds of randomly cracked 3D media. Each model realization differs in the porosity of the medium and is performed for dry and fluid-saturated pores. The synthetic results are compared with the predictions of the well known Gassmann equation and the Biot velocity relations. Although we have a very low porosity in our models, the numerical calculations showed that the Gassmann equation cannot be applied for isolated pores (thin penny-shaped cracks). For Fontainebleau sandstone we observe with our dynamic finite-difference approach the exact same elastic properties as with a static finite-element approach. For this case the Gassmann equation can be checked successfully. Additionally, we show that so-called open-cell Gaussian random field models are an useful tool to study wave propagation in fluid-saturated fractured media. For all synthetic models considered in this study the high-frequency limit of the Biot velocity relations is very close to the predictions of the Gassmann equation. However, using synthetic rock models saturated with artificial “heavy” water we can roughly estimate the corresponding tortuosity parameter.

Geometric control of earthquake magnitudes by fluid injections in rocks

Pore fluids in rocks and pore-pressure relaxation can induce aftershocks of earthquakes. Sometimes rock stimulations by fluid injections into geothermal boreholes are able to trigger perceptible or even potentially damaging earthquakes. This seems to be not the case by hydraulic fracturing of hydrocarbon reservoirs. Reasons of such a difference and factors defining magnitudes of induced earthquakes (triggered tectonically as well as induced artificially) remain unclear. Here we show that one of the main factors limiting the probability to induce a large-magnitude event is the minimum principal axis of a fluid-stimulated rock volume. This geometrical scale controls the order of a largest possible magnitude. We analyze an impact of the geometry of a stimulated volume on the Gutenberg–Richter-type frequency–magnitude distribution of induced earthquakes. It seems that a rupture is only probable along a surface located mainly inside a stimulated rock volume. Therefore, monitoring of a spatial growth of seismicity in real time can help to constrain a risk of damaging induced earthquakes.

Numerical Study of Transmission Signatures of Gas Hydrate-Bearing Microstructures

This work is inspired by the observation, that gas hydrate bearing sediments have a high velocity in combination with high attenuation. We study numerically the influence of different gas hydrate locations within the pore space on transmitted p-waves. From the wave propagation simulations on the micro--scale it can be seen, that different positions of the gas hydrate in the pore space results in almost the same effective velocities and attenuation, as long as the gas hydrate had contact to the sediment grains. This changes in the case of a suspension, here the attenuation increases and the effective velocity decreases. The resulting p-wave versus gas hydrate saturation plot is in a qualitatively good agreement with experimental results obtained for the Mallik 2L-38 well.