Different forms of Gassmann equation and implications for the estimation of rock properties
Tran Trung Dung (1), Carl H.Sondergeld (2), Jean-Claude Roegiers (2) 1. Phu Quoc POC 2. University of Oklahoma Email: dungtt@phuquocpoc.vn

Summary

Three new equivalent forms of Gassmann equation are presented that are useful when the unknown parameters are the Biot- Willis coefficient, the dry bulk modulus, and/or the grain matrix bulk modulus. We apply these equations to several sets of laboratory measurements to determine the profiles of grain matrix bulk modulus and Biot-Willis coefficient as functions of applied pressure, and perform a Monte Carlo simulation to examine the effect of uncertainty and/or measurement errors on the calculated grain matrix bulk modulus and Biot-Willis coefficient. The results show that the calculated grain matrix bulk modulus is relatively constant with applied differential pressure (up to 50MPa). However, it is very sensitive to the uncertainty of dry and saturated bulk modulus values. Thus, the presented new forms of Gassmann equation can be used to effectively quantify the uncertainty of dry and saturated bulk modulus (and subsequently, the seismic velocities) in fluid identification, fluid substitution, or reservoir monitoring applications.

Key words: Gassmann equation, bulk moduli, Biot-Willis coefficient, sensitivity analysis.

1.       Introduction

The Gassmann equations [1] have been used extensively in the oil and gas industry for fluid identification and reservoir monitoring applications, despite its various assumptions [2, 3]. The first Gassmann equation provides the relationship between  the  saturated  bulk  modulus  of a rock and its dry frame bulk modulus, porosity, bulk modulus of the mineral matrix,  and  bulk  modulus  of  the pore-filling fluid. Whereas, the second Gassmann equation simply states that the shear modulus of the rock is independent of the presence of the saturating fluid:





Where α is the Biot-Willis coefficient [4]:

 

The moduli are related to the seismic velocities and density by:



Berryman [5] gave a concise derivation of Gassmann equations for an isotropic and homogeneous medium using   the   quasi-static   poroelastic   theory.   Other forms of Equation (1) can be found in Mavko et al. [6]. Zimmerman [7] presented an equivalent  form  in  terms of compressibility. However, Equation (1) is probably the most intuitive in describing the effect of fluid presence on the bulk modulus.

White and Castagna [8] argued that, since all input parameters for Gassmann equations carry some degrees of uncertainty, a fluid modulus inversion should be performed using a probabilistic approach. Artola and Alvarado [9] evaluated the effect of uncertainty of different input parameters and showed that the computed compressional velocity of a saturated rock is most sensitive to uncertainties in the rock bulk density and the dry bulk and shear moduli, while other parameters (porosity, the grain matrix and fluid bulk moduli) have negligible effects.

Note that the three parameters: dry frame modulus (Kdry), Biot-Willis coefficient (α), and grain matrix bulk modulus (Km) are related by Equation (3); in many instances they are unknown. The fluid saturated bulk modulus (Ksat) and fluid bulk modulus (Kf) can also be unknown (e.g. in fluid substitution problem). Asaresult, ad-hoc and empirical correlations have been proposed to address this problem. There are many instances Biot-Willis coefficient is assumedto be 1 due to the lack of a better estimate. For sandstone at high  differential  pressure  (40MPa),  Han and Batzle [10] proposed α to be a polynomial function of porosity:



In this study, we present 3 new equivalent forms of Gassmann equation that are useful for different scenarios of available data. We apply these equations to  several sets of laboratory measurements. A stochastic simulation is performed to examine the effect of uncertainty and/    or measurement errors on calculated grain matrix bulk modulus and Biot-Willis coefficient.

2.       The new equivalent Gassmann equations

When (K , K , K , and φ) are known

This is generally the case for laboratory measurements on dry and wet rock samples (e.g. dry and brine saturated acoustic velocities are measured as functions of differential pressure along with rock porosity). In this case we can rewrite Equation (1) as a function of Biot-Willis coefficient α (see Appendix A for a detailed derivation):



Equation (7) is a quadratic equation A2 + B + C = 0, where all coefficients can be readily calculated.



This simple quadratic equation has two solutions:

(11)

However,  Berryman and Milton(11) showed that α is physically bounded between 0 and 1. Equations (9) and (10) show that B is negative since Kdry < Ksat, and C is also negative since Kf < Kdry for consolidated rocks. This means

∆ = B2 – 4AC> B2 and thus, 1=…. Is the only.

  Therefore, instead of having two non-linear equations for two unknowns (α and Km), we reduce the problem to one simple quadratic equation, Equation (7), that always gives one physically realistic solution.

This provides an independent methodology to calculate the grainmatrixbul kmodulusofarockfrom SCAL laboratory acoustic  measurements of dry and saturated rock samples. Traditionally, the grain matrix bulk moduli are estimated from averages of the rock mineralogical composition (e.g. Voigt-Reuss-Hill average or Hashin- Shtrikman bounds). These bounds may carry large uncertainties since many minerals, especially clays, have a high variance in their bulk modulus values depending the measurement conditions [12, 13]. We can further  postulate that: (a) the grain matrix calculated from Gassmann equation (using Equations (8) to (12)) must lie between the two bounds obtained from mixture theory, and (b) the calculated grain matrix values are insensitive to the first order to the applied pressure. Equations (7) to (11) can also be used to verify the applicability of existing estimate Biot-Willis coefficient) for different rocks.

2.2. When (Ksat1, Ksat2, Kf1k Kf2, and ɸ)are known

This case can be encountered in the field. The same rock can be fully saturated with brine in one well; or it can have verying saturation in the same well. Acousstic logs, and density- porosity logs are available. In this case, K   , K  , and α  are unknown in a system of three non-linear equations (two Equation (1) for two different saturation fluids and  Equation (3)). Starting from Equation (7) instead, we end up with (see Appendix B for detaild derivations)

 (12)


(13)

We can write Equation (13) in a more convenient form for numerical calculations:

 possible solution since α2 is negative.

The corresponding grain matrix bulk modulus then can be calculated from Equation (3):



 
Kdry, α, and Km can now be calculated very quickly using a simple iteration using Equation (14) and Equation (7)    as follows:

-        Step 1: Make an initial guess, for example:



-        Step 2: Use guessed Kdry value in Equation (7) to find two Biot-Willis coefficients αf1 and αf2 (for two saturations):



-         Step 3: Take the average for a new Biot –Wills coefficient



-        Step 4: Use this new α in Equation (14) to find new Kdry.



-        Step 5: Repeat steps 2 to 4 until Kdry converges:



-        Step 6: Use Equations (7) and (12) to find corresponding α and Km.

Note that we have assumed there are no softening or hardening effects caused by the saturating fluids on the grain bulk modulus (Km is constant). The secondassumption  is  that  the  rock  dry  frame  is  stiffer  than both fluids, Kdry > max { Kf 1, Kf 2 }, so that Equation (7) still gives only one positive (physically realistic) root. This assumption is generally valid for consolidated sedimentary rocks.

2.3. When (Km, Ksat, Kf, and φ) are known

In this case Kdry and α are unknown. An instance for this case is that fluid data, fluid saturation, acoustic log (Vp, Vs) and density-porosity log are available while Km is estimated from the mineralogical composition of the rock (FTIR, XRD, thin section of rock cuttings, or mineralogy log). The Biot-Willis coefficient can be estimated directly from the following equivalent Gassmann equation (see Appendix C for detailed derivations):

(15)

 
And applying this α to Equation (3) gives Kdry. This is equivalent to the Kdry solution of Zhu and McMechan [14].

3.       Numerical applications

          Han and Batzle’s sandstone data

We applied Equation (7) to the pressure dependent sandstone sample published by Han and Batzle [10] (Figure using the density relaρtionshρip:  

 (16)
The calculated Biot-Willis coefficient and grain matrix modulus as functions of pressure are plotted in Figure 1. The relatively constant value of the grain bulk modulus (39GPa) as a function of pressure is a good indicator that Gassmann equation is applicable for this rock. The variation of grain bulk modulus at low confining pressure (< 10MPa) is possibly due to higher uncertainty in input values (i.e. higher noise-to-signal ratio from velocity signals).

The Biot-Willis coefficient profile is remarkably similar to the result measured on a 26% porosity Boise sandstone sample by Fatt  [15]. Note that Gassmann equation gave   a higher value (0.73 at 40MPa) than Han and Batzle’s Equation (6) (0.63).

          Coyner’s limestone data

We employed the iteration procedure using Equations (7) to (14) on water and benzene saturated Bedford limestone sample published by  Coyner  [16]  (Figure  2).  In his experiment at room temperature, the fluid pore pressure in both saturation cases was maintained at 10MPa. The porosity of the rock is 11.9%. The shear modulus  profile  is  almost  identical  for  all  vacuum  dry, water saturated, and benzene saturated cases, suggesting

The back-calculated dry bulk modulus is also plotted (as a dashed line) against the various measured moduli in Figure 2. The profile is consistently higher than the measured vacuum dry bulk modulus profile by approximately 2.5GPa (or 5 - 9%). This is another evidence supporting the argument that the vacuum dry measured bulk modulus is too dry and should not be used in Gassmann equation [6, 17]. We could have applied the measured vacuum dry and either water- or benzene-saturated moduli values on Equation (7), but  that approach would give unrealistically high grain matrix bulk modulus.

 

Figure 1. Grain bulk modulus and Biot-Willis coefficient of a sandstone sample as a func- tion of pressure calculated from its dry and brine saturated moduli [10] using Gassmann equation.

 

Figure 2. Bulk and shear moduli as a function of differential pressure for Bedford lime- stone [6]. The dashed line is the dry bulk modulus calculated from Gassmann equation, consistently higher than the vacuum dry measured data.

 

Figure 3. Calculated grain matrix bulk modulus and Biot-Willis coefficient of Bedford limestone sample as a function of pressure using Gassmann equation from its water - and benzene - saturated moduli [16].

The grain matrix bulk modulus and Biot-Willis coefficient profiles obtained from the rock water- and benzene-saturated moduli are shown in Figure 3. While the grain matrix bulk modulus is similar to Coyner’s reported value of 65GPa, the back-calculated Biot-Willis coefficient profile decreases from 0.6 to 0.53, significantly lower than the commonly assumed value of 1 while significantly higher than the estimated value of 0.34 obtained from Han and Batzle’s 2004 correlation.

          Effects of input data errors on calculated grain bulk modulus

Measured values always have some associated errors. Velocities, especially shear wave velocities may carry significant uncertainties. We would like to determine the effects of uncertainties from porosity, Kdry, Ksat, and Kf to the uncertainty of the predicted Km. Since the relationship in Equation (7) is not linear, a Monte Carlo (stochastic) simulation was used.

Table    1   summarises   the   input   parameter  values [18] and their estimated ranges of  uncertainties.  The  rock sample is a Berea sandstone sample with Voigt- Reuss-Hill average grain bulk modulus  of  39.6GPa from its mineralogical composition. All parameters were assumed to have a normal distribution with means being the measured values and the errors represent the 95% confidence interval. Thus, the relative error (uncertainty) of each parameter is defined as:

 (17)

 
Where s is the standard deviation of the parameter’s sample.

For each set of perturbed errors, 10,000 sets of (porosity, dry bulk modulus, wet bulk modulus, and fluid modulus) values were generated to compute 10,000 grain bulk moduli, which are then analysed for the mean value and standard deviation.

In our base case, porosity is assigned a 1% error,  Kdry and Ksat are each assigned a 3% error,  and Kf carries  a 10% uncertainty. The resulting Km is also a Gaussian distribution with a mean of 44.6GPa and a standard deviation of 3.45GPa. The 95% confidence interval is, therefore, from 37.7GPa to 51.5GPa (or 16% error). The Biot-Willis coefficient α is also a Gaussian distribution with a mean of 0.62 and a standard deviation of 0.03. The 95% confidence interval is from 0.56 to 0.68 (or 10% error).

 

Table 1. Mean (measured) values of a Berea sandstone sample [18] and ranges of uncertainties used in Monte Carlo simulations

 

Figure 4. The uncertainty of the computed grain matrix bulk modulus Kmusing Gassmann equation as functions of percent error in one input parameter (Kf, α, Kdry, or Ksat), while the remaining input parameters carry the same errors as the base case. Errors from Ksat and Kdry have the largest impacts on the uncertainty of calculated Km. Porosity and fluid bulk modulus, on the other hand, show negligible effects.

Figure 4 shows the uncertainty of the computed grain matrix bulk modulus Km as functions of percent error in one input parameter (K , φ, K , or K ), while the remaining input parameters carry the same uncertainties as of the base case. Errors from Ksat and Kdry have the largest effects on the uncertainty of Km. Minus errors in Kdry and Ksat (even within laboratory measurement standard) can result in a large error in the estimated value of Km. Porosity and fluid bulk modulus, on the other hand, show negligible influence. This result is not surprising, as Kf, φ and Km should be uncorrelated parameters.

4.       Conclusions

Three equivalent forms of Gassmann equation were presented that can be useful for the determination of Biot-Willis  coefficient,  dry  bulk  modulus,  and/or  grain matrix bulk modulus of a rock. We demonstrated the applicability of these equations using several sets of published laboratory measurements and the implications of the results for other estimations of rock properties. Astochastic simulation was also performed  to  examine  the effect of uncertainty and/or measurement errors on calculated grain matrix bulk modulus. The results showed that the calculated grain matrix bulk modulus is relatively constant with applied differential pressure (up to 50MPa) for sedimentary rocks. However, the estimation is very sensitive to the uncertainty of dry and saturated bulk modulus values. Our new forms of Gassmann equation can be used to effectively quantify  the  uncertainty  of  dry and saturated bulk modulus (and subsequently, the seismic velocities) in fluid identification, fluid substitution, or reservoir monitoring applications.

NOMENCLATURE

K: bulk modulus (GPa or psi)

 Ksat: saturated bulk modulus (GPa or psi)

 Kdry: dry (frame) bulk modulus (GPa or psi) Km: grain (matrix) bulk modulus (GPa or psi) Kf: fluid bulk modulus (GPa or psi)

G: shear modulus (GPa or psi)

Gsat: saturated shear modulus (GPa or psi)

Gdry: dry (frame) shear modulus (GPa or psi)

α: Biot-Willis coefficient (dimensionless)

φ: porosity (dimensionless)

ƿ: density (g/cc)

Vp: compressional wave velocity (km/s)

Vs: shear wave velocity (km/s)

APPENDIX A: Derivation of Equation (7)

From Equation (3) we can write:
(A.1)

 Rewriting Equation (1) as a function of α gives:

(A2)

(A3)

Multipplying both side with



 which is Equation (7).

APPENDIX B: Derivation of Equation (13)

If the same rock is subjected to two different saturation fluids, then we have two equations in the form of Equation (7):



equation and rearrange Equation (1) as a function of α, the Biot-Willis coefficient only:



  which is Equation (15).

APPENDIX C: Derivation of Equation (15)

If Km  value can  be obtained (e.g using mixture theory), then one can substitute K = (1 − α) K into Gassmann sandstones. Elsevier Science. 1991.


 References

1.       Fritz Gassmann. Uber die Elastizität Poröser Medien. Vierteljahrschrift der Naturforschenden Gesellschaftin Zürich. 1951; 96: p. 1 - 23.

2.       Tad M.Smith, Carl H.Sondergeld, Chandra S.Rai. Gassmannfluid substitution. A tutorial. Geophysics. 2003; 68(2): p. 430 - 440.

3. Ludmila Adam, Michael Batzle, Ivar Brevik.Gassmannfluid substitution and shear modulus variabilityin  Geophysics. 2006; 71(6): p. F173 - F183.

4.Murice Anthony Biot, David G.Willis. The elastic coefficients of the theory of consolidation. Journal of Applied Machanic. 1957; 24: p.594-601.

5.       James G.Berryman. Origin of Gassmann’s equations. Geophysics. 1999; 64(5): p. 1627 - 1629.

6.       Gary Mavko, Tapan Mukerji, Jack Dvorkin. The rock physics handbook: Tools for seismic analysis in porous media. Cambridge University Press, Cambridge. 1998.

7.       Robert        W.Zimmerman.    Compressibility    of which is Equation (13). Equation (14) then can be readily obtained by multiplying both sides by (Ksat1 × Ksat2).


8.       Luther White, John Castagna. Stochastic fluid modulus inversion. Geophysics. 2002; 67(6): p. 1835 - 1843.

9.       Fredy A.V.Artola, Vladimir Alvarado. Sensitivity analysis of Gassmann's fluid substitution equations: Some implications in feasibility studies of time-lapse seismic reservoir monitoring. Jounal of Applied Geophysics. 2006; 59(1): p. 47-62.

10.     De-Hua Han, Michael L.Batzle. Gassmann’s equation and fluid-saturation effects on seismic velocities. Geophysics. 2004; 69(2): p. 398 - 405.

11.     James G.Berryman, Graeme W.Milton. Exact results for generalized Gassmann’s equations in composite porous media with two constituents. Geophysics. 1991; 56: p. 1950 - 1960.

12.     Keith W.Katahara. Clay mineral elastic properties. SEG Technical Program Expanded Abstracts. 1996: p. 1691 - 1694.

13.     Zhijing Jee Wang, Hui Wang, Michael E.Cates. Elastic properties of solid clays. SEG Technical Program Expanded Abstracts. 1998: p. 1045 - 1048.

14.     Xianhuai Zhu, George A.McMechan. Direct estimation of the bulk modulu softheframeinfluidsaturated elastic medium by Biot theory. SEG Technical Program Expanded Abstracts. 1990: p. 787 - 790.

15.    I.Fatt. The Biot-Willis elastic coefficients for a sandstone. Journal of Applied Mechanics. 1959; 26: p. 296 - 297.

16.     Karl B.Coyner. Effects of stress, pore pressure, and pore fluids on bulk strain, velocity, and permeability in rocks. Ph.D. dissertation, Massachusetts Institute of Technology, Cambridge, Massachusetts. 1984.

17.     Virginia    A.Clark,     Bernhard     R.Tittmann, Terry  W.Spencer.  Effect  of  volatiles  on   attenuation  (Q-1)   and   velocity   in    sedimentary    rocks.    Journal    of        Geophysical        Research.        1980;         85(B10):  p. 5190 - 5198.

18.     Tran Trung Dung, Chandra S.Rai, Carl H.Sondergeld. Changes in crack aspect-ratio concentration from heat treatment: A comparison between velocity inversion and experimental data. Geophysics. 2008; 73(4): p. E123 - E132.