esys.downunder.forwardmodels.magnetotelluric2d Package¶
Forward models for 2D MT (TE and TM mode)
Classes¶
-
class
esys.downunder.forwardmodels.magnetotelluric2d.
ForwardModel
¶ Bases:
object
An abstract forward model that can be plugged into a cost function. Subclasses need to implement
getDefect()
,getGradient()
, and possiblygetArguments()
and ‘getCoordinateTransformation’.-
__init__
()¶ Initialize self. See help(type(self)) for accurate signature.
-
getArguments
(x)¶
-
getCoordinateTransformation
()¶
-
getDefect
(x, *args)¶
-
getGradient
(x, *args)¶
-
-
class
esys.downunder.forwardmodels.magnetotelluric2d.
MT2DBase
(domain, omega, x, Z, eta=None, w0=1.0, mu=1.2566370614359173e-06, sigma0=0.01, airLayerLevel=None, fixAirLayer=False, coordinates=None, tol=1e-08, saveMemory=False, directSolver=True)¶ Bases:
esys.downunder.forwardmodels.base.ForwardModel
Base class for 2D MT forward models. See
MT2DModelTEMode
andMT2DModelTMMode
for actual implementations.-
__init__
(domain, omega, x, Z, eta=None, w0=1.0, mu=1.2566370614359173e-06, sigma0=0.01, airLayerLevel=None, fixAirLayer=False, coordinates=None, tol=1e-08, saveMemory=False, directSolver=True)¶ initializes a new forward model.
Parameters: - domain (
Domain
) – domain of the model - omega (positive
float
) – frequency - x (
list
oftuple
withfloat
) – coordinates of measurements - Z (
list
ofcomplex
) – measured impedance (possibly scaled) - eta (positive
float
orlist
of positivefloat
) – spatial confidence radius - w0 (
None
or a list of positivefloat
) – confidence factors for meassurements. - mu (
float
) – permeability - sigma0 (
float
) – background conductivity - airLayerLevel (
float
orNone
) – position of the air layer from to bottom of the domain. If not set the air layer starts at the top of the domain - fixAirLayer (
bool
) – fix air layer (TM mode) - coordinates (
ReferenceSystem
orSpatialCoordinateTransformation
) – 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 the PDE to minimize memory use. This will require more compute time as the matrix needs to be reallocated at each iteration. - directSolver (
bool
) – if true a direct solver (rather than an iterative solver) will be used to solve the PDE
- domain (
-
getArguments
(x)¶ Returns precomputed values shared by
getDefect()
andgetGradient()
. Needs to be implemented in subclasses.
-
getCoordinateTransformation
()¶ returns the coordinate transformation being used
Return type: CoordinateTransformation
-
getDefect
(x, Ex, dExdz)¶ Returns the defect value. Needs to be implemented in subclasses.
-
getDomain
()¶ Returns the domain of the forward model.
Return type: Domain
-
getGradient
(x, Ex, dExdz)¶ Returns the gradient. Needs to be implemented in subclasses.
-
getWeightingFactor
(x, wx0, x0, eta)¶ returns the weighting factor
-
setUpPDE
()¶ Return the underlying PDE.
Return type: LinearPDE
-
-
class
esys.downunder.forwardmodels.magnetotelluric2d.
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)¶ Bases:
esys.downunder.forwardmodels.magnetotelluric2d.MT2DBase
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()
andgetGradient()
.Parameters: sigma ( Data
of shape (2,)) – conductivityReturns: E_x, E_z Return type: Data
of shape (2,)
-
getCoordinateTransformation
()¶ returns the coordinate transformation being used
Return type: CoordinateTransformation
-
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
- sigma (
-
getDomain
()¶ Returns the domain of the forward model.
Return type: Domain
-
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
- sigma (
-
getWeightingFactor
(x, wx0, x0, eta)¶ returns the weighting factor
-
setUpPDE
()¶ Return the underlying PDE.
Return type: LinearPDE
-
class
esys.downunder.forwardmodels.magnetotelluric2d.
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)¶ Bases:
esys.downunder.forwardmodels.magnetotelluric2d.MT2DBase
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()
andgetGradient()
.Parameters: rho ( Data
of shape (2,)) – resistivityReturns: Hx, grad(Hx) Return type: tuple
ofData
-
getCoordinateTransformation
()¶ returns the coordinate transformation being used
Return type: CoordinateTransformation
-
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
- rho (
-
getDomain
()¶ Returns the domain of the forward model.
Return type: Domain
-
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
- rho (
-
getWeightingFactor
(x, wx0, x0, eta)¶ returns the weighting factor
-
setUpPDE
()¶ Return the underlying PDE.
Return type: LinearPDE