There will be a webinar on pyGIMLi hosted by SEG on March 19, 2024 at 4 pm CET. Register for free here.

pygimli.physics.gravimetry#

Solve gravimetric and magneto static problems in 2D and 3D analytically

Overview#

Functions

BZPoly(pnts, poly, mag[, openPoly])

Vertical magnetic gradient for polygone.

BaZCylinderHoriz(pnts, R, pos, M)

Magnetic anomaly for a horizontal cylinder.

BaZSphere(pnts, R, pos, M)

Magnetic anomaly for a sphere.

SolveGravMagHolstein(*args, **kwargs)

gradGZCylinderHoriz(r, a, rho[, pos])

TODO WRITEME.

gradGZHalfPlateHoriz(pnts, t, rho[, pos])

TODO WRITEME.

gradGZSphere(r, rad, rho[, pos])

TODO WRITEME.

gradUCylinderHoriz(r, a, rho[, pos])

2D Gradient of gravimetric potential of horizontal cylinder.

gradUHalfPlateHoriz(pnts, t, rho[, pos])

Gravitational field od a horizontal half plate.

gradUSphere(r, rad, rho[, pos])

Gravitational field of a sphere.

solveGravimetry(mesh[, dDensity, pnts, complete])

Solve gravimetric response.

uCylinderHoriz(pnts, rad, rho[, pos])

Gravitational potential of horizonzal cylinder.

uSphere(r, rad, rho[, pos])

Gravitational potential of a sphere.

Classes

GravityModelling(mesh, points[, cmp, foot])

Magnetics modelling operator using Holstein (2007).

GravityModelling2D([points])

Gravimetry modelling operator.

MagManager([data])

Magnetics Manager.

MagneticsModelling([mesh, points, cmp, ...])

Magnetics modelling operator using Holstein (2007).

Functions#

pygimli.physics.gravimetry.BZPoly(pnts, poly, mag, openPoly=False)[source]#

Vertical magnetic gradient for polygone.

Parameters:
  • pnts (list) – Measurement points [[p1x, p1z], [p2x, p2z],…]

  • poly (list) – Polygon [[p1x, p1z], [p2x, p2z],…]

  • mag ([M_x, M_y, M_z]) – Magnetization = [M_x, M_y, M_z]

pygimli.physics.gravimetry.BaZCylinderHoriz(pnts, R, pos, M)[source]#

Magnetic anomaly for a horizontal cylinder.

Calculate the vertical component of the anomalous magnetic field Bz for a buried horizontal cylinder at position pos with radius R for a given magnetization M at measurement points pnts.

TODO .. only 2D atm

Parameters:
  • pnts ([[x,z], ]) – measurement points – array[x,y,z]

  • R (float) – radius

  • pos ([float, float]) – [x,z] – sphere center

  • M ([float, float]) – [Mx, Mz] – magnetization

pygimli.physics.gravimetry.BaZSphere(pnts, R, pos, M)[source]#

Magnetic anomaly for a sphere.

Calculate the vertical component of the anomalous magnetic field Bz for a buried sphere at position pos with radius R for a given magnetization M at measurement points pnts.

Parameters:
  • pnts ([[x,y,z], ]) – measurement points – array[x,y,z]

  • R (float) – radius

  • pos ([float, float, float]) – [x,y,z] – sphere center

  • M ([float, float, float]) – [Mx, My, Mz] – magnetization

pygimli.physics.gravimetry.SolveGravMagHolstein(*args, **kwargs)#
pygimli.physics.gravimetry.gradGZCylinderHoriz(r, a, rho, pos=(0.0, 0.0))[source]#

TODO WRITEME.

\[g = -grad u(r), with r = [x,z], |r| = \sqrt{x^2+z^2}\]
Parameters:
  • r (list[[x, z]]) – Observation positions

  • a (float) – Cylinder radius in [meter]

  • rho – Density in [kg/m^3]

Return type:

grad gz, [gz_x, gz_z]

Examples using pygimli.physics.gravimetry.gradGZCylinderHoriz

Semianalytical Gravimetry and Geomagnetics in 2D

Semianalytical Gravimetry and Geomagnetics in 2D
pygimli.physics.gravimetry.gradGZHalfPlateHoriz(pnts, t, rho, pos=(0.0, 0.0))[source]#

TODO WRITEME.

\[g = -\nabla u\]
Parameters:
  • pnts (array (\(n\times 2\))) – n 2 dimensional measurement points

  • t (float) – Plate thickness in \([\text{m}]\)

  • rho (float) – Density in \([\text{kg}/\text{m}^3]\)

Returns:

gz – Gradient of z-component of g \(\nabla(\frac{\partial u}{\partial\vec{r}}_z)\)

Return type:

array

Examples using pygimli.physics.gravimetry.gradGZHalfPlateHoriz

Semianalytical Gravimetry and Geomagnetics in 2D

Semianalytical Gravimetry and Geomagnetics in 2D
pygimli.physics.gravimetry.gradGZSphere(r, rad, rho, pos=(0.0, 0.0, 0.0))[source]#

TODO WRITEME.

\[g = -\nabla u\]
Parameters:
Return type:

[d g_z /dx, d g_z /dy, d g_z /dz]

pygimli.physics.gravimetry.gradUCylinderHoriz(r, a, rho, pos=(0.0, 0.0))[source]#

2D Gradient of gravimetric potential of horizontal cylinder.

\[g = -G[m^3/(kg s^2)] * dM[kg/m] * 1/r[1/m] * grad(r)[1/1] = [m^3/(kg s^2)] * [kg/m] * 1/m * [1/1] == m/s^2\]
Parameters:
  • r (list[[x, z]]) – Observation positions

  • a (float) – Cylinder radius in [meter]

  • pos ([x,z]) – Center position of cylinder.

  • rho (float) – Delta density in [kg/m^3]

Returns:

g – Gradient of gravimetry potential [mGal].

Return type:

[dudx, dudz]

Examples using pygimli.physics.gravimetry.gradUCylinderHoriz

Semianalytical Gravimetry and Geomagnetics in 2D

Semianalytical Gravimetry and Geomagnetics in 2D

Gravimetry in 2D - Part I

Gravimetry in 2D - Part I
pygimli.physics.gravimetry.gradUHalfPlateHoriz(pnts, t, rho, pos=(0.0, 0.0))[source]#

Gravitational field od a horizontal half plate.

\[g = -grad u,\]
Parameters:
  • pnts

  • t

  • rho – Density in [kg/m^3]

Returns:

z-component of g .. math:: nabla(partial u/partialvec{r})_z

Return type:

gz

Examples using pygimli.physics.gravimetry.gradUHalfPlateHoriz

Semianalytical Gravimetry and Geomagnetics in 2D

Semianalytical Gravimetry and Geomagnetics in 2D
pygimli.physics.gravimetry.gradUSphere(r, rad, rho, pos=(0.0, 0.0, 0.0))[source]#

Gravitational field of a sphere.

\[g = -G[m^3/(kg s^2)] * dM[kg] * 1/r^2 1/m^2] * \grad(r)[1/1] = [m^3/(kg s^2)] * [kg] * 1/m^2 * [1/1] == m/s^2\]
Parameters:
Returns:

[gx, gy, gz] – gravitational acceleration (note that gz points negative)

Return type:

[float*3]

pygimli.physics.gravimetry.solveGravimetry(mesh, dDensity=None, pnts=None, complete=False)[source]#

Solve gravimetric response.

2D with pygimli.physics.gravimetry.lineIntegralZ_WonBevis

3D with pygimli.physics.gravimetry.gravMagBoundarySinghGup

Parameters:
  • mesh (GIMLI::Mesh) – 2d or 3d mesh with or without cells.

  • dDensity (float | array) –

    Density difference.

    • float – solve for positive boundary marker only.

      Assuming one inhomogeneity.

    • [[int, float]] – solve for multiple positive boundaries TOIMPL

    • array – solve for one delta density value per cell

    • None – return per cell kernel matrix G TOIMPL

  • pnts ([[x_i, y_i]]) – List of measurement positions.

  • complete (bool [False]) – If True return whole solution or matrix for [dgx, dgy, dgz]

Returns:

  • dg (array OR)

  • dz, dgz (arrays (if complete))

Examples using pygimli.physics.gravimetry.solveGravimetry

Semianalytical Gravimetry and Geomagnetics in 2D

Semianalytical Gravimetry and Geomagnetics in 2D

Gravimetry in 2D - Part I

Gravimetry in 2D - Part I
pygimli.physics.gravimetry.uCylinderHoriz(pnts, rad, rho, pos=[0.0, 0.0])[source]#

Gravitational potential of horizonzal cylinder.

Parameters:
  • pnts (iterable) – measuring point locations

  • rad (float) – radius of the cylinder

  • rho (float) – density contrast in kg/m^3

  • pos ([float, float]) – x,z position of the cylinder axis

Return type:

gravimetric potential at the given points

pygimli.physics.gravimetry.uSphere(r, rad, rho, pos=None)[source]#

Gravitational potential of a sphere.

Gravitational potential of a sphere with radius and density at a given position.

\[u = -G * dM * \frac{1}{r}\]
Parameters:

Classes#

class pygimli.physics.gravimetry.GravityModelling(mesh, points, cmp=['gz'], foot=None)[source]#

Bases: MeshModelling

Magnetics modelling operator using Holstein (2007).

__init__(mesh, points, cmp=['gz'], foot=None)[source]#

Setup forward operator.

Parameters:
  • mesh (pygimli:mesh) – tetrahedral or hexahedral mesh

  • points (list|array of (x, y, z)) – measuring points

  • cmp (list of str) – component of: gx, gy, gz, TFA, Bx, By, Bz, Bxy, Bxz, Byy, Byz, Bzz

createJacobian(model)[source]#

Do nothing as this is a linear problem.

createKernel()[source]#

Create computational kernel.

response(model)[source]#

Compute forward response.

setMesh(mesh, ignoreRegionManager=False)[source]#

Set the mesh.

class pygimli.physics.gravimetry.GravityModelling2D(points=None, **kwargs)[source]#

Bases: MeshModelling

Gravimetry modelling operator.

__init__(points=None, **kwargs)[source]#

Initialize forward operator, optional with mesh and points.

You can specify both the mesh and the measuring points, or set the latter after the mesh has been set.

Parameters:
  • mesh (pg.Mesh) – mesh for forward computation

  • points (array[x,y]) – measuring points

calcMatrix()[source]#

Create Jacobian matrix (density-independent).

createJacobian(model)[source]#

Create Jacobian.

createStartmodel(*args)[source]#

Create the default starting model.

response(model)[source]#

Calculate response for a given density distribution.

setSensorPositions(pnts)[source]#

Set measurement locations. [[x,y,z],…].

class pygimli.physics.gravimetry.MagManager(data=None, **kwargs)[source]#

Bases: MeshMethodManager

Magnetics Manager.

__init__(data=None, **kwargs)[source]#

Create Magnetics Manager instance.

createForwardOperator(**kwargs)[source]#

Create forward operator (computationally extensive!).

createGrid(dx=50, depth=800, bnd=0)[source]#

Create a grid.

createMesh(bnd=0, area=100000.0, depth=800, quality=1.3, addPLC=None, addPoints=True)[source]#

Create an unstructured mesh.

inversion(noise_level=2, noisify=False, **kwargs)[source]#

Run Inversion (requires mesh and FOP).

Parameters:
  • noise_level (float|array) – absolute noise level (absoluteError)

  • noisify (bool) – add noise before inversion

  • relativeError (float|array [0.01]) – relative error to stabilize very low data

  • depthWeighting (bool [True]) – apply depth weighting after Li&Oldenburg (1996)

  • z0 (float) – skin depth for depth weighting

  • mul (array) – multiply constraint weight with

  • arguments (standard inversion keyword) –

  • ....................................

  • C( (int|Matrix|[float, float, float]) – constraint order, matrix or correlation lengths

  • cType) (int|Matrix|[float, float, float]) – constraint order, matrix or correlation lengths

  • limits ([float, float]) – lower and upper parameter limits

  • symlogThreshold (float [0]) – threshold for symlog data trans (0 = linear)

  • startModel (float|array) – starting model (typically homogeneous)

  • lam (float) – regularization strength

  • robustData (bool) – L1 norm on data misfit and model roughness

  • blockyModel (bool) – L1 norm on data misfit and model roughness

  • maxIter (int) – maximum iteration number

Returns:

model – model vector (also saved in self.inv.model)

Return type:

array

show3DModel(label=None, trsh=0.025, synth=None, invert=False, position='yz', elevation=10, azimuth=25, zoom=1.2, **kwargs)[source]#

Standard 3D view.

showData(cmp=None, **kwargs)[source]#

Show data.

showDataFit()[source]#

Show data, model response and misfit.

class pygimli.physics.gravimetry.MagneticsModelling(mesh=None, points=None, cmp=['TFA'], igrf=[50, 13], foot=None)[source]#

Bases: MeshModelling

Magnetics modelling operator using Holstein (2007).

__init__(mesh=None, points=None, cmp=['TFA'], igrf=[50, 13], foot=None)[source]#

Setup forward operator.

Parameters:
  • mesh (pygimli:mesh) – tetrahedral or hexahedral mesh

  • points (list|array of (x, y, z)) – measuring points

  • cmp (list of str) – component of: gx, gy, gz, TFA, Bx, By, Bz, Bxy, Bxz, Byy, Byz, Bzz

  • igrf (list|array of size 3 or 7) –

    international geomagnetic reference field, either [D, I, H, X, Y, Z, F] - declination, inclination, horizontal field,

    X/Y/Z components, total field OR

    [X, Y, Z] - X/Y/Z components [lat, lon] - latitude, longitude (automatic IGRF)

computeKernel()[source]#

Compute the kernel.

createJacobian(model)[source]#

Do nothing as this is a linear problem.

response(model)[source]#

Compute forward response.