esys.downunder.forwardmodels Package

Classes

class esys.downunder.forwardmodels.AcousticWaveForm(domain, omega, w, data, F, coordinates=None, fixAtBottom=False, tol=1e-10, saveMemory=True, scaleF=True)

Forward Model for acoustic waveform inversion in the frequency domain. It defines a cost function:

Math:

defect = 1/2 integrate( ( w * ( a * u - data ) ) ** 2 )

where w are weighting factors, data are the measured data (as a 2-comp vector of real and imaginary part) for real frequency omega, and u is the corresponding result produced by the forward model. u (as a 2-comp vector) is the solution of the complex Helmholtz equation for frequency omega, source F and complex, inverse, squared p-velocity sigma:

Math:

-u_{ii} - omega**2 * sigma * u = F

It is assumed that the exact scale of source F is unknown and the scaling factor a of F is calculated by minimizing the defect.

__init__(domain, omega, w, data, F, coordinates=None, fixAtBottom=False, tol=1e-10, saveMemory=True, scaleF=True)

initializes a new forward model with acoustic wave form inversion.

Parameters:

  • domain (Domain) – domain of the model

  • w (Scalar) – weighting factors

  • data (escript.Data of shape (2,)) – real and imaginary part of data

  • F (escript.Data of shape (2,)) – real and imaginary part of source given at Dirac points, on surface or at volume.

  • coordinates (ReferenceSystem or SpatialCoordinateTransformation) – defines coordinate system to be used (not supported yet)

  • tol (positive float) – tolerance of underlying PDE

  • saveMemory (bool) – if true stiffness matrix is deleted after solution of PDE to minimize memory requests. This will require more compute time as the matrix needs to be reallocated.

  • scaleF (bool) – if true source F is scaled to minimize defect.

  • fixAtBottom (bool) – if true pressure is fixed to zero at the bottom of the domain

getArguments(sigma)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

sigma (escript.Data of shape (2,)) – a suggestion for complex 1/V**2

Returns:

solution, uTar, uTai, uTu

Return type:

escript.Data of shape (2,), 3 x float

getCoordinateTransformation()

returns the coordinate transformation being used

Return type:

CoordinateTransformation

getDefect(sigma, u, uTar, uTai, uTu)

Returns the defect value.

Parameters:

  • sigma (escript.Data of shape (2,)) – a suggestion for complex 1/V**2

  • u (escript.Data of shape (2,)) – a u vector

  • uTar (float) – equals integrate( w  * (data[0]*u[0]+data[1]*u[1]))

  • uTai – equals integrate( w  * (data[1]*u[0]-data[0]*u[1]))

  • uTu (float) – equals integrate( w  * (u,u))

Return type:

float

getDomain()

Returns the domain of the forward model.

Return type:

Domain

getGradient(sigma, u, uTar, uTai, uTu)

Returns the gradient of the defect with respect to density.

Parameters:

  • sigma (escript.Data of shape (2,)) – a suggestion for complex 1/V**2

  • u (escript.Data of shape (2,)) – a u vector

  • uTar (float) – equals integrate( w  * (data[0]*u[0]+data[1]*u[1]))

  • uTai – equals integrate( w  * (data[1]*u[0]-data[0]*u[1]))

  • uTu (float) – equals integrate( w  * (u,u))

getSourceScaling(u)

returns the scaling factor s required to rescale source F to minimize defect |s * u- data|^2

Parameters:

u (escript.Data of shape (2,)) – value of pressure solution (real and imaginary part)

Return type:

complex

getSurvey(index=None)

Returns the pair (data, weight)

If argument index is ignored.

rescaleWeights(scale=1.0, sigma_scale=1.0)

rescales the weights such that

Math:

integrate( ( w omega**2 * sigma_scale * data * ((1/L_j)**2)**-1) +1 )/(data*omega**2 * ((1/L_j)**2)**-1) * sigma_scale )=scale

Parameters:

  • scale (positive float) – scale of data weighting factors

  • sigma_scale (Scalar) – scale of 1/vp**2 velocity.

setUpPDE()

Creates and returns the underlying PDE.

Return type:

lpde.LinearPDE

class esys.downunder.forwardmodels.DcRes(domain, locator, delphiIn, sampleTags, phiPrimary, sigmaPrimary, w=1.0, coordinates=None, tol=1e-08, saveMemory=True, b=None)

Forward Model for DC resistivity, with a given source pair. The cost function is defined as:

Math:

defect = 1/2 (sum_s sum_pq w_pqs * ((phi_sp-phi_sq)-v_pqs)**2

__init__(domain, locator, delphiIn, sampleTags, phiPrimary, sigmaPrimary, w=1.0, coordinates=None, tol=1e-08, saveMemory=True, b=None)

setup new forward model

Parameters:

  • domain – the domain of the model

  • locator – contains locator to the measurement pairs

  • sampleTags (list of tuples) – tags of measurement points from which potential differences will be calculated.

  • phiPrimary (Scalar) – primary potential.

Type:

escript domain

Type:

list of Locator

Param:

delphiIn: this is v_pq, the potential difference for the current source and a set of measurement pairs. A list of measured potential differences is expected. Note this should be the secondary potential only.

getArguments(sigma)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

sigma (Data of shape (1,)) – conductivity

Returns:

phi

Return type:

Data of shape (1,)

getCoordinateTransformation()

returns the coordinate transformation being used

Return type:

CoordinateTransformation

getDefect(sigma, phi, loc_phi)

Returns the defect value.

Parameters:

  • sigma (Data of shape (1,)) – a suggestion for conductivity

  • phi (Data of shape (1,)) – potential field

Return type:

float

getDomain()

Returns the domain of the forward model.

Return type:

Domain

getGradient(sigma, phi, loc_phi)

Returns the gradient of the defect with respect to density.

Parameters:

  • sigma (Data of shape (1,)) – a suggestison for conductivity

  • phi (Data of shape (1,)) – potential field

getPrimaryPotential()

returns the primary potential :rtype: Data

setUpPDE()

Return the underlying PDE.

Return type:

LinearPDE

class esys.downunder.forwardmodels.ForwardModel

An abstract forward model that can be plugged into a cost function. Subclasses need to implement getDefect(), getGradient(), and possibly getArguments() and ‘getCoordinateTransformation’.

__init__()
getArguments(x)
getCoordinateTransformation()
getDefect(x, *args)
getGradient(x, *args)
class esys.downunder.forwardmodels.ForwardModelWithPotential(domain, w, data, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Base class for a forward model using a potential such as magnetic or gravity. It defines a cost function:

defect = 1/2 sum_s integrate( ( weight_i[s] * ( r_i - data_i[s] ) )**2 )

where s runs over the survey, weight_i are weighting factors, data_i are the data, and r_i are the results produced by the forward model. It is assumed that the forward model is produced through postprocessing of the solution of a potential PDE.

__init__(domain, w, data, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

initializes a new forward model with potential.

Parameters:

  • domain (Domain) – domain of the model

  • w (Vector or list of Vector) – data weighting factors

  • data (Vector or list of Vector) – data

  • coordinates (ReferenceSystem or SpatialCoordinateTransformation) – defines coordinate system to be used

  • fixPotentialAtBottom (bool) – if true potential is fixed to zero at the bottom of the domain in addition to the top.

  • tol (positive float) – tolerance of underlying PDE

getCoordinateTransformation()

returns the coordinate transformation being used

Return type:

CoordinateTransformation

getData()

Returns the data

Return type:

list of Data

getDataFunctionSpace()

Returns the FunctionSpace of the data

Return type:

FunctionSpace

getDefectGradient(result)
getDomain()

Returns the domain of the forward model.

Return type:

Domain

getMisfitWeights()

Returns the weights of the misfit function

Return type:

list of Data

getPDE()

Return the underlying PDE.

Return type:

LinearPDE

getSurvey(index=None)

Returns the pair (data_index, weight_index), where data_i is the data of survey i, weight_i is the weighting factor for survey i. If index is None, all surveys will be returned in a pair of lists.

class esys.downunder.forwardmodels.GravityModel(domain, w, g, gravity_constant=6.6742e-11, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Forward Model for gravity inversion as described in the inversion cookbook.

__init__(domain, w, g, gravity_constant=6.6742e-11, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Creates a new gravity model on the given domain with one or more surveys (w, g).

Parameters:

  • domain (Domain) – domain of the model

  • w (Vector or list of Vector) – data weighting factors

  • g (Vector or list of Vector) – gravity anomaly data

  • coordinates (ReferenceSystem` or SpatialCoordinateTransformation) – defines coordinate system to be used

  • tol (positive float) – tolerance of underlying PDE

  • fixPotentialAtBottom (bool) – if true potential is fixed to zero at the base of the domain in addition to the top

Note:

It is advisable to call rescaleWeights() to rescale weights before starting the inversion.

getArguments(rho)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

rho (Scalar) – a suggestion for the density distribution

Returns:

gravity potential and corresponding gravity field.

Return type:

Scalar, Vector

getDefect(rho, phi, gravity_force)

Returns the value of the defect

Parameters:

  • rho (Scalar) – density distribution

  • phi (Scalar) – corresponding potential

  • gravity_force (Vector) – gravity force

Return type:

float

getGradient(rho, phi, gravity_force)

Returns the gradient of the defect with respect to density.

Parameters:

  • rho (Scalar) – density distribution

  • phi (Scalar) – corresponding potential

  • gravity_force (Vector) – gravity force

Return type:

Scalar

getPotential(rho)

Calculates the gravity potential for a given density distribution.

Parameters:

rho (Scalar) – a suggestion for the density distribution

Returns:

gravity potential

Return type:

Scalar

rescaleWeights(scale=1.0, rho_scale=1.0)

rescales the weights such that

sum_s integrate( ( w_i[s] *g_i[s]) (w_j[s]*1/L_j) * L**2 * 4*pi*G*rho_scale )=scale

Parameters:

  • scale (positive float) – scale of data weighting factors

  • rho_scale (Scalar) – scale of density.

class esys.downunder.forwardmodels.IsostaticPressure(domain, p0=0.0, level0=0, gravity0=-9.81, background_density=2670.0, gravity_constant=6.6742e-11, coordinates=None, tol=1e-08)

class to calculate isostatic pressure field correction due to gravity forces

__init__(domain, p0=0.0, level0=0, gravity0=-9.81, background_density=2670.0, gravity_constant=6.6742e-11, coordinates=None, tol=1e-08)
Parameters:

  • domain (Domain) – domain of the model

  • p0 (scalar Data or float) – pressure at level0

  • background_density (float) – defines background_density in kg/m^3

  • coordinates (ReferenceSystem` or SpatialCoordinateTransformation) – defines coordinate system to be used

  • tol (positive float) – tolerance of underlying PDE

  • level0 (float) – pressure for z>=`level0` is set to zero.

  • gravity0 (float) – vertical background gravity at level0

getPressure(g=None, rho=None)

return the pressure for gravity force anomaly g and density anomaly rho

Parameters:

  • g (Vector) – gravity anomaly data

  • rho (Scalar) – gravity anomaly data

Returns:

pressure distribution

Return type:

Scalar

class esys.downunder.forwardmodels.MT2DModelTEMode(domain, omega, x, Z_XY, eta=None, w0=1.0, mu=1.2566370614359173e-06, sigma0=0.01, airLayerLevel=None, coordinates=None, Ex_top=1, fixAtTop=False, tol=1e-08, saveMemory=False, directSolver=True)

Forward Model for two dimensional MT model in the TE mode for a given frequency omega. It defines a cost function:

  • defect = 1/2 integrate( sum_s w^s * ( E_x/H_y - Z_XY^s ) ) ** 2 *

where E_x is the horizontal electric field perpendicular to the YZ-domain, horizontal magnetic field H_y=1/(i*omega*mu) * E_{x,z} with complex unit i and permeability mu. The weighting factor w^s is set to

  • w^s(X) = w_0^s *

if length(X-X^s) <= eta and zero otherwise. X^s is the location of impedance measurement Z_XY^s, w_0^s is the level of confidence (eg. 1/measurement error) and eta is level of spatial confidence.

E_x is given as solution of the PDE

  • -E_{x,ii} - i omega * mu * sigma * E_x = 0

where E_x at top and bottom is set to solution for background field. Homogeneous Neuman conditions are assumed elsewhere.

__init__(domain, omega, x, Z_XY, eta=None, w0=1.0, mu=1.2566370614359173e-06, sigma0=0.01, airLayerLevel=None, coordinates=None, Ex_top=1, fixAtTop=False, tol=1e-08, saveMemory=False, directSolver=True)

initializes a new forward model. See base class for a description of the arguments.

getArguments(sigma)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

sigma (Data of shape (2,)) – conductivity

Returns:

E_x, E_z

Return type:

Data of shape (2,)

getDefect(sigma, Ex, dExdz)

Returns the defect value.

Parameters:

  • sigma (Data of shape ()) – a suggestion for conductivity

  • Ex (Data of shape (2,)) – electric field

  • dExdz (Data of shape (2,)) – vertical derivative of electric field

Return type:

float

getGradient(sigma, Ex, dExdz)

Returns the gradient of the defect with respect to density.

Parameters:

  • sigma (Data of shape ()) – a suggestion for conductivity

  • Ex (Data of shape (2,)) – electric field

  • dExdz (Data of shape (2,)) – vertical derivative of electric field

class esys.downunder.forwardmodels.MT2DModelTMMode(domain, omega, x, Z_YX, eta=None, w0=1.0, mu=1.2566370614359173e-06, sigma0=0.01, airLayerLevel=None, coordinates=None, tol=1e-08, saveMemory=False, directSolver=True)

Forward Model for two-dimensional MT model in the TM mode for a given frequency omega. It defines a cost function:

  • defect = 1/2 integrate( sum_s w^s * ( rho*H_x/Hy - Z_YX^s ) ) ** 2 *

where H_x is the horizontal magnetic field perpendicular to the YZ-domain, horizontal magnetic field H_y=1/(i*omega*mu) * E_{x,z} with complex unit i and permeability mu. The weighting factor w^s is set to

  • w^s(X) = w_0^s *

if length(X-X^s) <= eta and zero otherwise. X^s is the location of impedance measurement Z_XY^s, w_0^s is the level of confidence (eg. 1/measurement error) and eta is level of spatial confidence.

H_x is given as solution of the PDE

  • -(rho*H_{x,i})_{,i} + i omega * mu * H_x = 0

where H_x at top and bottom is set to solution for background field. Homogeneous Neuman conditions are assumed elsewhere.

__init__(domain, omega, x, Z_YX, eta=None, w0=1.0, mu=1.2566370614359173e-06, sigma0=0.01, airLayerLevel=None, coordinates=None, tol=1e-08, saveMemory=False, directSolver=True)

initializes a new forward model. See base class for a description of the arguments.

getArguments(rho)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

rho (Data of shape (2,)) – resistivity

Returns:

Hx, grad(Hx)

Return type:

tuple of Data

getDefect(rho, Hx, g_Hx)

Returns the defect value.

Parameters:

  • rho (Data of shape ()) – a suggestion for resistivity

  • Hx (Data of shape (2,)) – magnetic field

  • g_Hx (Data of shape (2,2)) – gradient of magnetic field

Return type:

float

getGradient(rho, Hx, g_Hx)

Returns the gradient of the defect with respect to resistivity.

Parameters:

  • rho (Data of shape ()) – a suggestion for resistivity

  • Hx (Data of shape (2,)) – magnetic field

  • g_Hx (Data of shape (2,2)) – gradient of magnetic field

class esys.downunder.forwardmodels.MagneticIntensityModel(domain, w, b, background_magnetic_flux_density, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Forward Model for magnetic intensity inversion as described in the inversion cookbook.

__init__(domain, w, b, background_magnetic_flux_density, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Creates a new magnetic intensity model on the given domain with one or more surveys (w, b).

Parameters:

  • domain (Domain) – domain of the model

  • w (Scalar or list of Scalar) – data weighting factors

  • b (Scalar or list of Scalar) – magnetic intensity field data

  • tol (positive float) – tolerance of underlying PDE

  • background_magnetic_flux_density (Vector or list of float) – background magnetic flux density (in Tesla) with components (B_east, B_north, B_vertical)

  • coordinates (None) – defines coordinate system to be used

  • fixPotentialAtBottom (bool) – if true potential is fixed to zero at the bottom of the domain in addition to the top

getArguments(k)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

k (Scalar) – susceptibility

Returns:

scalar magnetic potential and corresponding magnetic field

Return type:

Scalar, Vector

getDefect(k, phi, magnetic_flux_density)

Returns the value of the defect.

Parameters:

  • k (Scalar) – susceptibility

  • phi (Scalar) – corresponding potential

  • magnetic_flux_density (Vector) – magnetic field

Return type:

float

getGradient(k, phi, magnetic_flux_density)

Returns the gradient of the defect with respect to susceptibility.

Parameters:

  • k (Scalar) – susceptibility

  • phi (Scalar) – corresponding potential

  • magnetic_flux_density (Vector) – magnetic field

Return type:

Scalar

getPotential(k)

Calculates the magnetic potential for a given susceptibility.

Parameters:

k (Scalar) – susceptibility

Returns:

magnetic potential

Return type:

Scalar

rescaleWeights(scale=1.0, k_scale=1.0)

rescales the weights such that

sum_s integrate( ( w_i[s] *B_i[s]) (w_j[s]*1/L_j) * L**2 * (background_magnetic_flux_density_j[s]*1/L_j) * k_scale )=scale

Parameters:

  • scale (positive float) – scale of data weighting factors

  • k_scale (Scalar) – scale of susceptibility.

class esys.downunder.forwardmodels.MagneticModel(domain, w, B, background_magnetic_flux_density, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Forward Model for magnetic inversion as described in the inversion cookbook.

__init__(domain, w, B, background_magnetic_flux_density, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Creates a new magnetic model on the given domain with one or more surveys (w, B).

Parameters:

  • domain (Domain) – domain of the model

  • w (Vector or list of Vector) – data weighting factors

  • B (Vector or list of Vector) – magnetic field data

  • tol (positive float) – tolerance of underlying PDE

  • background_magnetic_flux_density (Vector or list of float) – background magnetic flux density (in Tesla) with components (B_east, B_north, B_vertical)

  • coordinates (ReferenceSystem or SpatialCoordinateTransformation) – defines coordinate system to be used

  • fixPotentialAtBottom (bool) – if true potential is fixed to zero at the bottom of the domain in addition to the top

getArguments(k)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

k (Scalar) – susceptibility

Returns:

scalar magnetic potential and corresponding magnetic field

Return type:

Scalar, Vector

getDefect(k, phi, magnetic_flux_density)

Returns the value of the defect.

Parameters:

  • k (Scalar) – susceptibility

  • phi (Scalar) – corresponding potential

  • magnetic_flux_density (Vector) – magnetic field

Return type:

float

getGradient(k, phi, magnetic_flux_density)

Returns the gradient of the defect with respect to susceptibility.

Parameters:

  • k (Scalar) – susceptibility

  • phi (Scalar) – corresponding potential

  • magnetic_flux_density (Vector) – magnetic field

Return type:

Scalar

getPotential(k)

Calculates the magnetic potential for a given susceptibility.

Parameters:

k (Scalar) – susceptibility

Returns:

magnetic potential

Return type:

Scalar

rescaleWeights(scale=1.0, k_scale=1.0)

rescales the weights such that

sum_s integrate( ( w_i[s] *B_i[s]) (w_j[s]*1/L_j) * L**2 * (background_magnetic_flux_density_j[s]*1/L_j) * k_scale )=scale

Parameters:

  • scale (positive float) – scale of data weighting factors

  • k_scale (Scalar) – scale of susceptibility.

class esys.downunder.forwardmodels.SelfDemagnetizationModel(domain, w, B, background_magnetic_flux_density, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Forward Model for magnetic inversion with self-demagnetization as described in the inversion cookbook.

__init__(domain, w, B, background_magnetic_flux_density, coordinates=None, fixPotentialAtBottom=False, tol=1e-08)

Creates a new magnetic model on the given domain with one or more surveys (w, B).

Parameters:

  • domain (Domain) – domain of the model

  • w (Vector or list of Vector) – data weighting factors

  • B (Vector or list of Vector) – magnetic field data

  • background_magnetic_flux_density (Vector or list of float) – background magnetic flux density (in Tesla) with components (B_east, B_north, B_vertical)

  • coordinates (ReferenceSystem or SpatialCoordinateTransformation) – defines coordinate system to be used

  • fixPotentialAtBottom (bool) – if true potential is fixed to zero at the bottom of the domain in addition to the top

  • tol (positive float) – tolerance of underlying PDE

getArguments(k)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

k (Scalar) – susceptibility

Returns:

scalar magnetic potential and corresponding magnetic field

Return type:

Scalar, Vector

getDefect(k, phi, grad_phi, magnetic_flux_density)

Returns the value of the defect.

Parameters:

  • k (Scalar) – susceptibility

  • phi (Scalar) – corresponding potential

  • magnetic_flux_density (Vector) – magnetic field

Return type:

float

getGradient(k, phi, grad_phi, magnetic_flux_density)

Returns the gradient of the defect with respect to susceptibility.

Parameters:

  • k (Scalar) – susceptibility

  • phi (Scalar) – corresponding potential

  • magnetic_flux_density (Vector) – magnetic field

Return type:

Scalar

getPotential(k)

Calculates the magnetic potential for a given susceptibility.

Parameters:

k (Scalar) – susceptibility

Returns:

magnetic potential

Return type:

Scalar

rescaleWeights(scale=1.0, k_scale=1.0)

rescales the weights such that

sum_s integrate( ( w_i[s] *B_i[s]) (w_j[s]*1/L_j) * L**2 * (background_magnetic_flux_density_j[s]*1/L_j) * k_scale )=scale

Parameters:

  • scale (positive float) – scale of data weighting factors

  • k_scale (Scalar) – scale of susceptibility.

class esys.downunder.forwardmodels.Subsidence(domain, w, d, lam, mu, coordinates=None, tol=1e-08)

Forward Model for subsidence inversion minimizing integrate( (inner(w,u)-d)**2) where u is the surface displacement due to a pressure change P

__init__(domain, w, d, lam, mu, coordinates=None, tol=1e-08)

Creates a new subsidence on the given domain

Parameters:

  • domain (Domain) – domain of the model

  • w (Vector with FunctionOnBoundary) – data weighting factors and direction

  • d (Scalar with FunctionOnBoundary) – displacement measured at surface

  • lam (Scalar with Function) – 1st Lame coefficient

  • lam – 2st Lame coefficient/Shear modulus

  • coordinates (ReferenceSystem or SpatialCoordinateTransformation) – defines coordinate system to be used (not supported yet))

  • tol (positive float) – tolerance of underlying PDE

getArguments(P)

Returns precomputed values shared by getDefect() and getGradient().

Parameters:

P (Scalar) – pressure

Returns:

displacement u

Return type:

Vector

getDefect(P, u)

Returns the value of the defect.

Parameters:

  • P (Scalar) – pressure

  • u (Vector) – corresponding displacement

Return type:

float

getGradient(P, u)

Returns the gradient of the defect with respect to susceptibility.

Parameters:

  • P (Scalar) – pressure

  • u (Vector) – corresponding displacement

Return type:

Scalar

rescaleWeights(scale=1.0, P_scale=1.0)

rescales the weights

Parameters:

  • scale (positive float) – scale of data weighting factors

  • P_scale (Scalar) – scale of pressure increment

Functions

Others

Packages