wassname/simpeg
1
0
Fork 0
mirror of https://github.com/wassname/simpeg.git synced 2026-06-27 20:53:38 +08:00
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #313 from simpeg/em/dev

Browse Source 
Em/dev: Naming conventions
This commit is contained in:
Lindsey
2016-05-25 10:53:43 -07:00
parent e5ec512517 2c87a50d29
commit d98eef2560
14 changed files with 235 additions and 192 deletions
Download Patch File Download Diff File Expand all files Collapse all files
SimPEG/EM/Base.py
+18
View File
@@ -61,6 +61,15 @@ class BaseEMProblem(Problem.BaseProblem):
self._Me = self.mesh.getEdgeInnerProduct()
return self._Me
@property
def MeI(self):
"""
Edge inner product matrix
"""
if getattr(self, '_MeI', None) is None:
self._MeI = self.mesh.getEdgeInnerProduct(invMat=True)
return self._MeI
@property
def Mf(self):
"""
@@ -70,6 +79,15 @@ class BaseEMProblem(Problem.BaseProblem):
self._Mf = self.mesh.getFaceInnerProduct()
return self._Mf
@property
def MfI(self):
"""
Face inner product matrix
"""
if getattr(self, '_MfI', None) is None:
self._MfI = self.mesh.getFaceInnerProduct(invMat=True)
return self._MfI
# ----- Magnetic Permeability ----- #
@property
SimPEG/EM/FDEM/FieldsFDEM.py
+8 -8
View File
@@ -160,9 +160,9 @@ class Fields(SimPEG.Problem.Fields):
return self._jDeriv_u(src, v, adjoint), self._jDeriv_m(src, v, adjoint)
return np.array(self._jDeriv_u(src, du_dm_v, adjoint) + self._jDeriv_m(src, v, adjoint), dtype = complex)
class Fields_e(Fields):
class Fields3D_e(Fields):
"""
Fields object for Problem_e.
Fields object for Problem3D_e.
:param Mesh mesh: mesh
:param Survey survey: survey
@@ -426,9 +426,9 @@ class Fields_e(Fields):
class Fields_b(Fields):
class Fields3D_b(Fields):
"""
Fields object for Problem_b.
Fields object for Problem3D_b.
:param Mesh mesh: mesh
:param Survey survey: survey
@@ -693,9 +693,9 @@ class Fields_b(Fields):
return Zero()
class Fields_j(Fields):
class Fields3D_j(Fields):
"""
Fields object for Problem_j.
Fields object for Problem3D_j.
:param Mesh mesh: mesh
:param Survey survey: survey
@@ -988,9 +988,9 @@ class Fields_j(Fields):
return 1./(1j * omega(src.freq)) * VI * (self._aveE2CCV * ( s_mDeriv(v) - self._edgeCurl.T * ( self._MfRhoDeriv(jSolution) * v ) ) )
class Fields_h(Fields):
class Fields3D_h(Fields):
"""
Fields object for Problem_h.
Fields object for Problem3D_h.
:param Mesh mesh: mesh
:param Survey survey: survey
SimPEG/EM/FDEM/FDEM.py → SimPEG/EM/FDEM/ProblemFDEM.py
+16 -17
View File
@@ -1,7 +1,7 @@
from SimPEG import Problem, Utils, np, sp, Solver as SimpegSolver
from scipy.constants import mu_0
from SurveyFDEM import Survey as SurveyFDEM
from FieldsFDEM import Fields, Fields_e, Fields_b, Fields_h, Fields_j
from FieldsFDEM import Fields, Fields3D_e, Fields3D_b, Fields3D_h, Fields3D_j
from SimPEG.EM.Base import BaseEMProblem
from SimPEG.EM.Utils import omega
@@ -17,8 +17,8 @@ class BaseFDEMProblem(BaseEMProblem):
\mathbf{C} \mathbf{e} + i \omega \mathbf{b} = \mathbf{s_m} \\\\
{\mathbf{C}^{\\top} \mathbf{M_{\mu^{-1}}^f} \mathbf{b} - \mathbf{M_{\sigma}^e} \mathbf{e} = \mathbf{s_e}}
if using the E-B formulation (:code:`Problem_e`
or :code:`Problem_b`). Note that in this case, :math:`\mathbf{s_e}` is an integrated quantity.
if using the E-B formulation (:code:`Problem3D_e`
or :code:`Problem3D_b`). Note that in this case, :math:`\mathbf{s_e}` is an integrated quantity.
If we write Maxwell's equations in terms of
\\\(\\\mathbf{h}\\\) and current density \\\(\\\mathbf{j}\\\)
@@ -28,7 +28,7 @@ class BaseFDEMProblem(BaseEMProblem):
\mathbf{C}^{\\top} \mathbf{M_{\\rho}^f} \mathbf{j} + i \omega \mathbf{M_{\mu}^e} \mathbf{h} = \mathbf{s_m} \\\\
\mathbf{C} \mathbf{h} - \mathbf{j} = \mathbf{s_e}
if using the H-J formulation (:code:`Problem_j` or :code:`Problem_h`). Note that here, :math:`\mathbf{s_m}` is an integrated quantity.
if using the H-J formulation (:code:`Problem3D_j` or :code:`Problem3D_h`). Note that here, :math:`\mathbf{s_m}` is an integrated quantity.
The problem performs the elimination so that we are solving the system for \\\(\\\mathbf{e},\\\mathbf{b},\\\mathbf{j} \\\) or \\\(\\\mathbf{h}\\\)
"""
@@ -87,7 +87,7 @@ class BaseFDEMProblem(BaseEMProblem):
du_dm_v = Ainv * ( - dA_dm_v + dRHS_dm_v )
for rx in src.rxList:
df_dmFun = getattr(f, '_%sDeriv'%rx.projField, None)
df_dmFun = getattr(f, '_{0}Deriv'.format(rx.projField), None)
df_dm_v = df_dmFun(src, du_dm_v, v, adjoint=False)
Jv[src, rx] = rx.evalDeriv(src, self.mesh, f, df_dm_v)
Ainv.clean()
@@ -125,7 +125,7 @@ class BaseFDEMProblem(BaseEMProblem):
for rx in src.rxList:
PTv = rx.evalDeriv(src, self.mesh, f, v[src, rx], adjoint=True) # wrt f, need possibility wrt m
df_duTFun = getattr(f, '_%sDeriv'%rx.projField, None)
df_duTFun = getattr(f, '_{0}Deriv'.format(rx.projField), None)
df_duT, df_dmT = df_duTFun(src, None, PTv, adjoint=True)
ATinvdf_duT = ATinv * df_duT
@@ -137,10 +137,9 @@ class BaseFDEMProblem(BaseEMProblem):
df_dmT = df_dmT + du_dmT
# TODO: this should be taken care of by the reciever?
real_or_imag = rx.projComp
if real_or_imag is 'real':
if rx.component is 'real':
Jtv += np.array(df_dmT, dtype=complex).real
elif real_or_imag is 'imag':
elif rx.component is 'imag':
Jtv += - np.array(df_dmT, dtype=complex).real
else:
raise Exception('Must be real or imag')
@@ -177,7 +176,7 @@ class BaseFDEMProblem(BaseEMProblem):
################################ E-B Formulation #########################################
##########################################################################################
class Problem_e(BaseFDEMProblem):
class Problem3D_e(BaseFDEMProblem):
"""
By eliminating the magnetic flux density using
@@ -199,7 +198,7 @@ class Problem_e(BaseFDEMProblem):
_solutionType = 'eSolution'
_formulation = 'EB'
fieldsPair = Fields_e
fieldsPair = Fields3D_e
def __init__(self, mesh, **kwargs):
BaseFDEMProblem.__init__(self, mesh, **kwargs)
@@ -288,7 +287,7 @@ class Problem_e(BaseFDEMProblem):
return C.T * (MfMui * s_mDeriv(v)) -1j * omega(freq) * s_eDeriv(v)
class Problem_b(BaseFDEMProblem):
class Problem3D_b(BaseFDEMProblem):
"""
We eliminate :math:`\mathbf{e}` using
@@ -310,7 +309,7 @@ class Problem_b(BaseFDEMProblem):
_solutionType = 'bSolution'
_formulation = 'EB'
fieldsPair = Fields_b
fieldsPair = Fields3D_b
def __init__(self, mesh, **kwargs):
BaseFDEMProblem.__init__(self, mesh, **kwargs)
@@ -436,7 +435,7 @@ class Problem_b(BaseFDEMProblem):
##########################################################################################
class Problem_j(BaseFDEMProblem):
class Problem3D_j(BaseFDEMProblem):
"""
We eliminate \\\(\\\mathbf{h}\\\) using
@@ -458,7 +457,7 @@ class Problem_j(BaseFDEMProblem):
_solutionType = 'jSolution'
_formulation = 'HJ'
fieldsPair = Fields_j
fieldsPair = Fields3D_j
def __init__(self, mesh, **kwargs):
BaseFDEMProblem.__init__(self, mesh, **kwargs)
@@ -577,7 +576,7 @@ class Problem_j(BaseFDEMProblem):
class Problem_h(BaseFDEMProblem):
class Problem3D_h(BaseFDEMProblem):
"""
We eliminate \\\(\\\mathbf{j}\\\) using
@@ -596,7 +595,7 @@ class Problem_h(BaseFDEMProblem):
_solutionType = 'hSolution'
_formulation = 'HJ'
fieldsPair = Fields_h
fieldsPair = Fields3D_h
def __init__(self, mesh, **kwargs):
BaseFDEMProblem.__init__(self, mesh, **kwargs)
SimPEG/EM/FDEM/RxFDEM.py
+126
View File
@@ -0,0 +1,126 @@
import SimPEG
from SimPEG import sp
class BaseRx(SimPEG.Survey.BaseRx):
"""
Frequency domain receiver base class
:param numpy.ndarray locs: receiver locations (ie. :code:`np.r_[x,y,z]`)
:param string orientation: receiver orientation 'x', 'y' or 'z'
:param string component: real or imaginary component 'real' or 'imag'
"""
def __init__(self, locs, orientation=None, component=None):
assert(orientation in ['x','y','z']), "Orientation %s not known. Orientation must be in 'x', 'y', 'z'. Arbitrary orientations have not yet been implemented."%orientation
assert(component in ['real', 'imag']), "'component' must be 'real' or 'imag', not %s"%component
self.projComp = orientation
self.component = component
SimPEG.Survey.BaseRx.__init__(self, locs, rxType=None) #TODO: remove rxType from baseRx
def projGLoc(self, u):
"""Grid Location projection (e.g. Ex Fy ...)"""
return u._GLoc(self.projField) + self.projComp
def eval(self, src, mesh, f):
"""
Project fields to recievers to get data.
:param Source src: FDEM source
:param Mesh mesh: mesh used
:param Fields f: fields object
:rtype: numpy.ndarray
:return: fields projected to recievers
"""
P = self.getP(mesh, self.projGLoc(f))
f_part_complex = f[src, self.projField]
f_part = getattr(f_part_complex, self.component) # get the real or imag component
return P*f_part
def evalDeriv(self, src, mesh, f, v, adjoint=False):
"""
Derivative of projected fields with respect to the inversion model times a vector.
:param Source src: FDEM source
:param Mesh mesh: mesh used
:param Fields f: fields object
:param numpy.ndarray v: vector to multiply
:rtype: numpy.ndarray
:return: fields projected to recievers
"""
P = self.getP(mesh, self.projGLoc(f))
if not adjoint:
Pv_complex = P * v
Pv = getattr(Pv_complex, self.component)
elif adjoint:
Pv_real = P.T * v
if self.component == 'imag':
Pv = 1j*Pv_real
elif self.component == 'real':
Pv = Pv_real.astype(complex)
else:
raise NotImplementedError('must be real or imag')
return Pv
class Point_e(BaseRx):
"""
Electric field FDEM receiver
:param numpy.ndarray locs: receiver locations (ie. :code:`np.r_[x,y,z]`)
:param string orientation: receiver orientation 'x', 'y' or 'z'
:param string component: real or imaginary component 'real' or 'imag'
"""
def __init__(self, locs, orientation=None, component=None):
self.projField = 'e'
super(Point_e, self).__init__(locs, orientation, component)
class Point_b(BaseRx):
"""
Magnetic flux FDEM receiver
:param numpy.ndarray locs: receiver locations (ie. :code:`np.r_[x,y,z]`)
:param string orientation: receiver orientation 'x', 'y' or 'z'
:param string component: real or imaginary component 'real' or 'imag'
"""
def __init__(self, locs, orientation=None, component=None):
self.projField = 'b'
super(Point_b, self).__init__(locs, orientation, component)
class Point_h(BaseRx):
"""
Magnetic field FDEM receiver
:param numpy.ndarray locs: receiver locations (ie. :code:`np.r_[x,y,z]`)
:param string orientation: receiver orientation 'x', 'y' or 'z'
:param string component: real or imaginary component 'real' or 'imag'
"""
def __init__(self, locs, orientation=None, component=None):
self.projField = 'h'
super(Point_h, self).__init__(locs, orientation, component)
class Point_j(BaseRx):
"""
Current density FDEM receiver
:param numpy.ndarray locs: receiver locations (ie. :code:`np.r_[x,y,z]`)
:param string orientation: receiver orientation 'x', 'y' or 'z'
:param string component: real or imaginary component 'real' or 'imag'
"""
def __init__(self, locs, orientation=None, component=None):
self.projField = 'j'
super(Point_j, self).__init__(locs, orientation, component)
SimPEG/EM/FDEM/SrcFDEM.py
+31 -21
View File
@@ -9,8 +9,14 @@ class BaseSrc(Survey.BaseSrc):
"""
freq = None
# rxPair = RxFDEM
integrate = True
integrate = False
_ePrimary = None
_bPrimary = None
_hPrimary = None
_jPrimary = None
def __init__(self, rxList, **kwargs):
Survey.BaseSrc.__init__(self, rxList, **kwargs)
def eval(self, prob):
"""
@@ -51,7 +57,9 @@ class BaseSrc(Survey.BaseSrc):
:rtype: numpy.ndarray
:return: primary magnetic flux density
"""
return Zero()
if self._bPrimary is None:
return Zero()
return self._bPrimary
def hPrimary(self, prob):
"""
@@ -61,7 +69,9 @@ class BaseSrc(Survey.BaseSrc):
:rtype: numpy.ndarray
:return: primary magnetic field
"""
return Zero()
if self._hPrimary is None:
return Zero()
return self._hPrimary
def ePrimary(self, prob):
"""
@@ -71,7 +81,9 @@ class BaseSrc(Survey.BaseSrc):
:rtype: numpy.ndarray
:return: primary electric field
"""
return Zero()
if self._ePrimary is None:
return Zero()
return self._ePrimary
def jPrimary(self, prob):
"""
@@ -81,7 +93,9 @@ class BaseSrc(Survey.BaseSrc):
:rtype: numpy.ndarray
:return: primary current density
"""
return Zero()
if self._jPrimary is None:
return Zero()
return self._jPrimary
def s_m(self, prob):
"""
@@ -136,15 +150,14 @@ class RawVec_e(BaseSrc):
:param list rxList: receiver list
:param float freq: frequency
:param numpy.array s_e: electric source term
:param bool integrate: Integrate the source term (multiply by Me) [True]
:param bool integrate: Integrate the source term (multiply by Me) [False]
"""
def __init__(self, rxList, freq, s_e, integrate=True): #, ePrimary=None, bPrimary=None, hPrimary=None, jPrimary=None):
def __init__(self, rxList, freq, s_e, **kwargs):
self._s_e = np.array(s_e, dtype=complex)
self.freq = float(freq)
self.integrate = integrate
BaseSrc.__init__(self, rxList)
BaseSrc.__init__(self, rxList, **kwargs)
def s_e(self, prob):
"""
@@ -166,15 +179,14 @@ class RawVec_m(BaseSrc):
:param float freq: frequency
:param rxList: receiver list
:param numpy.array s_m: magnetic source term
:param bool integrate: Integrate the source term (multiply by Me) [True]
:param bool integrate: Integrate the source term (multiply by Me) [False]
"""
def __init__(self, rxList, freq, s_m, integrate=True): #ePrimary=Zero(), bPrimary=Zero(), hPrimary=Zero(), jPrimary=Zero()):
def __init__(self, rxList, freq, s_m, **kwargs): #ePrimary=Zero(), bPrimary=Zero(), hPrimary=Zero(), jPrimary=Zero()):
self._s_m = np.array(s_m, dtype=complex)
self.freq = float(freq)
self.integrate = integrate
BaseSrc.__init__(self, rxList)
BaseSrc.__init__(self, rxList, **kwargs)
def s_m(self, prob):
"""
@@ -197,14 +209,13 @@ class RawVec(BaseSrc):
:param float freq: frequency
:param numpy.array s_m: magnetic source term
:param numpy.array s_e: electric source term
:param bool integrate: Integrate the source term (multiply by Me) [True]
:param bool integrate: Integrate the source term (multiply by Me) [False]
"""
def __init__(self, rxList, freq, s_m, s_e, integrate=True):
def __init__(self, rxList, freq, s_m, s_e, **kwargs):
self._s_m = np.array(s_m, dtype=complex)
self._s_e = np.array(s_e, dtype=complex)
self.freq = float(freq)
self.integrate = integrate
BaseSrc.__init__(self, rxList)
BaseSrc.__init__(self, rxList, **kwargs)
def s_m(self, prob):
"""
@@ -278,14 +289,13 @@ class MagDipole(BaseSrc):
:param float mu: background magnetic permeability
"""
def __init__(self, rxList, freq, loc, orientation='Z', moment=1., mu=mu_0):
def __init__(self, rxList, freq, loc, orientation='Z', moment=1., mu=mu_0, **kwargs):
self.freq = float(freq)
self.loc = loc
self.orientation = orientation
assert orientation in ['X','Y','Z'], "Orientation (right now) doesn't actually do anything! The methods in SrcUtils should take care of this..."
self.moment = moment
self.mu = mu
self.integrate = False
BaseSrc.__init__(self, rxList)
def bPrimary(self, prob):
@@ -543,7 +553,7 @@ class CircularLoop(BaseSrc):
if not prob.mesh.isSymmetric:
# TODO ?
raise NotImplementedError('Non-symmetric cyl mesh not implemented yet!')
a = MagneticDipoleVectorPotential(self.loc, gridY, 'y', moment=self.radius, mu=self.mu)
a = MagneticLoopVectorPotential(self.loc, gridY, 'y', moment=self.radius, mu=self.mu)
else:
srcfct = MagneticDipoleVectorPotential
SimPEG/EM/FDEM/SurveyFDEM.py
+2 -119
View File
@@ -4,126 +4,9 @@ from SimPEG.EM.Base import BaseEMSurvey
from scipy.constants import mu_0
from SimPEG.Utils import Zero, Identity
import SrcFDEM as Src
import RxFDEM as Rx
from SimPEG import sp
####################################################
# Receivers
####################################################
class Rx(SimPEG.Survey.BaseRx):
"""
Frequency domain receivers
:param numpy.ndarray locs: receiver locations (ie. :code:`np.r_[x,y,z]`)
:param string rxType: reciever type from knownRxTypes
"""
knownRxTypes = {
'exr':['e', 'x', 'real'],
'eyr':['e', 'y', 'real'],
'ezr':['e', 'z', 'real'],
'exi':['e', 'x', 'imag'],
'eyi':['e', 'y', 'imag'],
'ezi':['e', 'z', 'imag'],
'bxr':['b', 'x', 'real'],
'byr':['b', 'y', 'real'],
'bzr':['b', 'z', 'real'],
'bxi':['b', 'x', 'imag'],
'byi':['b', 'y', 'imag'],
'bzi':['b', 'z', 'imag'],
'jxr':['j', 'x', 'real'],
'jyr':['j', 'y', 'real'],
'jzr':['j', 'z', 'real'],
'jxi':['j', 'x', 'imag'],
'jyi':['j', 'y', 'imag'],
'jzi':['j', 'z', 'imag'],
'hxr':['h', 'x', 'real'],
'hyr':['h', 'y', 'real'],
'hzr':['h', 'z', 'real'],
'hxi':['h', 'x', 'imag'],
'hyi':['h', 'y', 'imag'],
'hzi':['h', 'z', 'imag'],
}
radius = None
def __init__(self, locs, rxType):
SimPEG.Survey.BaseRx.__init__(self, locs, rxType)
@property
def projField(self):
"""Field Type projection (e.g. e b ...)"""
return self.knownRxTypes[self.rxType][0]
@property
def projComp(self):
"""Component projection (real/imag)"""
return self.knownRxTypes[self.rxType][2]
def projGLoc(self, u):
"""Grid Location projection (e.g. Ex Fy ...)"""
return u._GLoc(self.rxType[0]) + self.knownRxTypes[self.rxType][1]
def eval(self, src, mesh, f):
"""
Project fields to recievers to get data.
:param Source src: FDEM source
:param Mesh mesh: mesh used
:param Fields f: fields object
:rtype: numpy.ndarray
:return: fields projected to recievers
"""
# projGLoc = u._GLoc(self.knownRxTypes[self.rxType][0])
# projGLoc += self.knownRxTypes[self.rxType][1]
P = self.getP(mesh, self.projGLoc(f))
f_part_complex = f[src, self.projField]
# get the real or imag component
real_or_imag = self.projComp
f_part = getattr(f_part_complex, real_or_imag)
return P*f_part
def evalDeriv(self, src, mesh, f, v, adjoint=False):
"""
Derivative of projected fields with respect to the inversion model times a vector.
:param Source src: FDEM source
:param Mesh mesh: mesh used
:param Fields f: fields object
:param numpy.ndarray v: vector to multiply
:rtype: numpy.ndarray
:return: fields projected to recievers
"""
P = self.getP(mesh, self.projGLoc(f))
if not adjoint:
Pv_complex = P * v
real_or_imag = self.projComp
Pv = getattr(Pv_complex, real_or_imag)
elif adjoint:
Pv_real = P.T * v
real_or_imag = self.projComp
if real_or_imag == 'imag':
Pv = 1j*Pv_real
elif real_or_imag == 'real':
Pv = Pv_real.astype(complex)
else:
raise NotImplementedError('must be real or imag')
return Pv
####################################################
# Survey
####################################################
class Survey(BaseEMSurvey):
"""
Frequency domain electromagnetic survey
@@ -132,7 +15,7 @@ class Survey(BaseEMSurvey):
"""
srcPair = Src.BaseSrc
rxPair = Rx
rxPair = Rx.BaseRx
def __init__(self, srcList, **kwargs):
# Sort these by frequency
SimPEG/EM/FDEM/__init__.py
+5 -3
View File
@@ -1,3 +1,5 @@
from SurveyFDEM import Rx, Src, Survey
from FDEM import BaseFDEMProblem, Problem_e, Problem_b, Problem_j, Problem_h
from FieldsFDEM import *
from SurveyFDEM import Survey
import SrcFDEM as Src
import RxFDEM as Rx
from ProblemFDEM import Problem3D_e, Problem3D_b, Problem3D_j, Problem3D_h
from FieldsFDEM import Fields3D_e, Fields3D_b, Fields3D_j, Fields3D_h
SimPEG/EM/Utils/testingUtils.py
+17 -12
View File
@@ -20,56 +20,61 @@ def getFDEMProblem(fdemType, comp, SrcList, freq, useMu=False, verbose=False):
mesh = Mesh.TensorMesh([hx,hy,hz],['C','C','C'])
if useMu is True:
mapping = [('sigma', Maps.ExpMap(mesh)), ('mu', Maps.IdentityMap(mesh))]
mapping = [('sigma', Maps.ExpMap(mesh)), ('mu', Maps.IdentityMap(mesh))]
else:
mapping = Maps.ExpMap(mesh)
x = np.array([np.linspace(-5.*cs,-2.*cs,3),np.linspace(5.*cs,2.*cs,3)]) + cs/4. #don't sample right by the source, slightly off alignment from either staggered grid
XYZ = Utils.ndgrid(x,x,np.linspace(-2.*cs,2.*cs,5))
Rx0 = EM.FDEM.Rx(XYZ, comp)
Rx0 = getattr(EM.FDEM.Rx, 'Point_' + comp[0])
if comp[2] == 'r':
real_or_imag = 'real'
elif comp[2] == 'i':
real_or_imag = 'imag'
rx0 = Rx0(XYZ, comp[1], 'imag')
Src = []
for SrcType in SrcList:
if SrcType is 'MagDipole':
Src.append(EM.FDEM.Src.MagDipole([Rx0], freq=freq, loc=np.r_[0.,0.,0.]))
Src.append(EM.FDEM.Src.MagDipole([rx0], freq=freq, loc=np.r_[0.,0.,0.]))
elif SrcType is 'MagDipole_Bfield':
Src.append(EM.FDEM.Src.MagDipole_Bfield([Rx0], freq=freq, loc=np.r_[0.,0.,0.]))
Src.append(EM.FDEM.Src.MagDipole_Bfield([rx0], freq=freq, loc=np.r_[0.,0.,0.]))
elif SrcType is 'CircularLoop':
Src.append(EM.FDEM.Src.CircularLoop([Rx0], freq=freq, loc=np.r_[0.,0.,0.]))
Src.append(EM.FDEM.Src.CircularLoop([rx0], freq=freq, loc=np.r_[0.,0.,0.]))
elif SrcType is 'RawVec':
if fdemType is 'e' or fdemType is 'b':
S_m = np.zeros(mesh.nF)
S_e = np.zeros(mesh.nE)
S_m[Utils.closestPoints(mesh,[0.,0.,0.],'Fz') + np.sum(mesh.vnF[:1])] = 1e-3
S_e[Utils.closestPoints(mesh,[0.,0.,0.],'Ez') + np.sum(mesh.vnE[:1])] = 1e-3
Src.append(EM.FDEM.Src.RawVec([Rx0], freq, S_m, S_e))
Src.append(EM.FDEM.Src.RawVec([rx0], freq, S_m, mesh.getEdgeInnerProduct()*S_e))
elif fdemType is 'h' or fdemType is 'j':
S_m = np.zeros(mesh.nE)
S_e = np.zeros(mesh.nF)
S_m[Utils.closestPoints(mesh,[0.,0.,0.],'Ez') + np.sum(mesh.vnE[:1])] = 1e-3
S_e[Utils.closestPoints(mesh,[0.,0.,0.],'Fz') + np.sum(mesh.vnF[:1])] = 1e-3
Src.append(EM.FDEM.Src.RawVec([Rx0], freq, S_m, S_e))
Src.append(EM.FDEM.Src.RawVec([rx0], freq, mesh.getEdgeInnerProduct()*S_m, S_e))
if verbose:
print ' Fetching %s problem' % (fdemType)
if fdemType == 'e':
survey = EM.FDEM.Survey(Src)
prb = EM.FDEM.Problem_e(mesh, mapping=mapping)
prb = EM.FDEM.Problem3D_e(mesh, mapping=mapping)
elif fdemType == 'b':
survey = EM.FDEM.Survey(Src)
prb = EM.FDEM.Problem_b(mesh, mapping=mapping)
prb = EM.FDEM.Problem3D_b(mesh, mapping=mapping)
elif fdemType == 'j':
survey = EM.FDEM.Survey(Src)
prb = EM.FDEM.Problem_j(mesh, mapping=mapping)
prb = EM.FDEM.Problem3D_j(mesh, mapping=mapping)
elif fdemType == 'h':
survey = EM.FDEM.Survey(Src)
prb = EM.FDEM.Problem_h(mesh, mapping=mapping)
prb = EM.FDEM.Problem3D_h(mesh, mapping=mapping)
else:
raise NotImplementedError()
@@ -90,7 +95,7 @@ def crossCheckTest(SrcList, fdemType1, fdemType2, comp, addrandoms = False, useM
prb1 = getFDEMProblem(fdemType1, comp, SrcList, freq, useMu, verbose)
mesh = prb1.mesh
print 'Cross Checking Forward: %s, %s formulations - %s' % (fdemType1, fdemType2, comp)
logsig = np.log(np.ones(mesh.nC)*CONDUCTIVITY)
mu = np.ones(mesh.nC)*MU
SimPEG/Examples/EM_FDEM_1D_Inversion.py
+3 -3
View File
@@ -42,8 +42,8 @@ def run(plotIt=True):
ax.grid(color='k', alpha=0.5, linestyle='dashed', linewidth=0.5)
rxOffset=10.
bzi = EM.FDEM.Rx(np.array([[rxOffset, 0., 1e-3]]), 'bzi')
rxOffset=10.
bzi = EM.FDEM.Rx.Point_b(np.array([[rxOffset, 0., 1e-3]]), orientation='z', component='imag')
freqs = np.logspace(1,3,10)
srcLoc = np.array([0., 0., 10.])
@@ -51,7 +51,7 @@ def run(plotIt=True):
srcList = [EM.FDEM.Src.MagDipole([bzi],freq, srcLoc,orientation='Z') for freq in freqs]
survey = EM.FDEM.Survey(srcList)
prb = EM.FDEM.Problem_b(mesh, mapping=mapping)
prb = EM.FDEM.Problem3D_b(mesh, mapping=mapping)
try:
from pymatsolver import MumpsSolver
SimPEG/Examples/EM_Schenkel_Morrison_Casing.py
+1 -1
View File
@@ -215,7 +215,7 @@ def run(plotIt=True):
# ------------ Problem and Survey ---------------
survey = FDEM.Survey(sg_p + dg_p)
mapping = [('sigma', Maps.IdentityMap(mesh))]
problem = FDEM.Problem_h(mesh, mapping=mapping)
problem = FDEM.Problem3D_h(mesh, mapping=mapping)
problem.pair(survey)
# ------------- Solve ---------------------------
SimPEG/MT/BaseMT.py
+1 -1
View File
@@ -1,5 +1,5 @@
from SimPEG import SolverLU as SimpegSolver, PropMaps, Utils, mkvc, sp, np
from SimPEG.EM.FDEM.FDEM import BaseFDEMProblem
from SimPEG.EM.FDEM.ProblemFDEM import BaseFDEMProblem
from SurveyMT import Survey, Data
from FieldsMT import BaseMTFields
docs/examples/DC_Forward_PseudoSection.rst
+1 -1
View File
@@ -22,7 +22,7 @@ radi = Radius of spheres [r1,r2]
param = Conductivity of background and two spheres [m0,m1,m2]
stype = survey type "pdp" (pole dipole) or "dpdp" (dipole dipole)
dtype = Data type "appr" (app res) | "appc" (app cond) | "volt" (potential)
Created by @fourndo on Mon Feb 01 19:28:06 2016
Created by @fourndo
tests/em/fdem/forward/test_FDEM_analytics.py
+4 -4
View File
@@ -28,12 +28,12 @@ class FDEM_analyticTests(unittest.TestCase):
x = np.linspace(-10,10,5)
XYZ = Utils.ndgrid(x,np.r_[0],np.r_[0])
rxList = EM.FDEM.Rx(XYZ, 'exi')
rxList = EM.FDEM.Rx.Point_e(XYZ, orientation='x', component='imag')
Src0 = EM.FDEM.Src.MagDipole([rxList],loc=np.r_[0.,0.,0.], freq=freq)
survey = EM.FDEM.Survey([Src0])
prb = EM.FDEM.Problem_b(mesh, mapping=mapping)
prb = EM.FDEM.Problem3D_b(mesh, mapping=mapping)
prb.pair(survey)
try:
@@ -125,8 +125,8 @@ class FDEM_analyticTests(unittest.TestCase):
mapping = [('sigma', Maps.IdentityMap(mesh)),('mu', Maps.IdentityMap(mesh))]
prbe = EM.FDEM.Problem_h(mesh, mapping=mapping)
prbm = EM.FDEM.Problem_e(mesh, mapping=mapping)
prbe = EM.FDEM.Problem3D_h(mesh, mapping=mapping)
prbm = EM.FDEM.Problem3D_e(mesh, mapping=mapping)
prbe.pair(surveye) # pair problem and survey
prbm.pair(surveym)
tests/em/fdem/forward/test_FDEM_forwardHB.py
+2 -2
View File
@@ -12,7 +12,7 @@ testBH = True
verbose = False
TOLEJHB = 1 # averaging and more sensitive to boundary condition violations (ie. the impact of violating the boundary conditions in each case is different.)
#TODO: choose better testing parameters to lower this
#TODO: choose better testing parameters to lower this
SrcList = ['RawVec', 'MagDipole_Bfield', 'MagDipole', 'CircularLoop']
@@ -125,4 +125,4 @@ class FDEM_CrossCheck(unittest.TestCase):
self.assertTrue(crossCheckTest(SrcList, 'b', 'h', 'hzi', verbose=verbose, TOL=TOLEJHB))
if __name__ == '__main__':
unittest.main()
unittest.main()