Gravimetry in 2D - Part IIΒΆ

Simple gravimetric field solution with Integration after [WB87].

../../_images/sphx_glr_plot_mod-gravimetry-integration-2d_001.png
import numpy as np
import pygimli as pg
from pygimli.meshtools import createCircle
from pygimli.physics.gravimetry import solveGravimetry
from pygimli.physics.gravimetry import gradUCylinderHoriz, gradGZCylinderHoriz
from pygimli.physics.gravimetry import gradUHalfPlateHoriz
from pygimli.physics.gravimetry import gradGZHalfPlateHoriz

radius = 2.
depth = 5.
rho = 1000.0

x = np.arange(-20, 30, 1)
pnts = np.zeros((len(x), 2))
pnts[:, 0] = x
pos = [0, -depth]


def plot(x, a1, ga, gza, a2, g, gz):
    a1.plot(x, ga[:, 0],  label=r'Analytical $\frac{\partial u}{\partial x}$')
    a1.plot(x, ga[:, 1],  label=r'Analytical $\frac{\partial u}{\partial z}$')

    a1.plot(x, g[:, 0], label=r'Won & Bevis: $\frac{\partial u}{\partial x}$',
            marker='o', linewidth=0)
    a1.plot(x, g[:, 2], label=r'Won & Bevis: $\frac{\partial u}{\partial z}$',
            marker='o', linewidth=0)

    a2.plot(x, gza[:, 0],
            label=r'Analytical $\frac{\partial^2 u}{\partial z,x}$')
    a2.plot(x, gza[:, 1],
            label=r'Analytical $\frac{\partial^2 u}{\partial z,z}$')

    a2.plot(x, gz[:, 0], marker='o', linestyle='',
            label=r'Won & Bevis: $\frac{\partial^2 u}{\partial z,x}$')
    a2.plot(x, gz[:, 2], marker='o', linestyle='',
            label=r'Won & Bevis: $\frac{\partial^2 u}{\partial z,z}$')
    a1.set_xlabel('$x$-coordinate [m]')
    a1.set_ylabel(r'$\frac{\partial u}{\partial (x,z)}$ [mGal]')
    a1.legend(loc='best')

    a2.set_xlabel('$x$-coordinate [m]')
    a2.legend(loc='best')


fig = pg.plt.figure(figsize=(8,8))
ax = [fig.add_subplot(2, 2, i) for i in range(1, 5)]

# Horizontal cylinder

ga = gradUCylinderHoriz(pnts, radius, rho, pos=pos)
gza = gradGZCylinderHoriz(pnts, radius, rho, pos=pos)

circ = createCircle([0, -depth], radius=radius, marker=2, area=0.1,
                    segments=32)
g, gz = solveGravimetry(circ, rho, pnts, complete=True)

plot(x, ax[0], ga, gza, ax[1], g, gz)

# Half plate

thickness = 0.1

# mesh = pg.createGrid(x=[-2,2], y=[-2,2], z=[-3,-7])
mesh = pg.createGrid(x=np.linspace(0, 5000, 2),
                     y=[-depth-thickness/2.0, -depth+thickness/2.0])

ga = gradUHalfPlateHoriz(pnts, thickness, rho, pos=[0, -depth])
gza = gradGZHalfPlateHoriz(pnts, thickness, rho, pos=[0, -depth])
g, gz = solveGravimetry(mesh, rho, pnts, complete=True)

plot(x, ax[2], ga, gza, ax[3], g, gz)

Total running time of the script: ( 0 minutes 0.377 seconds)

Generated by Sphinx-Gallery

News (Aug 2017):
pyGIMLi Paper out now! Download PDF here.


© 2017 - GIMLi Development Team
Created using Bootstrap and Sphinx. Last updated on Aug 15, 2017.