Computation of wave attenuation and dispersion, by using quasi-static finite difference modeling method in frequency domain

Computation of wave attenuation and dispersion, by using quasi-static finite difference modeling method in frequency domain Qazi Adnan Ahmad1,*, Guochen Wu2, Wu Jianlu3 1 China University of Petroleum (East China), School of Geoscienc, Qingdao, Shandong, China 2 China University of Petroleum (East China), College of Geo Resources and Information, Department of Geophysics Dong Ying, Shan Dong, China 3 China University of Petroleum Huadong, Dongying, Shandong, China


Introduction
The main challenge we are facing in oil and gas exploration is the estimation of effects on wave characteristics due to the presence of multiscale subsurface heterogeneities.During wave propagation, the heterogeneous nature of earth subsurface results into pressure gradient at different spatial scale.Consequently, the resulted pressured gradient will accelerate the fluid to flow (at different scales) from high pressured zone to comparatively low-pressure zone.The wave induced fluid flow (WIFF) at wavelength scale is named as macroscopic, whereas, at pore scale, it is termed as microscopic, on the other hand, at a scale much smaller than wavelength but larger than the pores scale it is named as mesoscopic fluid flow [Muller et al. 2010].
The influences of fluids on seismic responses remain a key topic among researchers.Biot's [M. A. Biot 1956, M. A. Biot 1956] classical theory of poroelasticity led the investigation regarding wave attenuation and velocity dispersion in a saturated porous medium.In his theory, macroscopic scale heterogeneities due to presence of single phase fluid were outlined as predicted cause of wave attenuation and velocity dispersion.Wave propagation through such macroscopic heterogeneous media creates pressure gradient at wavelength scale, equilibration of this pressure gradient results into a loss in wave energy.Two compressional and one shear wave was predicted in Biot's theory, experimental results also validate the presence of slow P-wave [Plona 1980, Kelder et al. 1997].Meanwhile, another promising feature in Biot's theory is the low-frequency limiting velocities are identical with Gassmann predicted velocities.Although Biot's research opens new horizons in the field of geophysical prospecting, however, the outcomes of Biot's theory deals significantly at higher frequencies subsequently, it exceedingly under estimate wave attenuation and velocity dispersion within low-frequency limits [Vogelaar et al. 2007].Accordingly, extrinsic attenuation (scattering) dominates in Biot's theory, when the wavelengths are of the order of grain scale (microscopic scale).Therefore, some new techniques are required to overcome these existing challenges.
In order to compute wave attenuation and velocity dispersion at grain (microscopic) scale Mavko and Nur [Mavko et al. 1975, Mavko et al. 1979, Dvorkin et al. 1995] presented the concept of fluid flow from soft crack to stiff pores, that eventually causes wave attenuation and velocity dispersion.Several analytical solutions are also put forward to estimate wave attenuation and velocity dispersion due to WIFF at microscopic scale [Mavko et al. 1975, Mavko et al. 1979, Budiansky et al. 1976, Palmer et al. 1980, Dvorkin et al. 1995, Chapman et al. 2002, Pride et al. 2003, Gurevich et al. 2009].Apart from the ability to precisely predict wave attenuation and velocity dispersion at sonic frequency frequencies, investigations reveal that the microscopic scale theories lack to predict wave attenuation and velocity dispersion within seismic frequency band.This inadequacy motivates the researchers in opening new avenues to estimate the influences of WIFF at mesoscopic scale.
The study of wave propagation through saturated porous media demands to estimate inevitable wave attenuation and velocity dispersion within seismic frequency band.Wave propagation through saturated porous rock exhibiting heterogeneity at mesoscopic scale (which is larger than microscale but much smaller than the macro scale) causes significant attenuation and velocity dispersion within a seismic frequency band [Muller et al. 2010].Pressure gradient arises at mesoscopic scale, due to presence of heterogeneity within solid (i.e. between complaint pores/cracks and stiff pores) and fluids (i.e.due to compressibility difference between saturating fluids) at mesoscopic scale.This heterogeneity give rise to fluid flow at mesoscopic scale that eventually leads to wave attenuation and velocity dispersion [Pride 2004].
Within a porous rock presence of two or more immiscible fluids at mesoscopic scale arises patchy satura-tion.Exploring patchy saturation in a spherical manner, White [White 1975] conducted a significant work and declared that mesoscopic fluid flow plays an influential role in wave attenuation and velocity dispersion.White's concept was further extended to patchy saturation in layered form [White et al. 1975].In subsequent studies, Dutta and Ode [Dutta et al. 1979] further investigated gas-water patchy saturation in terms of Biot's poroelastic theory [M. A. Biot 1956].Johnson [Johnson 2001] analyzed consequences of arbitrary shape patchy saturation on wave characteristics and figured out the low and high-frequency limits for p-wave velocity.Additionally, with the growing importance of mesoscopic scale theory, upcoming researchers analyzed wave attenuation and velocity dispersion due to the presence of three-dimensional heterogeneities at mesoscopic scale [Müller et al. 2005, Toms et al. 2007].
In recent years, there has been a growing interest in analyzing the effects of mesoscopic heterogeneity on wave characteristics.During seismic and acoustic propagation through saturated porous rock with mesoscopic heterogeneity, significant wave attenuation and velocity dispersion can be produced due to WIFF [Muller et al. 2010] within seismic frequency band.In order to analyze the effects of meso-scopic fluid flow on wave characteristics, seismic modeling is proved as an effective tool which is frequently applied in recent studies concerning estimation of wave attenuation and velocity dispersion [Masson et al. 2006, Rubino et al. 2011, Rubino et al. 2012, Quintal et al. 2011, Milani et al. 2015].Masson and Pride [Masson et al. 2007] proposed an efficient method to estimate wave attenuation and velocity dispersion due to wave induced fluid flow (WIFF) at mesoscopic scale.According to their proposed quasi-static creep test, they used finite difference method in the time domain for the solution of Biot's quasi-static equations for wave propagation in poroelastic media.Furthermore, Rubino [Rubino et al. 2009] and Quintal [Quintal et al. 2011] have put forward two different quasi-static numerical strategies, in which they used the finite element method in frequency and time domain, for the computation of seismic attenuation due to wave-induced fluid flow.Exploiting the quasi-static finite element method in frequency domain, Carcione [Carcione et al. 2012] and Rubino [Quintal et al. 2014] have discussed the effects of fracture connectivity and anisotropy on the seismic attenuation and dispersion.Lately, the effects on seismic characteristics of fractured magmatic geothermal reservoirs were numerically modeled [Grab et al. 2017] in a simple and efficient way.
Despite the fact that, the application of seismic modeling to efficiently estimate wave attenuation and velocity dispersion due to wave induced fluid flow (WIFF) at mesoscopic scale has been improved in recent years.Nonetheless, it is possible to further improve the efficiency of the methods by put forwarding more flexible methods for accurate prediction of wave attenuation and velocity dispersion at mesoscopic scale.To achieve this goal, recent work seeks to propose a simple and more effective method to compute P-wave modulus for a wide range of frequencies compared with prior to the finite difference in time domain methods.
In this article, a new numerical technique is suggested in which finite difference in frequency domain is used to solve Biot's [Biot 1941] quasi-static equations of consolidation based on the finite difference in time domain approach, presented by Masson and pride [Masson et al. 2007].Our suggested methodology eventually improves the efficiency to estimate seismic wave attenuation and dispersion in a 1D saturated medium exhibiting mesoscopic heterogeneity.We further predicted that seismic wave attenuation is intensively sensitive to variation in patch size in a one-dimensional saturated porous medium.Both, the implementation and the computational procedure for measuring the complex P-wave modulus of our suggested methodology is much easier.At the same time, our proposed method can be further extended to two-dimensional case, in a flexible way.
The remainder of the paper is organized as follows: The outlines of relationships of 1D Biot's poroelastic theory are given first.Afterwards, the presentation of poroelastic relations of Biot's theory and the outlines of quasi-static finite difference in frequency domain are given, which will assist in computing wave attenuation and velocity dispersion in frequency domain.Then the implementation of proposed quasi-static numerical method in frequency domain and its comparison with analytical results of 1D White's model for rock saturated with two immiscible fluids is described.Then the numerically estimated results for different patch sizes having identical gas saturation are given.After that, the numerical measurements of 1D model saturated with three different kinds of fluids are presented.The end of the paper summarizes the results and draws the conclusions from our work.

1D Biot's poroelastic relationships
We have implemented a quasi-static finite difference modeling method in frequency domain for the estimation of frequency dependent complex P-wave modulus of a mesoscopic heterogeneous porous media.We have used the Biot's [Maurice A. Biot 1941] quasi-static equations for linear consolidation in a sat-urated porous media in 1D frequency domain, which can be written as (1) (2) (3) (4) In above relations τ zz and P defines total stress and fluid pressure, where U s and U f represents the solid and fluid displacement.The bulk modulus and shear modulus of the dry frame is indicated as k and m, where η, ω and k 0 represents viscosity, angular frequency and permeability of the saturated porous medium, , α and M are given as (5) In above equations, K s and K f represents the bulk modulus of the solid matrix and the fluid phase respectively φ represents the porosity of the media.

The outline of quasi-static finite difference in frequency domain
The aim of this study is to compute complex and frequency dependent P-wave modulus for a certain frequency band in case of a saturated porous media exhibiting mesoscopic heterogeneity.To achieve this goal, many researchers have proposed various methods, including finite element method in time/frequency domain and the finite difference method in time domain.To estimate wave attenuation and velocity dispersion due to WIFF at the mesoscopic scale, we have implemented traditional frequency domain finite difference method to compute complex and frequency dependent P-wave modulus.Then the boundary conditions are adjusted in such a manner that the applied stress positions, pressure, and solid/fluid displacement are remained fixed on the same grid point.In order to ensure the undrained condition, we have applied the normal stress P 0 on both sides of 1D synthetic sample and kept the fluid displacement as zero at both ends of the sample(Figure 1).By doing so, the solid displacement at both ends (U s top (ω) and U s bottom (ω)) of the sample can be measured easily.These measurements will assist in measuring the total strain (ω) in the synthetic sample.
After estimating the strain and knowing the applied stress, the frequency dependent complex P-wave modulus M (ω) can be easily computed by using equations (7-8).Eventually, seismic wave attenuation Q -1 (ω) and P-wave velocity V p (ω) can be measured by using equation [9-10) [Masson and Pride 2007]. ( Where L indicate the total length of the described model, is the effective density which can be obtained by using the formula, given below (11) Where ρ s and ρ f are the solid and fluid densities.Furthermore, the first-space derivative is obtained by applying second order space-operator of differencing approximation is given by ( 12)

1D White's model
We have proposed a quasi-static finite difference numerical method to quantify seismic wave attenuation and dispersion in frequency domain.To verify our quasistatic finite difference method, we implemented our   ε proposed method for the calculation of seismic wave attenuation and velocity dispersion of the White's 1D model (Figure 2a) exhibiting mesoscopic heterogeneity due to alternating layers of gas and water within a porous media.The statistical representative volume element (REV) enclosing a pair of alternating layers of saturated fluid, is given in Figure 2b.We implemented our proposed method on the described model having petrophysical properties given in table 1, where Dz = 0.1cm is the grid size and10cm as thickness of each layer.The results of total stress, pressure, solid and fluid displacement are given in Figure 3.
The comparison of our numerical results with the White's [White et al. 1975] 1D layer model's analytical results of seismic wave attenuation and dispersion is given in Figure 4.The results validate the accuracy of our proposed method over a wide frequency range (10-2-105).

Numerical experiments with different patch sizes for the same gas saturation
In order to compute seismic wave attenuation and velocity dispersion due to WIFF at mesoscopic scale, a novel technique is presented in this manuscript.The supposed caused of WIFF is due to the presence of alternating layers of gas and water within a porous media creating heterogeneity at mesoscopic scale.We further investigated the effect of variation in the patch size by keeping constant gas saturation.
The detail of our 8 different cases about the variation in patch size is given in Figure 5, the description of mesoscopic heterogeneity in first four cases indicate a gradual decrease in patch size where the random patch sizes are selected in other cases.The regions saturated with gas (brighter) and water (darker) are all divided into two parts in accordance with the size of 2 to 1   respectively.The low-frequency limit of Gassmann-Wood and high-frequency limit of Gassmann-Hill [Müller et al. 2010] are also calculated by using the physical properties of the solid and the saturating fluids given in Table 1.
Figure 6 shows the visual representation of frequency dependence of wave attenuation and velocity dispersion with a change in the size of the heterogeneity.The results for the case of a decrease in patch sizes shows that the velocity dispersion and attenuation peaks shifted towards higher frequencies with the decrease in size of heterogeneity.From case 5 to 8 patch sizes are randomly selected, there is a little difference in results of seismic wave attenuation and velocity dispersion.
Infect the results of rock sample of the case 5 and case 6 have almost identical values of wave attenuation and velocity dispersion, also the rock samples having patch size as shown in case 7 and case 8 have the same behavior.It can be concluded from results that, the change in the spatial size of the saturating fluids having lesser compressibility values have fewer effects on seismic wave attenuation and velocity dispersion than the fluids having higher compressibility values.

Application of quasi-static finite difference in frequency domain, on a 1D model, saturated with alternating layers of three different fluids
We designed a new quasi-static finite difference technique for the estimation of seismic wave attenuation and velocity dispersion in frequency domain.
Although several analytical/numerical studies have been carried out to estimate wave attenuation and velocity dispersion for rock saturated with two distinct fluid, little attention has been paid to predict wave attenuation and velocity dispersion for rock saturated with three different fluids.The problem with the analytical prediction of wave attenuation and velocity dispersion in rock saturated with three different fluids is that it is not easy to predict the complex P-wave modulus analytically.This part of the article demonstrates the feasibility of estimating the complex P-wave modulus and eventually the wave attenuation and velocity dispersion in rock saturated with three different fluids.
We implemented the proposed quasi-static finite difference modeling in frequency domain to a 1D model saturated with alternating layers of three different kinds of fluids (Figure 7).Numerical simulations of the investigated sample are demonstrated in Figure 8 and the estimated results of wave attenuation and velocity dispersion are visualized in Figure 9.According to our findings, this is the first study in which rock saturated with three different sorts of fluids is investigated and interestingly, there are two peaks in the attenuation curve at seismic frequencies and also the results are stable over a wide frequency range.The attenuation peak appears at the low frequency is because of the reason that at the interface between gas and oil the complex P-wave modulus gets increased whereas the peak appears at high frequency is indicating the increase in complex P-wave modulus at the interface between water and gas.Velocity dispersion at the same frequencies where there is high modulus dispersion validates the concept of higher values of complex P-wave modulus at the interfaces between different fluids.Because of little difference between modulus of oil and water, attenuation at the interface between oil and water is disappeared in the attenuation curve.

Effects of changing the order of saturating fluids
Presence of fluids within the subsurface rocks creates heterogeneity at a different scale which significantly influences the characteristics of propagating wave.Accurate computation of these influences on wave characteristics is an essential need of oil and gas exploration.
We have proposed a quasi-static finite difference method in frequency domain to compute the influence of three different fluids within a porous rock.Figure 10 gives the pictorial view of numerical simulation of investigated sample with change in fluid order (i.e.10a Gas, Oil, Water and 10b Gas, Water Oil) and the results of wave attenuation and velocity dispersion due to the change in the order of saturating fluids are demonstrated in Figure 11 (i.e.11a wave attenuation and 11b Velocity dispersion).
The black solid line shows the wave attenuation and velocity dispersion when the fluids are in order of Gas, Oil, and water while blue circle shows the order of Gas, Water, and Oil.It is interesting to notice that, there is a good match between Gas, Oil, Water and Gas, Water, Oil saturation.One possible reason for this situation could be, because of having little difference in compressibilities of water and oil, the pore pressure difference created by the passing wave in both cases is same.So the footprints of pore pressure equilibration at the interface between Gas/Oil and Gas/water occurs at the same frequencies.

Effects of changing the layer thickness of saturating fluids
For a visual representation of the dependence of layer thickness of three-phase fluids in a porous rock, the reader is referred to Figure 12.From the figure, it can be seen that with decreasing the layer thickness of saturating fluids, the phenomenon of pore pressure equilibration starts to occur at higher frequencies.The possible reason for the shift towards higher frequencies is that at a wavelength higher than layer thickness it's AHMAD ET AL. 8 easy for the wave to pass through the layer, while it gets trapped (attenuate) at higher frequencies when the layer thickness becomes equal to the wavelength.The result of velocity dispersion shows the same trend of shifting towards higher frequencies while the upper and lower velocity limits remain the same.

Effects on wave attenuation and velocity dispersion due to change in thickness of layer's having Gas saturation
We started by investigating the effects of mesoscopic heterogeneity (due to the presence of fluids) on wave propagation.Several numerical experiments were carried out to compute the influence of varying fluids properties, on wave propagation.In results of Figure 13, we analyzed the influence of a change in thickness of a gas layer on wave characteristics, by keeping the total gas saturation constant.The analysis of numerical findings indicates that variation in thickness of gas layer at top and bottom significantly influences wave characteristics.The increase in thickness of gas layer at the top cause increase in attenuation at low frequencies due to presence of highly compressible fluid (gas), while ANALYSIS OF SEISMIC RESPONSES FROM MICRO STRUCTURES   this attenuation peak gradually become flat with a decrease in thickness of gas layer at the top.Maximum attenuation obtained when the thickness of gas layer becomes identical at top and bottom of the representative element.

Conclusions
This study is related to estimation of wave attenuation and velocity dispersion in fluid-saturated porous media.Our work, suggest a new technique that improves the ability to estimate complex P-wave modulus and eventually wave attenuation and velocity dispersion.In our proposed technique, by solving Biot's quasi-static equations, a 1D finite difference scheme in frequency domain is used for the deduction of complex P-wave modulus of a mesoscopic heterogeneous porous media.In order to investigate the validity and accuracy of our proposed numerical scheme, we have compared our results with the analytical results of White's 1D model.Our results are in good agreement with the analytical results of White's 1D model exhibiting mesoscopic heterogeneity due to saturation of two distinct fluids in a porous media.Moreover, we also investigated the effect of, variation in mesoscopic patches, on wave attenuation and velocity dispersion.The findings of our research are quite convincing, and thus the following conclusions can be drawn: The proposed scheme of quasi-static finite difference in frequency domain, assist well in estimating more precise results of frequency dependent P-wave modulus over a wide frequency range (10-2-105).Attenuation and velocity dispersion is much more sensitive to change in the patch size of fluid having higher compressibility values.Also, attenuation peak shift towards higher frequency with the decrease in patch size.Furthermore, variation in patch size shows negligible effect on high and low-frequency limit in case of identical gas saturation.

Figure 1 .
Figure 1.The outline of quasi-static finite difference modeling in frequency domain.

Figure 6 .
Figure 6.Geometries of the partially saturated rock samples with different patchy sizes under the same condition of gas saturation, grey and black phases are saturated with gas and water respectively.

Figure 7 .
Figure 7. Rock sample saturated with three kinds of fluid.

Figure 8 .
Figure 8. Numerical simulations of the sample saturated with three kinds of fluid Stress (b) Fluid pressure (c) Solid displacement (d) Fluid displacement.

Figure 9 .
Figure 9. Seismic attenuation (a) and velocity dispersion (b) of the sample saturated with three kinds of fluid.(a) (b)

Figure 10 .
Figure 10.Numerical simulations of the sample saturated with three kinds of fluid in different order (a) Gas, Oil and Water (b) Gas, Water and Oil.

Figure 13 .
Figure 13.Seismic attenuation (a) and velocity dispersion (b) of the sample saturated with three kinds of fluid.

Figure 12 .
Figure 12.Seismic attenuation (a) and velocity dispersion (b) of the sample saturated with three kinds of fluid.

Figure 11 .
Figure 11.Seismic attenuation (a) and velocity dispersion (b) of the sample saturated with three kinds of fluid in different order.(a)(b)

Table 1 .
Properties of the skeleton and fluids.