bruges.rockphysics package#

Submodules#

bruges.rockphysics.anisotropy module#

bruges.rockphysics.anisotropy.backus(vp, vs, rho, lb, dz)[source]#

Backus averaging. Using Liner’s algorithm (2014; see Notes).

Parameters
  • vp (ndarray) – P-wave interval velocity.

  • vs (ndarray) – S-wave interval velocity.

  • rho (ndarray) – Bulk density.

  • lb (float) – The Backus averaging length in m.

  • dz (float) – The depth sample interval in m.

Returns

the smoothed logs: vp, vs, plus rho. Useful for computing

other elastic parameters at a seismic scale.

Return type

namedtuple

Notes

Liner, C (2014), Long-wave elastic attenuation produced by horizontal layering. The Leading Edge, June 2014, p 634-638.

bruges.rockphysics.anisotropy.backus_parameters(vp, vs, rho, lb, dz)[source]#

Intermediate parameters for Backus averaging. This is expected to be a private function. You probably want backus() and not this.

Parameters
  • vp (ndarray) – P-wave interval velocity.

  • vs (ndarray) – S-wave interval velocity.

  • rho (ndarray) – Bulk density.

  • lb (float) – The Backus averaging length in m.

  • dz (float) – The depth sample interval in m.

Returns

Liner’s 5 intermediate parameters: A, C, F, L and M.

Return type

tuple

Notes

Liner, C (2014), Long-wave elastic attenuation produced by horizontal layering. The Leading Edge, June 2014, p 634-638.

bruges.rockphysics.anisotropy.backus_quality_factor(vp, vs, rho, lb, dz)[source]#

Compute Qp and Qs from Liner (2014) equation 10.

Parameters
  • vp (ndarray) – P-wave interval velocity.

  • vs (ndarray) – S-wave interval velocity.

  • rho (ndarray) – Bulk density.

  • lb (float) – The Backus averaging length in m.

  • dz (float) – The depth sample interval in m.

Returns

Qp and Qs.

Return type

namedtuple

bruges.rockphysics.anisotropy.blangy(vp1, vs1, rho1, d1, e1, vp0, vs0, rho0, d0, e0, theta)[source]#

Blangy, JP, 1994, AVO in transversely isotropic media-An overview. Geophysics 59 (5), 775-781. DOI: 10.1190/1.1443635

Provide Vp, Vs, rho, delta, epsilon for the upper and lower intervals, and theta, the incidence angle.

Parameters
  • vp1 – The p-wave velocity of the upper medium.

  • vs1 – The s-wave velocity of the upper medium.

  • rho1 – The density of the upper medium.

  • d1 – Thomsen’s delta of the upper medium.

  • e1 – Thomsen’s epsilon of the upper medium.

  • vp0 – The p-wave velocity of the lower medium.

  • vs0 – The s-wave velocity of the lower medium.

  • rho0 – The density of the lower medium.

  • d0 – Thomsen’s delta of the lower medium.

  • e0 – Thomsen’s epsilon of the lower medium.

  • theta – A scalar [degrees].

Returns

the isotropic and anisotropic reflectivities in a tuple. The isotropic result is equivalent to Aki-Richards.

TODO

Use rocks.

bruges.rockphysics.anisotropy.crack_density(porosity, aspect)[source]#

Returns crack density from porosity and aspect ratio, phi and alpha respectively in the unnumbered equation between 15.40 and 15.41 in Dvorkin et al. 2014.

Parameters
  • porosity (float) – Fractional porosity.

  • aspect (float) – Aspect ratio.

Returns

Crack density.

Return type

float

bruges.rockphysics.anisotropy.dispersion_parameter(qp)[source]#

Kjartansson (1979). Journal of Geophysical Research, 84 (B9), 4737-4748. DOI: 10.1029/JB084iB09p04737.

bruges.rockphysics.anisotropy.hudson_delta_G(porosity, aspect, mu, lam=None, pmod=None)[source]#

The approximate reduction in shear modulus G (or mu) in the direction normal to a set of aligned cracks. Eqn 15.42 in Dvorkin et al (2014).

Parameters
  • porosity (float) – Fractional porosity, phi.

  • aspect (float) – Aspect ratio, alpha.

  • mu (float) – Shear modulus, sometimes called G.

  • lam (float) – Lame’s first parameter, lambda.

  • pmod (float) – Compressional modulus, M.

Returns

M_inf - M_0 = Delta c_11.

Return type

float

bruges.rockphysics.anisotropy.hudson_delta_M(porosity, aspect, mu, lam=None, pmod=None)[source]#

The approximate reduction in compressional modulus M in the direction normal to a set of aligned cracks. Eqn 15.40 in Dvorkin et al (2014).

Parameters
  • porosity (float) – Fractional porosity, phi.

  • aspect (float) – Aspect ratio, alpha.

  • mu (float) – Shear modulus, sometimes called G.

  • lam (float) – Lame’s first parameter.

  • pmod (float) – Compressional modulus, M.

Returns

M_inf - M_0 = Delta c_11.

Return type

float

bruges.rockphysics.anisotropy.hudson_inverse_Q_ratio(mu=None, pmod=None, pr=None, vp=None, vs=None, aligned=True)[source]#

Dvorkin et al. (2014), Eq 15.44 (aligned) and 15.48 (not aligned). You must provide one of the following: pr, or vp and vs, or mu and pmod.

Parameters
  • mu (float) – Shear modulus, sometimes called G.

  • pmod (float) – Compressional modulus, M.

  • pr (float) – Poisson’s ratio, somtimes called v.

  • vp (ndarray) – P-wave interval velocity.

  • vs (ndarray) – S-wave interval velocity.

  • aligned (bool) – Either treats cracks as alligned (Default, True) or assumes defects are randomly oriented (False)

Returns

2Q_s^-1

Return type

float

bruges.rockphysics.anisotropy.hudson_quality_factor(porosity, aspect, mu, lam=None, pmod=None)[source]#

Returns Q_p and Q_s for cracked media. Equations 15.41 and 15.43 in Dvorkin et al. (2014).

Parameters
  • porosity (float) – Fractional porosity, phi.

  • aspect (float) – Aspect ratio, alpha.

  • mu (float) – Shear modulus, sometimes called G.

  • lam (float) – Lame’s first parameter, lambda.

  • pmod (float) – Compressional modulus, M.

Returns

Q_p float: Q_s

Return type

float

bruges.rockphysics.anisotropy.ruger(vp1, vs1, rho1, d1, e1, vp2, vs2, rho2, d2, e2, theta)[source]#

Coded by Alessandro Amato del Monte and (c) 2016 by him https://github.com/aadm/avo_explorer/blob/master/avo_explorer_v2.ipynb

Rüger, A., 1997, P -wave reflection coefficients for transversely isotropic models with vertical and horizontal axis of symmetry: Geophysics, v. 62, no. 3, p. 713–722.

Provide Vp, Vs, rho, delta, epsilon for the upper and lower intervals, and theta, the incidence angle.

Parameters
  • vp1 – The p-wave velocity of the upper medium.

  • vs1 – The s-wave velocity of the upper medium.

  • rho1 – The density of the upper medium.

  • d1 – Thomsen’s delta of the upper medium.

  • e1 – Thomsen’s epsilon of the upper medium.

  • vp0 – The p-wave velocity of the lower medium.

  • vs0 – The s-wave velocity of the lower medium.

  • rho0 – The density of the lower medium.

  • d0 – Thomsen’s delta of the lower medium.

  • e0 – Thomsen’s epsilon of the lower medium.

  • theta – A scalar [degrees].

Returns

anisotropic reflectivity.

bruges.rockphysics.anisotropy.thomsen_parameters(vp, vs, rho, lb, dz)[source]#

Liner, C, and T Fei (2006). Layer-induced seismic anisotropy from full-wave sonic logs: Theory, application, and validation. Geophysics 71 (6), p D183–D190. DOI:10.1190/1.2356997

Parameters
  • vp (ndarray) – P-wave interval velocity.

  • vs (ndarray) – S-wave interval velocity.

  • rho (ndarray) – Bulk density.

  • lb (float) – The Backus averaging length in m.

  • dz (float) – The depth sample interval in m.

Returns

delta, epsilon and gamma.

Return type

namedtuple

bruges.rockphysics.bounds module#

bruges.rockphysics.bounds.hashin_shtrikman(f, k, mu, modulus='bulk')[source]#

Hashin-Shtrikman bounds for a mixture of two constituents. The best bounds for an isotropic elastic mixture, which give the narrowest possible range of elastic modulus without specifying anything about the geometries of the constituents.

Parameters
  • f – list or array of volume fractions (must sum to 1.00 or 100%).

  • k – bulk modulus of constituents (list or array).

  • mu – shear modulus of constituents (list or array).

  • modulus – A string specifying whether to return either the ‘bulk’ or ‘shear’ HS bound.

Returns

The Hashin Shtrikman (lower, upper) bounds.

Return type

namedtuple

Source

Berryman, J.G., 1993, Mixture theories for rock properties Mavko, G., 1993, Rock Physics Formulas.

: Written originally by Xingzhou ‘Frank’ Liu, in MATLAB : modified by Isao Takahashi, 4/27/99, : Translated into Python by Evan Bianco

bruges.rockphysics.bounds.hill_average(f, m)[source]#

The Hill average effective elastic modulus, mh of a mixture of N material phases. This is defined as the simple average of the Reuss (lower) and Voigt (upper) bounds.

Parameters
  • f – list or array of N volume fractions (must sum to 1 or 100).

  • m – elastic modulus of N constituents (list or array).

Returns

Hill average.

Return type

mh

bruges.rockphysics.bounds.reuss_bound(f, m)[source]#

The lower bound on the effective elastic modulus of a mixture of N material phases. This is defined at the harmonic average of the constituents. Same as Wood’s equation for homogeneous mixed fluids.

Parameters
  • f – list or array of N volume fractions (must sum to 1 or 100).

  • m – elastic modulus of N constituents (list or array).

Returns

Reuss lower bound.

Return type

mr

bruges.rockphysics.bounds.voigt_bound(f, m)[source]#

The upper bound on the effective elastic modulus, mv of a mixture of N material phases. This is defined at the arithmetic average of the constituents.

Parameters
  • f – list or array of N volume fractions (must sum to 1 or 100).

  • m – elastic modulus of N constituents (list or array).

Returns

Voigt upper bound.

Return type

mv

bruges.rockphysics.elastic module#

bruges.rockphysics.elastic.elastic_impedance(vp, vs, rho, theta1, k=None, normalize=False, constants=None, use_sin=False, rho_term=False)[source]#

Returns the elastic impedance as defined by Connolly, 1999; we are using the formulation reported in Whitcombe et al. (2001). Inputs should be shape m x 1, angles should be n x 1. The result will be m x n.

Parameters
  • vp (ndarray) – The P-wave velocity scalar or 1D array.

  • vs (ndarray) – The S-wave velocity scalar or 1D array.

  • rho (ndarray) – The bulk density scalar or 1D array.

  • theta1 (ndarray) – Incident angle(s) in degrees, scalar or array.

  • k (float) – A constant, see Connolly (1999). Default is None, which computes it from Vp and Vs.

  • normalize (bool) – if True, returns the normalized form of Whitcombe.

  • constants (tuple) – A sequence of 3 constants to use for normalization. If you don’t provide this, then normalization just uses the means of the inputs. If these are scalars, the result will be the acoustic impedance (see Whitcombe et al., 2001).

  • use_sin (bool) – If True, use sin^2 for the first term, instead of tan^2 (see Connolly).

  • rho_term (bool) – If True, provide alternative form, with Vp factored out; use in place of density in generating synthetics in other software (see Connolly). In other words, the result can be multipled with Vp to get the elastic impedance.

Returns

The elastic impedance log at the specficied angle or angles.

Return type

ndarray

bruges.rockphysics.fluids module#

bruges.rockphysics.fluids.rho_brine(temperature, pressure, salinity)[source]#

The density of NaCl brine, given temperature, pressure, and salinity. The density of pure water is computed from rho_water(). Implements equation 27b from Batzle and Wang (1992).

Use scalars or arrays; if you use arrays, they must be the same size.

Parameters
  • temperature (array) – The temperature in degrees Celsius.

  • pressure (array) – The pressure in pascals.

  • salinity (array) – The weight fraction of NaCl, e.g. 35e-3 for 35 parts per thousand, or 3.5% (the salinity of seawater).

Returns

The density in kg/m3.

Return type

array

bruges.rockphysics.fluids.rho_gas(temperature, pressure, molecular_weight)[source]#

Gas density, given temperature (in deg C), pressure (in Pa), and molecular weight.

Parameters
  • temperature (array) – The temperature in degrees Celsius.

  • pressure (array) – The pressure in pascals.

  • molecular_weight (array) – The molecular weight.

Returns

The density in kg/m3.

Return type

array

bruges.rockphysics.fluids.rho_water(temperature, pressure)[source]#

The density of pure water, as a function of temperature and pressure. Implements equation 27a from Batzle and Wang (1992).

Use scalars or arrays; if you use arrays, they must be the same size.

Parameters
  • temperature (array) – The temperature in degrees Celsius.

  • pressure (array) – The pressure in pascals.

Returns

The density in kg/m3.

Return type

array

bruges.rockphysics.fluids.v_brine(temperature, pressure, salinity)[source]#

The acoustic velocity of brine, as a function of temperature (deg C), pressure (Pa), and salinity (weight fraction). Implements equation 29 from Batzle and Wang (1992).

Note that this function does not work at pressures above about 100 MPa.

Use scalars or arrays; if you use arrays, they must be the same size.

Parameters
  • temperature (array) – The temperature in degrees Celsius.

  • pressure (array) – The pressure in pascals.

  • salinity (array) – The weight fraction of NaCl, e.g. 35e-3 for 35 parts per thousand, or 3.5% (the salinity of seawater).

Returns

The velocity in m/s.

Return type

array

bruges.rockphysics.fluids.v_water(temperature, pressure)[source]#

The acoustic velocity of pure water, as a function of temperature and pressure. Implements equation 28 from Batzle and Wang (1992), using the coefficients in Table 1.

Note that this function does not work at pressures above about 100 MPa.

Use scalars or arrays; if you use arrays, they must be the same size.

Parameters
  • temperature (array) – The temperature in degrees Celsius.

  • pressure (array) – The pressure in pascals.

Returns

The velocity in m/s.

Return type

array

bruges.rockphysics.fluids.wood(Kf1, Kf2, Sf1)[source]#

Wood’s equation, per equation 35b in Batzle and Wang (1992).

bruges.rockphysics.fluidsub module#

bruges.rockphysics.fluidsub.avseth_fluidsub(vp, vs, rho, phi, rhof1, rhof2, kmin, kf1, kf2)[source]#

Naive fluid substitution from Avseth et al, section 1.3.1. Bulk modulus of fluid 1 and fluid 2 are already known, and the bulk modulus of the dry frame, Kmin, is known. Use SI units.

Parameters
  • vp (float) – P-wave velocity

  • vs (float) – S-wave velocity

  • rho (float) – bulk density

  • phi (float) – porosity (volume fraction, e.g. 0.20)

  • rhof1 (float) – bulk density of original fluid (base case)

  • rhof2 (float) – bulk density of substitute fluid (subbed case)

  • kmin (float) – bulk modulus of solid mineral (s)

  • kf1 (float) – bulk modulus of original fluid

  • kf2 (float) – bulk modulus of substitute fluid

Returns

Vp, Vs, and rho for the substituted case

Return type

Tuple

bruges.rockphysics.fluidsub.avseth_gassmann(ksat1, kf1, kf2, k0, phi)[source]#

Applies the inverse Gassmann’s equation to calculate the rock bulk modulus saturated with fluid with bulk modulus. Requires as inputs the insitu fluid bulk modulus, the insitu saturated rock bulk modulus, the mineral matrix bulk modulus and the porosity.

Parameters
  • ksat1 (float) – initial saturated rock bulk modulus.

  • kf1 (float) – initial fluid bulk modulus.

  • kf2 (float) – final fluid bulk modulus.

  • k0 (float) – mineral bulk modulus.

  • phi (float) – porosity.

Returns

final saturated rock bulk modulus.

Return type

ksat2 (float)

bruges.rockphysics.fluidsub.rhogas(gravity, temp, pressure)[source]#

From http://www.spgindia.org/geohorizon/jan_2006/dhananjay_paper.pdf

bruges.rockphysics.fluidsub.rhosat(phi, sw, rhomin, rhow, rhohc)[source]#

Density of partially saturated rock.

bruges.rockphysics.fluidsub.smith_fluidsub(vp, vs, rho, phi, rhow, rhohc, sw, swnew, kw, khc, kclay, kqtz, vclay, rhownew=None, rhohcnew=None, kwnew=None, khcnew=None)[source]#

Naive fluid substitution from Smith et al. 2003. No pressure/temperature correction. Only works for SI units right now.

Parameters
  • vp (float) – P-wave velocity

  • vs (float) – S-wave velocity

  • rho (float) – bulk density

  • phi (float) – porosity (v/v, fraction)

  • rhow (float) – density of water

  • rhohc (float) – density of HC

  • sw (float) – water saturation in base case

  • swnew (float) – water saturation in subbed case

  • kw (float) – bulk modulus of water

  • khc (float) – bulk modulus of hydrocarbon

  • kclay (float) – bulk modulus of clay

  • kqtz (float) – bulk modulus of quartz

  • vclay (float) – Vclay, v/v

  • rhownew (float) – density of water in subbed case (optional)

  • rhohcnew (float) – density of HC in subbed case (optional)

  • kwnew (float) – bulk modulus of water in subbed case (optional)

  • khcnew (float) – bulk modulus of hydrocarbon in subbed case (optional)

Returns

Vp, Vs, and rho for the substituted case.

Return type

Tuple

bruges.rockphysics.fluidsub.smith_gassmann(kdry, k0, kf, phi)[source]#

Applies the direct Gassmann’s equation to calculate the saturated rock bulk modulus from porosity and the dry-rock, mineral and fluid bulk moduli.

Parameters
  • kdry (float) – dry-rock bulk modulus.

  • k0 (float) – mineral bulk modulus.

  • kf (float) – fluid bulk modulus.

  • phi (float) – Porosity.

Returns

saturated rock bulk modulus.

Return type

ksat (float)

bruges.rockphysics.fluidsub.vels(Kdry, Gdry, K0, D0, Kf, Df, phi)[source]#

Calculate velocities and densities of saturated rock via Gassmann equation. Provide all quantities in SI units.

Parameters
  • Kdry (float or array_like) – Dry rock bulk & shear modulus in Pa.

  • Gdry (float or array_like) – Dry rock bulk & shear modulus in Pa.

  • K0 (float or array_like) – Mineral bulk & shear modulus in Pa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in Pa.

  • Kf (float or array_like) – Fluid bulk modulus in Pa and density in kg/m^3.

  • Df (float or array_like) – Fluid bulk modulus in Pa and density in kg/m^3.

  • phi (float or array_like) – Porosity, v/v.

Returns

  • vp, vs (float or array_like) – Saturated rock P- and S-wave velocities in m/s.

  • rho (float or array_like) – Saturated rock density in kg/m^3.

  • K (float or array_like) – Saturated rock bulk modulus in Pa.

bruges.rockphysics.fluidsub.vrh(kclay, kqtz, vclay)[source]#

Voigt-Reuss-Hill average to find Kmatrix from clay and qtz components. Works for any two components, they don’t have to be clay and quartz.

From Smith et al, Geophysics 68(2), 2003.

Parameters
  • kclay (float) – K_clay.

  • kqtz (float) – K_quartz.

  • vclay (float) – V_clay.

Returns

Kvrh, also known as Kmatrix.

bruges.rockphysics.moduli module#

bruges.rockphysics.moduli.bulk(vp=None, vs=None, rho=None, mu=None, lam=None, youngs=None, pr=None, pmod=None)[source]#

Computes bulk modulus given either Vp, Vs, and rho, or any two elastic moduli (e.g. lambda and mu, or Young’s and P moduli). SI units only.

Parameters
  • vp

  • vs

  • rho (and) –

  • lam (or any 2 from) –

  • mu

  • youngs

  • pr

  • pmod (and) –

Returns

Bulk modulus in pascals, Pa

bruges.rockphysics.moduli.lam(vp=None, vs=None, rho=None, pr=None, mu=None, youngs=None, bulk=None, pmod=None)[source]#

Computes lambda given either Vp, Vs, and rho, or any two elastic moduli (e.g. bulk and mu, or Young’s and P moduli). SI units only.

Parameters
  • vp

  • vs

  • rho (and) –

  • bulk (or any 2 from) –

  • mu

  • youngs

  • pr

  • pmod (and) –

Returns

Lambda in pascals, Pa

bruges.rockphysics.moduli.moduli_dict(vp, vs, rho)[source]#

Computes elastic moduli given Vp, Vs, and rho. SI units only.

Parameters
  • Vp

  • Vs

  • rho (and) –

Returns

A dict of elastic moduli, plus P-wave impedance.

bruges.rockphysics.moduli.mu(vp=None, vs=None, rho=None, pr=None, lam=None, youngs=None, bulk=None, pmod=None)[source]#

Computes shear modulus given either Vp, Vs, and rho, or any two elastic moduli (e.g. lambda and bulk, or Young’s and P moduli). SI units only.

Parameters
  • vp

  • vs

  • rho (and) –

  • lam (or any 2 from) –

  • bulk

  • youngs

  • pr

  • pmod (and) –

Returns

Shear modulus in pascals, Pa

bruges.rockphysics.moduli.pmod(vp=None, vs=None, rho=None, pr=None, mu=None, lam=None, youngs=None, bulk=None)[source]#

Computes P-wave modulus given either Vp, Vs, and rho, or any two elastic moduli (e.g. lambda and mu, or Young’s and bulk moduli). SI units only.

Parameters
  • vp

  • vs

  • rho (and) –

  • lam (or any 2 from) –

  • mu

  • youngs

  • pr

  • bulk (and) –

Returns

P-wave modulus in pascals, Pa

bruges.rockphysics.moduli.pr(vp=None, vs=None, rho=None, mu=None, lam=None, youngs=None, bulk=None, pmod=None)[source]#

Computes Poisson ratio given either Vp, Vs, and rho, or any two elastic moduli (e.g. lambda and mu, or Young’s and P moduli). SI units only.

Parameters
  • vp

  • vs

  • rho (and) –

  • lam (or any 2 from) –

  • mu

  • youngs

  • bulk

  • pmod (and) –

Returns

Poisson’s ratio, dimensionless

bruges.rockphysics.moduli.vp(youngs=None, vs=None, rho=None, mu=None, lam=None, bulk=None, pr=None, pmod=None)[source]#

Computes Vp given bulk density and any two elastic moduli (e.g. lambda and mu, or Young’s and P moduli). SI units only.

Parameters
  • lam (Any 2 from) –

  • mu

  • youngs

  • pr

  • pmod

  • bulk

  • Rho

Returns

Vp in m/s

bruges.rockphysics.moduli.vs(youngs=None, vp=None, rho=None, mu=None, lam=None, bulk=None, pr=None, pmod=None)[source]#

Computes Vs given bulk density and shear modulus. SI units only.

Parameters
  • Mu

  • Rho

Returns

Vs in m/s

bruges.rockphysics.moduli.youngs(vp=None, vs=None, rho=None, mu=None, lam=None, bulk=None, pr=None, pmod=None)[source]#

Computes Young’s modulus given either Vp, Vs, and rho, or any two elastic moduli (e.g. lambda and mu, or bulk and P moduli). SI units only.

Parameters
  • vp

  • vs

  • rho (and) –

  • lam (or any 2 from) –

  • mu

  • bulk

  • pr

  • pmod (and) –

Returns

Young’s modulus in pascals, Pa

bruges.rockphysics.rockphysicsmodels module#

bruges.rockphysics.rockphysicsmodels.constant_cement(K0, G0, phi, phi_cem=0.38, phi_c=0.4, Cn=8.6, Kc=37, Gc=45, scheme=2)[source]#

Constant cement model, Avseth et al. (2000). This model describes the elastic behaviour (K=bulk and G=shear moduli) of high porosity dry sand with a certain initial cementation by interpolating with the lower Hashin-Shtrikman bound the two end members at zero porosity and critical porosity. The zero porosity end member has K and G equal to mineral. The high porosity end member has K and G given by the Contact Cement model. It is assumed that the porosity reduction is due to non-cementing material filling in the available pore space.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • phi_cem (float, optional) – Porosity at initial cementation. Default: 0.38.

  • phi_c (float, optional) – Critical porosity. Default: 0.4.

  • Cn (float, optional) – Coordination number. Default: 8.6.

  • Kc (float, optional) – Cement bulk & shear modulus in GPa. Default: 37, 45.

  • Gc (float, optional) – Cement bulk & shear modulus in GPa. Default: 37, 45.

  • scheme (int, optional) – Cementation scheme, can be either 1 or 2: 1: cement deposited at grain contacts 2: cement as uniform layer around grains. Default: 2.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Dvorkin et al. (2014), Seismic Reflections of Rock Properties, Cambridge University Press (p.30-31)

bruges.rockphysics.rockphysicsmodels.contact_cement(K0, G0, phi, phi_c=0.4, Cn=8.6, Kc=37, Gc=45, scheme=2)[source]#

Contact cement or cemented sand model,. This model describes the elastic behaviour (K=bulk and G=shear moduli) of dry sand where cement is deposited at grain contacts. The cement properties can be modified as well as the type of cementation scheme.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • phi_c (float, optional) – Critical porosity. Default: 0.4

  • Cn (float, optional) – Coordination number Default: 8.6.

  • Kc (float, optional) – Cement bulk & shear modulus in GPa. Default: 37, 45.

  • Gc (float, optional) – Cement bulk & shear modulus in GPa. Default: 37, 45.

  • scheme (int, optional) – Cementation scheme, can be either 1 or 2: 1: cement deposited at grain contacts 2: cement as uniform layer around grains. Default: 2.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Dvorkin-Nur (1996), Elasticity of High-Porosity Sandstones: Theory for Two North Sea Data Sets. Geophysics 61, no. 5 (1996). Mavko et al. (2009), The Rock Physics Handbook, Cambridge University Press (p.255)

bruges.rockphysics.rockphysicsmodels.critical_porosity(K0, G0, phi, phi_c=0.4)[source]#

Critical porosity model. This model describes the elastic behaviour (K=bulk and G=shear moduli) of dry sand for porosities below the critical porosity (phi_c). Above phi_c, the fluid phase is load-bearing, below phi_c the solid phase (mineral grains) is load-bearing. The equations here describe the variation of K and G for porosities below phi_c as straight lines. Critical porosity is usually 0.4 for sandstone, 0.7 for chalk, 0.02-0.03 for granites.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • phi_c (float, optional) – Critical porosity. Default: 0.4

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Mavko et al. (2009), The Rock Physics Handbook, Cambridge University Press (p.370)

bruges.rockphysics.rockphysicsmodels.hertz_mindlin(K0, G0, sigma, phi_c=0.4, Cn=8.6, f=1)[source]#

Hertz-Mindlin model. This model describes the elastic behaviour (K=bulk and G=shear moduli) of a dry pack of spheres subject to a hydrostatic confining pressure.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • sigma (float) – Effective stress in MPa.

  • phi_c (float, optional) – Critical porosity. Default: 0.4

  • Cn (float, optional) – Coordination number Default: 8.6.

  • f (float, optional) – Shear modulus correction factor, f=1 for dry pack with perfect adhesion between particles and f=0 for dry frictionless pack.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Mavko et al. (2009), The Rock Physics Handbook, Cambridge University Press (p.246)

bruges.rockphysics.rockphysicsmodels.increasing_cement(K0, G0, phi, phi_cem=0.38, phi_c=0.4, Cn=8.6, Kc=37, Gc=45, scheme=2)[source]#

Increasing cement model (Modified Hashin-Shtrikman upper bound). This model describes the elastic behaviour (K=bulk and G=shear moduli) of a dry sand with a certain initial cementation by interpolating with the upper Hashin-Shtrikman bound the two end members at zero porosity and critical porosity. The zero porosity end member has K and G equal to mineral. The high porosity end member has K and G given by the Contact Cement model. Probably best to avoid using if for porosities>phi_cem. Need to check references.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • phi_cem (float, optional) – Porosity at initial cementation. Default: 0.38.

  • phi_c (float, optional) – Critical porosity. Default: 0.4.

  • Cn (float, optional) – Coordination number. Default: 8.6.

  • Kc (float, optional) – Cement bulk & shear modulus in GPa. Default: 37, 45.

  • Gc (float, optional) – Cement bulk & shear modulus in GPa. Default: 37, 45.

  • scheme (int, optional) – Cementation scheme, can be either 1 or 2: 1: cement deposited at grain contacts 2: cement as uniform layer around grains. Default: 2.

Returns

Kdry, Gdry – dry rock bulk & shear modulus in GPa.

Return type

float or array_like

bruges.rockphysics.rockphysicsmodels.soft_sand(K0, G0, phi, sigma, phi_c=0.4, Cn=8.6, f=1)[source]#

Soft sand, or friable sand or uncemented sand model. This model describes the elastic behaviour (K=bulk and G=shear moduli) of poorly sorted dry sand by interpolating with the lower Hashin-Shtrikman bound the two end members at zero porosity and critical porosity. The zero porosity end member has K and G equal to mineral. The end member at critical porosity has K and G given by Hertz-Mindlin model.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • sigma (float) – Effective stress in MPa.

  • phi_c (float, optional) – Critical porosity. Default: 0.4

  • Cn (float, optional) – Coordination number Default: 8.6.

  • f (float, optional) – Shear modulus correction factor, f=1 for dry pack with perfect adhesion between particles and f=0 for dry frictionless pack.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Mavko et al. (2009), The Rock Physics Handbook, Cambridge University Press (p.258)

bruges.rockphysics.rockphysicsmodels.stiff_sand(K0, G0, phi, sigma, phi_c=0.4, Cn=8.6, f=1)[source]#

Stiff sand model. This model describes the elastic behaviour (K=bulk and G=shear moduli) of stiff dry sands by interpolating with the upper Hashin-Shtrikman bound the two end members at zero porosity and critical porosity. The zero porosity end member has K and G equal to mineral. The end member at critical porosity has K and G given by Hertz-Mindlin model.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • sigma (float) – Effective stress in MPa.

  • phi_c (float, optional) – Critical porosity. Default: 0.4

  • Cn (float, optional) – Coordination number Default: 8.6.

  • f (float, optional) – Shear modulus correction factor, f=1 for dry pack with perfect adhesion between particles and f=0 for dry frictionless pack.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa

Return type

float or array_like

References

Mavko et al. (2009), The Rock Physics Handbook, Cambridge University Press (p.260)

bruges.rockphysics.rockphysicsmodels.vernik_consol_sand(K0, G0, phi, sigma, b=10)[source]#

Vernik & Kachanov Consolidated Sand Model. This model describes the elastic behaviour (K=bulk and G=shear moduli) of consolidated dry sand subject to a hydrostatic confining pressure as a continuous solid containing pores and cracks.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • sigma (float) – Effective stress in MPa.

  • b (float, optional) – Slope parameter in pore shape empirical equation, range: 8-12. Default: 10.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Vernik & Kachanov (2010), Modeling elastic properties of siliciclastic rocks, Geophysics v.75 n.6

bruges.rockphysics.rockphysicsmodels.vernik_sand_diagenesis(K0, G0, phi, sigma, phi_c=0.36, phi_con=0.26, b=10, n=2.0, m=2.05)[source]#

Vernik & Kachanov Sandstone Diagenesis Model. This model describes the elastic behaviour (K=bulk and G=shear moduli) of dry sand modeled as a continuous solid containing pores and cracks for porosities below phi_con (consolidation porosity) using Vernik’s Consolidated Sand Model, and as a granular material for porosities above phi_con using Vernik’s Soft Sand Model 1.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • sigma (float) – Effective stress in MPa.

  • phi_c (float, optional) – Critical porosity, range 0.30-0.42. Default: 0.36.

  • phi_con (float, optional) – Consolidation porosity, range 0.22-0.30. Default: 0.26.

  • b (float, optional) – Slope parameter in pore shape empirical equation, range: 8-12. Default: 10.

  • n (float) – Empirical factors. Default: 2.00, 2.05.

  • m (float) – Empirical factors. Default: 2.00, 2.05.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Vernik & Kachanov (2010), Modeling elastic properties of siliciclastic rocks, Geophysics v.75 n.6

bruges.rockphysics.rockphysicsmodels.vernik_shale(vclay, phi, rhom=2.73, rhob=1, Mqz=96, c33_clay=33.4, A=0.00284)[source]#

Vernik & Kachanov Shale Model. This model describes the elastic behaviour in terms of velocities and density of inorganic shales.

Parameters
  • vclay (float or array_like) – Dry clay content volume fraction.

  • phi (float or array_like) – Porosity, maximum 0.40.

  • rhom (float, optional) – Shale matrix density in g/cc. Default: 2.73.

  • rhob (float, optional) – Brine density in g/cc. Default: 1.

  • Mqz (float, optional) – P-wave elastic modulus of remaining minerals in GPa Default: 96.

  • c33_clay (float, optional) – Anisotropic clay constant in GPa. Default: 33.4.

  • A (float, optional) – Empirical coefficient for Vs. Default is good for illite/smectite/chlorite, can be raised up to .006 for kaolinite-rich clays. Default: 0.00284.

Returns

vp, vs, density – P- and S-wave velocities in m/s, density in g/cc.

Return type

float or array_like

Notes

Shale matrix density (rhom) averages 2.73 +/- 0.03 g/cc at porosities below 0.25. It gradually varies with compaction and smectite-to-illite transition. A more accurate estimate can be calculated with this equation: rhom = 2.76+0.001*((rho-2)-230*np.exp(-4*(rho-2)))

References

Vernik & Kachanov (2010), Modeling elastic properties of siliciclastic rocks, Geophysics v.75 n.6

bruges.rockphysics.rockphysicsmodels.vernik_soft_sand_1(K0, G0, phi, sigma, phi_c=0.36, phi_con=0.26, b=10, n=2.0, m=2.05)[source]#

Vernik & Kachanov Soft Sand Model 1. This model describes the elastic behaviour (K=bulk and G=shear moduli) of dry sand modeled as a granular material. Only applicable for porosities between the low-porosity end-member (at the consolidation porosity phi_con) and the high-porosity end-member (at the critical porosity phi_c). The low-porosity end member is calculated with Vernik’s Consolidated Sand Model.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • sigma (float) – Effective stress in MPa.

  • phi_c (float, optional) – Critical porosity, range 0.30-0.42. Default: 0.36.

  • phi_con (float, optional) – Consolidation porosity, range 0.22-0.30. Default: 0.26.

  • b (float, optional) – Slope parameter in pore shape empirical equation, range: 8-12. Default: 10.

  • n (float) – Empirical factors. Default: 2.00, 2.05.

  • m (float) – Empirical factors. Default: 2.00, 2.05.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Vernik & Kachanov (2010), Modeling elastic properties of siliciclastic rocks, Geophysics v.75 n.6

bruges.rockphysics.rockphysicsmodels.vernik_soft_sand_2(K0, G0, phi, p=20, q=20)[source]#

Vernik & Kachanov Soft Sand Model 2. This model describes the elastic behaviour (K=bulk and G=shear moduli) of dry sand modeled as a granular material. Applicable in the entire porosity range.

Parameters
  • K0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • G0 (float or array_like) – Mineral bulk & shear modulus in GPa.

  • phi (float or array_like) – Porosity.

  • p (float, optional) – Pore shape factor for K and G, range: 10-45. Default: 20.

  • q (float, optional) – Pore shape factor for K and G, range: 10-45. Default: 20.

Returns

Kdry, Gdry – Dry rock bulk & shear modulus in GPa.

Return type

float or array_like

References

Vernik & Kachanov (2010), Modeling elastic properties of siliciclastic rocks, Geophysics v.75 n.6

bruges.rockphysics.rpm module#

bruges.rockphysics.rpm.contactcement(K0, G0, phi, phic=0.4, Cn=8.6, Kc=37, Gc=45, scheme=2)[source]#

Contact cement (cemented sand) model, Dvorkin-Nur (1996) written by Alessandro Amato del Monte (2015) from Mavko et al., Rock Physics Handbook, p.255

Parameters
  • K0 – mineral bulk & shear modulus in GPa

  • G0 – mineral bulk & shear modulus in GPa

  • phi – porosity

  • phic – critical porosity (default 0.4)

  • Cn – coordination nnumber (default 8.6)

  • Kc – cement bulk & shear modulus in GPa (default 37, 45 i.e. quartz)

  • Gc – cement bulk & shear modulus in GPa (default 37, 45 i.e. quartz)

  • scheme – 1=cement deposited at grain contacts, 2=in uniform layer around grains (default 2)

Returns

K_DRY, G_DRY: dry rock bulk & shear modulus in GPa

Return type

Tuple

bruges.rockphysics.rpm.critpor(K0, G0, phi, phic=0.4)[source]#

Critical porosity, Nur et al. (1991, 1995) written by Alessandro Amato del Monte (2015) from Mavko et al., Rock Physics Handbook, p.353.

Parameters
  • K0 – float or array_like Mineral bulk & shear modulus in GPa

  • G0 – float or array_like Mineral bulk & shear modulus in GPa

  • phi – float or array_like porosity

  • phic – float or array_like critical porosity (default 0.4)

Returns

K_DRY, G_DRY: dry rock bulk & shear modulus in GPa

Return type

Tuple

bruges.rockphysics.rpm.hertzmindlin(K0, G0, phi, phic=0.4, Cn=8.6, P=10, f=1)[source]#

Hertz-Mindlin model written by Alessandro Amato del Monte (2015) from Mavko et al., Rock Physics Handbook, p.246

Parameters
  • K0 – mineral bulk & shear modulus in GPa

  • G0 – mineral bulk & shear modulus in GPa

  • phi – porosity

  • phic – critical porosity (default 0.4)

  • Cn – coordination nnumber (default 8.6)

  • P – confining pressure in MPa (default 10)

  • f – shear modulus correction factor, f=1 for dry pack with perfect adhesion between particles and f=0 for dry frictionless pack.

Returns

K_DRY, G_DRY: dry rock bulk & shear modulus in GPa

Return type

Tuple

bruges.rockphysics.rpm.softsand(K0, G0, phi, phic=0.4, Cn=8.6, P=10, f=1)[source]#

Soft-sand (uncemented) model written by Alessandro Amato del Monte (2015) from Mavko et al., Rock Physics Handbook, p.258

Parameters
  • K0 – mineral bulk & shear modulus in GPa

  • G0 – mineral bulk & shear modulus in GPa

  • phi – porosity

  • phic – critical porosity (default 0.4)

  • Cn – coordination nnumber (default 8.6)

  • P – confining pressure in MPa (default 10)

  • f – shear modulus correction factor, f=1 for dry pack with perfect adhesion between particles and f=0 for dry frictionless pack

Returns

K_DRY, G_DRY: dry rock bulk & shear modulus in GPa

Return type

Tuple

bruges.rockphysics.rpm.stiffsand(K0, G0, phi, phic=0.4, Cn=8.6, P=10, f=1)[source]#

Stiff-sand model written by Alessandro Amato del Monte (2015) from Mavko et al., Rock Physics Handbook, p.260

Parameters
  • K0 – mineral bulk & shear modulus in GPa

  • G0 – mineral bulk & shear modulus in GPa

  • phi – porosity

  • phic – critical porosity (default 0.4)

  • Cn – coordination nnumber (default 8.6)

  • P – confining pressure in MPa (default 10)

  • f – shear modulus correction factor, f=1 for dry pack with perfect adhesion between particles and f=0 for dry frictionless pack

Returns

Tuple: K_DRY, G_DRY: dry rock bulk & shear modulus in GPa

bruges.rockphysics.rpt module#

bruges.rockphysics.rpt.rpt(model='soft', vsh=0.0, fluid='gas', phic=0.4, Cn=8, P=10, f=1, cement='quartz')[source]#
bruges.rockphysics.rpt.vrh(volumes, k, mu)[source]#

Calculates Voigt-Reuss-Hill bounds.

Parameters
  • volumes – array with volumetric fractions

  • k – array with bulk modulus

  • mu – array with shear modulus

Returns

upper (Voigt) and lower (Reuss) average for k mu_u, mu_l: upper (Voigt) and lower (Reuss) average for mu k0, mu0: Hill average of k and mu

Return type

k_u, k_l

Module contents#