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

simpegEM --> SimPEG.EM

Browse Source
This commit is contained in:
Rowan Cockett
2015-11-04 10:13:24 -08:00
parent e622dfb5cc
commit 4db756b060
6 changed files with 96 additions and 96 deletions
Download Patch File Download Diff File Expand all files Collapse all files
SimPEG/EM/Analytics/FDEMcasing.py
+20 -20
View File
@@ -1,6 +1,6 @@
from SimPEG import Utils, np
from scipy.constants import mu_0, epsilon_0
from simpegEM.Utils.EMUtils import k
from SimPEG.EM.Utils.EMUtils import k
def getKc(freq,sigma,a,b,mu=mu_0,eps=epsilon_0):
a = float(a)
@@ -16,43 +16,43 @@ def _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps
nobs = obsloc.shape[0]
dxyz = obsloc - np.c_[np.ones(nobs)]*np.r_[srcloc]
r2 = _r2(dxyz[:,:2])
sqrtr2z2 = np.sqrt(r2 + dxyz[:,2]**2)
k2 = k(freq,sigma[2],mu[2],eps)
k2 = k(freq,sigma[2],mu[2],eps)
return Kc1 * moment / (4.*np.pi) *np.exp(-1j*k2*sqrtr2z2) / sqrtr2z2
def _getCasingHertzMagDipoleDeriv_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
nobs = obsloc.shape[0]
dxyz = obsloc - np.c_[np.ones(nobs)]*np.r_[srcloc]
r2 = _r2(dxyz[:,:2])
sqrtr2z2 = np.sqrt(r2 + dxyz[:,2]**2)
k2 = k(freq,sigma[2],mu[2],eps)
k2 = k(freq,sigma[2],mu[2],eps)
return -HertzZ * np.sqrt(r2) / sqrtr2z2 * (1j*k2 + 1./ sqrtr2z2)
def _getCasingHertzMagDipoleDeriv_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
nobs = obsloc.shape[0]
dxyz = obsloc - np.c_[np.ones(nobs)]*np.r_[srcloc]
r2z2 = _r2(dxyz)
sqrtr2z2 = np.sqrt(r2z2)
k2 = k(freq,sigma[2],mu[2],eps)
k2 = k(freq,sigma[2],mu[2],eps)
return -HertzZ*dxyz[:,2] /sqrtr2z2 * (1j*k2 + 1./sqrtr2z2)
def _getCasingHertzMagDipole2Deriv_z_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
dHertzZdr = _getCasingHertzMagDipoleDeriv_r(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
nobs = obsloc.shape[0]
dxyz = obsloc - np.c_[np.ones(nobs)]*np.r_[srcloc]
@@ -60,14 +60,14 @@ def _getCasingHertzMagDipole2Deriv_z_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.o
r = np.sqrt(r2)
z = dxyz[:,2]
sqrtr2z2 = np.sqrt(r2 + z**2)
k2 = k(freq,sigma[2],mu[2],eps)
k2 = k(freq,sigma[2],mu[2],eps)
return dHertzZdr*(-z/sqrtr2z2)*(1j*k2+1./sqrtr2z2) + HertzZ*(z*r/sqrtr2z2**3)*(1j*k2 + 2./sqrtr2z2)
def _getCasingHertzMagDipole2Deriv_z_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
dHertzZdz = _getCasingHertzMagDipoleDeriv_z(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
nobs = obsloc.shape[0]
dxyz = obsloc - np.c_[np.ones(nobs)]*np.r_[srcloc]
@@ -75,10 +75,10 @@ def _getCasingHertzMagDipole2Deriv_z_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.o
r = np.sqrt(r2)
z = dxyz[:,2]
sqrtr2z2 = np.sqrt(r2 + z**2)
k2 = k(freq,sigma[2],mu[2],eps)
k2 = k(freq,sigma[2],mu[2],eps)
return (dHertzZdz*z + HertzZ)/sqrtr2z2*(-1j*k2 - 1./sqrtr2z2) + HertzZ*z/sqrtr2z2**3*(1j*k2*z + 2.*z/sqrtr2z2)
def getCasingEphiMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
return 1j * omega(freq) * mu * _getCasingHertzMagDipoleDeriv_r(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
@@ -95,4 +95,4 @@ def getCasingBrMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=eps
return mu_0 * getCasingHrMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
def getCasingBzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
return mu_0 * getCasingHzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
return mu_0 * getCasingHzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
SimPEG/EM/FDEM/FDEM.py
+44 -44
View File
@@ -2,33 +2,33 @@ from SimPEG import Survey, Problem, Utils, np, sp, Solver as SimpegSolver
from scipy.constants import mu_0
from SurveyFDEM import SurveyFDEM
from FieldsFDEM import FieldsFDEM, FieldsFDEM_e, FieldsFDEM_b, FieldsFDEM_h, FieldsFDEM_j
from simpegEM.Base import BaseEMProblem
from simpegEM.Utils.EMUtils import omega
from SimPEG.EM.Base import BaseEMProblem
from SimPEG.EM.Utils.EMUtils import omega
class BaseFDEMProblem(BaseEMProblem):
"""
We start by looking at Maxwell's equations in the electric
field \\\(\\\mathbf{e}\\\) and the magnetic flux
We start by looking at Maxwell's equations in the electric
field \\\(\\\mathbf{e}\\\) and the magnetic flux
density \\\(\\\mathbf{b}\\\)
.. math ::
\mathbf{C} \mathbf{e} + i \omega \mathbf{b} = \mathbf{s_m} \\\\
{\mathbf{C}^T \mathbf{M_{\mu^{-1}}^f} \mathbf{b} - \mathbf{M_{\sigma}^e} \mathbf{e} = \mathbf{M^e} \mathbf{s_e}}
if using the E-B formulation (:code:`ProblemFDEM_e`
or :code:`ProblemFDEM_b`) or the magnetic field
if using the E-B formulation (:code:`ProblemFDEM_e`
or :code:`ProblemFDEM_b`) or the magnetic field
\\\(\\\mathbf{h}\\\) and current density \\\(\\\mathbf{j}\\\)
.. math ::
.. math ::
\mathbf{C}^T \mathbf{M_{\\rho}^f} \mathbf{j} + i \omega \mathbf{M_{\mu}^e} \mathbf{h} = \mathbf{M^e} \mathbf{s_m} \\\\
\mathbf{C} \mathbf{h} - \mathbf{j} = \mathbf{s_e}
if using the H-J formulation (:code:`ProblemFDEM_j` or :code:`ProblemFDEM_h`).
The problem performs the elimination so that we are solving the system for \\\(\\\mathbf{e},\\\mathbf{b},\\\mathbf{j} \\\) or \\\(\\\mathbf{h}\\\)
The problem performs the elimination so that we are solving the system for \\\(\\\mathbf{e},\\\mathbf{b},\\\mathbf{j} \\\) or \\\(\\\mathbf{h}\\\)
"""
surveyPair = SurveyFDEM
@@ -36,7 +36,7 @@ class BaseFDEMProblem(BaseEMProblem):
def fields(self, m=None):
"""
Solve the forward problem for the fields.
Solve the forward problem for the fields.
"""
self.curModel = m
@@ -59,7 +59,7 @@ class BaseFDEMProblem(BaseEMProblem):
"""
if f is None:
f = self.fields(m)
f = self.fields(m)
self.curModel = m
@@ -67,7 +67,7 @@ class BaseFDEMProblem(BaseEMProblem):
for freq in self.survey.freqs:
dA_du = self.getA(freq) #
dA_duI = self.Solver(dA_du, **self.solverOpts)
dA_duI = self.Solver(dA_du, **self.solverOpts)
for src in self.survey.getSrcByFreq(freq):
ftype = self._fieldType + 'Solution'
@@ -97,7 +97,7 @@ class BaseFDEMProblem(BaseEMProblem):
return Utils.mkvc(Jv)
def Jtvec(self, m, v, f=None):
def Jtvec(self, m, v, f=None):
"""
Sensitivity transpose times a vector
"""
@@ -157,7 +157,7 @@ class BaseFDEMProblem(BaseEMProblem):
def getSourceTerm(self, freq):
"""
Evaluates the sources for a given frequency and puts them in matrix form
Evaluates the sources for a given frequency and puts them in matrix form
:param float freq: Frequency
:rtype: numpy.ndarray (nE or nF, nSrc)
@@ -165,11 +165,11 @@ class BaseFDEMProblem(BaseEMProblem):
"""
Srcs = self.survey.getSrcByFreq(freq)
if self._eqLocs is 'FE':
S_m = np.zeros((self.mesh.nF,len(Srcs)), dtype=complex)
S_m = np.zeros((self.mesh.nF,len(Srcs)), dtype=complex)
S_e = np.zeros((self.mesh.nE,len(Srcs)), dtype=complex)
elif self._eqLocs is 'EF':
S_m = np.zeros((self.mesh.nE,len(Srcs)), dtype=complex)
S_e = np.zeros((self.mesh.nF,len(Srcs)), dtype=complex)
S_e = np.zeros((self.mesh.nF,len(Srcs)), dtype=complex)
for i, src in enumerate(Srcs):
smi, sei = src.eval(self)
@@ -198,9 +198,9 @@ class ProblemFDEM_e(BaseFDEMProblem):
.. math ::
\\left(\mathbf{C}^T \mathbf{M_{\mu^{-1}}^f} \mathbf{C}+ i \omega \mathbf{M^e_{\sigma}} \\right)\mathbf{e} = \mathbf{C}^T \mathbf{M_{\mu^{-1}}^f}\mathbf{s_m} -i\omega\mathbf{M^e}\mathbf{s_e}
\\left(\mathbf{C}^T \mathbf{M_{\mu^{-1}}^f} \mathbf{C}+ i \omega \mathbf{M^e_{\sigma}} \\right)\mathbf{e} = \mathbf{C}^T \mathbf{M_{\mu^{-1}}^f}\mathbf{s_m} -i\omega\mathbf{M^e}\mathbf{s_e}
which we solve for \\\(\\\mathbf{e}\\\).
which we solve for \\\(\\\mathbf{e}\\\).
"""
_fieldType = 'e'
@@ -271,10 +271,10 @@ class ProblemFDEM_e(BaseFDEMProblem):
return - 1j * omega(freq) * S_eDerivv
else:
return None
else:
else:
S_mDerivv, S_eDerivv = S_mDeriv(v), S_eDeriv(v)
if S_mDerivv is not None and S_eDerivv is not None:
if S_mDerivv is not None and S_eDerivv is not None:
return C.T * (MfMui * S_mDerivv) -1j * omega(freq) * S_eDerivv
elif S_mDerivv is not None:
return C.T * (MfMui * S_mDerivv)
@@ -289,17 +289,17 @@ class ProblemFDEM_b(BaseFDEMProblem):
We eliminate \\\(\\\mathbf{e}\\\) using
.. math ::
\mathbf{e} = \mathbf{M^e_{\sigma}}^{-1} \\left(\mathbf{C}^T \mathbf{M_{\mu^{-1}}^f} \mathbf{b} - \mathbf{s_e}\\right)
and solve for \\\(\\\mathbf{b}\\\) using:
.. math ::
\\left(\mathbf{C} \mathbf{M^e_{\sigma}}^{-1} \mathbf{C}^T \mathbf{M_{\mu^{-1}}^f} + i \omega \\right)\mathbf{b} = \mathbf{s_m} + \mathbf{M^e_{\sigma}}^{-1}\mathbf{M^e}\mathbf{s_e}
\\left(\mathbf{C} \mathbf{M^e_{\sigma}}^{-1} \mathbf{C}^T \mathbf{M_{\mu^{-1}}^f} + i \omega \\right)\mathbf{b} = \mathbf{s_m} + \mathbf{M^e_{\sigma}}^{-1}\mathbf{M^e}\mathbf{s_e}
.. note ::
The inverse problem will not work with full anisotropy
The inverse problem will not work with full anisotropy
"""
_fieldType = 'b'
@@ -345,14 +345,14 @@ class ProblemFDEM_b(BaseFDEMProblem):
return MeSigmaIDeriv.T * (C.T * v)
if self._makeASymmetric is True:
return MfMui.T * ( C * ( MeSigmaIDeriv * v ) )
return C * ( MeSigmaIDeriv * v )
return MfMui.T * ( C * ( MeSigmaIDeriv * v ) )
return C * ( MeSigmaIDeriv * v )
def getRHS(self, freq):
"""
.. math ::
\mathbf{RHS} = \mathbf{s_m} + \mathbf{M^e_{\sigma}}^{-1}\mathbf{s_e}
\mathbf{RHS} = \mathbf{s_m} + \mathbf{M^e_{\sigma}}^{-1}\mathbf{s_e}
:param float freq: Frequency
:rtype: numpy.ndarray (nE, nSrc)
@@ -404,7 +404,7 @@ class ProblemFDEM_b(BaseFDEMProblem):
SrcDeriv = C * (self.MeSigmaI * S_eDeriv)
elif adjoint:
SrcDeriv = self.MeSigmaI.T * ( C.T * S_eDeriv)
else:
else:
SrcDeriv = None
if RHSderiv is not None and SrcDeriv is not None:
@@ -412,7 +412,7 @@ class ProblemFDEM_b(BaseFDEMProblem):
elif SrcDeriv is not None:
RHSderiv = SrcDeriv
if RHSderiv is not None:
if RHSderiv is not None:
if self._makeASymmetric is True and not adjoint:
return MfMui.T * RHSderiv
@@ -427,13 +427,13 @@ class ProblemFDEM_b(BaseFDEMProblem):
class ProblemFDEM_j(BaseFDEMProblem):
"""
We eliminate \\\(\\\mathbf{h}\\\) using
We eliminate \\\(\\\mathbf{h}\\\) using
.. math ::
.. math ::
\mathbf{h} = \\frac{1}{i \omega} \mathbf{M_{\mu}^e}^{-1} \\left(-\mathbf{C}^T \mathbf{M_{\\rho}^f} \mathbf{j} + \mathbf{M^e} \mathbf{s_m} \\right)
and solve for \\\(\\\mathbf{j}\\\) using
and solve for \\\(\\\mathbf{j}\\\) using
.. math ::
@@ -510,7 +510,7 @@ class ProblemFDEM_j(BaseFDEMProblem):
S_m, S_e = self.getSourceTerm(freq)
C = self.mesh.edgeCurl
MeMuI = self.MeMuI
MeMuI = self.MeMuI
RHS = C * (MeMuI * S_m) - 1j * omega(freq) * S_e
if self._makeASymmetric is True:
@@ -521,7 +521,7 @@ class ProblemFDEM_j(BaseFDEMProblem):
def getRHSDeriv_m(self, src, v, adjoint=False):
C = self.mesh.edgeCurl
MeMuI = self.MeMuI
MeMuI = self.MeMuI
S_mDeriv, S_eDeriv = src.evalDeriv(self, adjoint)
if adjoint:
@@ -538,10 +538,10 @@ class ProblemFDEM_j(BaseFDEMProblem):
return - 1j * omega(freq) * S_eDerivv
else:
return None
else:
else:
S_mDerivv, S_eDerivv = S_mDeriv(v), S_eDeriv(v)
if S_mDerivv is not None and S_eDerivv is not None:
if S_mDerivv is not None and S_eDerivv is not None:
RHSDeriv = C * (MeMuI * S_mDerivv) - 1j * omega(freq) * S_eDerivv
elif S_mDerivv is not None:
RHSDeriv = C * (MeMuI * S_mDerivv)
@@ -560,17 +560,17 @@ class ProblemFDEM_j(BaseFDEMProblem):
class ProblemFDEM_h(BaseFDEMProblem):
"""
We eliminate \\\(\\\mathbf{j}\\\) using
We eliminate \\\(\\\mathbf{j}\\\) using
.. math ::
\mathbf{j} = \mathbf{C} \mathbf{h} - \mathbf{s_e}
and solve for \\\(\\\mathbf{h}\\\) using
and solve for \\\(\\\mathbf{h}\\\) using
.. math ::
\\left(\mathbf{C}^T \mathbf{M_{\\rho}^f} \mathbf{C} + i \omega \mathbf{M_{\mu}^e}\\right) \mathbf{h} = \mathbf{M^e} \mathbf{s_m} + \mathbf{C}^T \mathbf{M_{\\rho}^f} \mathbf{s_e}
\\left(\mathbf{C}^T \mathbf{M_{\\rho}^f} \mathbf{C} + i \omega \mathbf{M_{\mu}^e}\\right) \mathbf{h} = \mathbf{M^e} \mathbf{s_m} + \mathbf{C}^T \mathbf{M_{\\rho}^f} \mathbf{s_e}
"""
@@ -612,7 +612,7 @@ class ProblemFDEM_h(BaseFDEMProblem):
"""
.. math ::
\mathbf{RHS} = \mathbf{M^e} \mathbf{s_m} + \mathbf{C}^T \mathbf{M_{\\rho}^f} \mathbf{s_e}
\mathbf{RHS} = \mathbf{M^e} \mathbf{s_m} + \mathbf{C}^T \mathbf{M_{\\rho}^f} \mathbf{s_e}
:param float freq: Frequency
:rtype: numpy.ndarray (nE, nSrc)
@@ -647,12 +647,12 @@ class ProblemFDEM_h(BaseFDEMProblem):
S_eDeriv = S_eDeriv(v)
if S_mDeriv is not None:
if RHSDeriv is not None:
if RHSDeriv is not None:
RHSDeriv += S_mDeriv(v)
else:
else:
RHSDeriv = S_mDeriv(v)
if S_eDeriv is not None:
if RHSDeriv is not None:
if RHSDeriv is not None:
RHSDeriv += C.T * (MfRho * S_e)
else:
RHSDeriv = C.T * (MfRho * S_e)
SimPEG/EM/FDEM/FieldsFDEM.py
+22 -22
View File
@@ -1,5 +1,5 @@
from SimPEG import Survey, Problem, Utils, np, sp
from simpegEM.Utils.EMUtils import omega
from SimPEG import Survey, Problem, Utils, np, sp
from SimPEG.EM.Utils.EMUtils import omega
class FieldsFDEM(Problem.Fields):
@@ -30,11 +30,11 @@ class FieldsFDEM_e(FieldsFDEM):
for i, src in enumerate(srcList):
ep = src.ePrimary(self.prob)
if ep is not None:
ePrimary[:,i] = ep
ePrimary[:,i] = ep
return ePrimary
def _eSecondary(self, eSolution, srcList):
return eSolution
return eSolution
def _e(self, eSolution, srcList):
return self._ePrimary(eSolution,srcList) + self._eSecondary(eSolution,srcList)
@@ -54,7 +54,7 @@ class FieldsFDEM_e(FieldsFDEM):
bPrimary[:,i] += bp
return bPrimary
def _bSecondary(self, eSolution, srcList):
def _bSecondary(self, eSolution, srcList):
C = self._edgeCurl
b = (C * eSolution)
for i, src in enumerate(srcList):
@@ -123,7 +123,7 @@ class FieldsFDEM_b(FieldsFDEM):
return bSolution
def _b(self, bSolution, srcList):
return self._bPrimary(bSolution, srcList) + self._bSecondary(bSolution, srcList)
return self._bPrimary(bSolution, srcList) + self._bSecondary(bSolution, srcList)
def _bDeriv_u(self, src, v, adjoint=False):
return None
@@ -142,7 +142,7 @@ class FieldsFDEM_b(FieldsFDEM):
def _eSecondary(self, bSolution, srcList):
e = self._MeSigmaI * ( self._edgeCurl.T * ( self._MfMui * bSolution))
for i,src in enumerate(srcList):
for i,src in enumerate(srcList):
_,S_e = src.eval(self.prob)
if S_e is not None:
e[:,i] += -self._MeSigmaI * S_e
@@ -150,7 +150,7 @@ class FieldsFDEM_b(FieldsFDEM):
def _eSecondaryDeriv_u(self, src, v, adjoint=False):
if not adjoint:
return self._MeSigmaI * ( self._edgeCurl.T * ( self._MfMui * v) )
return self._MeSigmaI * ( self._edgeCurl.T * ( self._MfMui * v) )
else:
return self._MfMui.T * (self._edgeCurl * (self._MeSigmaI.T * v))
@@ -168,7 +168,7 @@ class FieldsFDEM_b(FieldsFDEM):
if not adjoint:
de_dm = self._MeSigmaIDeriv(w) * v
elif adjoint:
elif adjoint:
de_dm = self._MeSigmaIDeriv(w).T * v
_, S_eDeriv = src.evalDeriv(self.prob, adjoint)
@@ -205,7 +205,7 @@ class FieldsFDEM_j(FieldsFDEM):
FieldsFDEM.__init__(self,mesh,survey,**kwargs)
def startup(self):
self.prob = self.survey.prob
self.prob = self.survey.prob
self._edgeCurl = self.survey.prob.mesh.edgeCurl
self._MeMuI = self.survey.prob.MeMuI
self._MfRho = self.survey.prob.MfRho
@@ -215,7 +215,7 @@ class FieldsFDEM_j(FieldsFDEM):
def _jPrimary(self, jSolution, srcList):
jPrimary = np.zeros_like(jSolution,dtype = complex)
for i, src in enumerate(srcList):
jp = src.jPrimary(self.prob)
jp = src.jPrimary(self.prob)
if jp is not None:
jPrimary[:,i] += jp
return jPrimary
@@ -238,11 +238,11 @@ class FieldsFDEM_j(FieldsFDEM):
for i, src in enumerate(srcList):
hp = src.hPrimary(self.prob)
if hp is not None:
hPrimary[:,i] = hp
hPrimary[:,i] = hp
return hPrimary
def _hSecondary(self, jSolution, srcList):
h = self._MeMuI * (self._edgeCurl.T * (self._MfRho * jSolution) )
h = self._MeMuI * (self._edgeCurl.T * (self._MfRho * jSolution) )
for i, src in enumerate(srcList):
h[:,i] *= -1./(1j*omega(src.freq))
S_m,_ = src.eval(self.prob)
@@ -251,7 +251,7 @@ class FieldsFDEM_j(FieldsFDEM):
return h
def _hSecondaryDeriv_u(self, src, v, adjoint=False):
if not adjoint:
if not adjoint:
return -1./(1j*omega(src.freq)) * self._MeMuI * (self._edgeCurl.T * (self._MfRho * v) )
elif adjoint:
return -1./(1j*omega(src.freq)) * self._MfRho.T * (self._edgeCurl * ( self._MeMuI.T * v))
@@ -264,10 +264,10 @@ class FieldsFDEM_j(FieldsFDEM):
MfRhoDeriv = self._MfRhoDeriv
Me = self._Me
if not adjoint:
hDeriv_m = -1./(1j*omega(src.freq)) * MeMuI * (C.T * (MfRhoDeriv(jSolution)*v ) )
if not adjoint:
hDeriv_m = -1./(1j*omega(src.freq)) * MeMuI * (C.T * (MfRhoDeriv(jSolution)*v ) )
elif adjoint:
hDeriv_m = -1./(1j*omega(src.freq)) * MfRhoDeriv(jSolution).T * ( C * (MeMuI.T * v ) )
hDeriv_m = -1./(1j*omega(src.freq)) * MfRhoDeriv(jSolution).T * ( C * (MeMuI.T * v ) )
S_mDeriv,_ = src.evalDeriv(self.prob, adjoint)
@@ -282,14 +282,14 @@ class FieldsFDEM_j(FieldsFDEM):
return hDeriv_m
def _h(self, jSolution, srcList):
def _h(self, jSolution, srcList):
return self._hPrimary(jSolution, srcList) + self._hSecondary(jSolution, srcList)
def _hDeriv_u(self, src, v, adjoint=False):
return self._hSecondaryDeriv_u(src, v, adjoint)
def _hDeriv_m(self, src, v, adjoint=False):
# assuming the primary doesn't depend on the model
# assuming the primary doesn't depend on the model
return self._hSecondaryDeriv_m(src, v, adjoint)
@@ -339,7 +339,7 @@ class FieldsFDEM_h(FieldsFDEM):
for i, src in enumerate(srcList):
jp = src.jPrimary(self.prob)
if jp is not None:
jPrimary[:,i] = jp
jPrimary[:,i] = jp
return jPrimary
def _jSecondary(self, hSolution, srcList):
@@ -370,5 +370,5 @@ class FieldsFDEM_h(FieldsFDEM):
return self._jSecondaryDeriv_u(src,v,adjoint)
def _jDeriv_m(self, src, v, adjoint=False):
# assuming the primary does not depend on the model
return self._jSecondaryDeriv_m(src,v,adjoint)
# assuming the primary does not depend on the model
return self._jSecondaryDeriv_m(src,v,adjoint)
SimPEG/EM/FDEM/SurveyFDEM.py
+7 -7
View File
@@ -1,10 +1,10 @@
from SimPEG import Survey, Problem, Utils, np, sp
from simpegEM.Utils import SrcUtils
from simpegEM.Utils.EMUtils import omega, e_from_j, j_from_e, b_from_h, h_from_b
from SimPEG.EM.Utils import SrcUtils
from SimPEG.EM.Utils.EMUtils import omega, e_from_j, j_from_e, b_from_h, h_from_b
from scipy.constants import mu_0
####################################################
# Receivers
# Receivers
####################################################
class RxFDEM(Survey.BaseRx):
@@ -99,13 +99,13 @@ class SrcFDEM(Survey.BaseSrc):
def eval(self, prob):
S_m = self.S_m(prob)
S_e = self.S_e(prob)
return S_m, S_e
return S_m, S_e
def evalDeriv(self, prob, v, adjoint=False):
return lambda v: self.S_mDeriv(prob,v,adjoint), lambda v: self.S_eDeriv(prob,v,adjoint)
def bPrimary(self, prob):
return None
return None
def hPrimary(self, prob):
return None
@@ -225,7 +225,7 @@ class SrcFDEM_RawVec(SrcFDEM):
return prob.Me * self._S_e
return self._S_e
class SrcFDEM_MagDipole(SrcFDEM):
#TODO: right now, orientation doesn't actually do anything! The methods in SrcUtils should take care of that
@@ -275,7 +275,7 @@ class SrcFDEM_MagDipole(SrcFDEM):
def S_m(self, prob):
b_p = self.bPrimary(prob)
return -1j*omega(self.freq)*b_p
return -1j*omega(self.freq)*b_p
def S_e(self, prob):
if all(np.r_[self.mu] == np.r_[prob.curModel.mu]):
SimPEG/EM/TDEM/BaseTDEM.py
+2 -2
View File
@@ -1,10 +1,10 @@
from SimPEG import Solver, Problem
from SimPEG.Problem import BaseTimeProblem
from simpegEM.Utils import SrcUtils
from SimPEG.EM.Utils import SrcUtils
from scipy.constants import mu_0
from SimPEG.Utils import sdiag, mkvc
from SimPEG import Utils, Mesh
from simpegEM.Base import BaseEMProblem
from SimPEG.EM.Base import BaseEMProblem
import numpy as np
SimPEG/EM/TDEM/SurveyTDEM.py
+1 -1
View File
@@ -1,6 +1,6 @@
from SimPEG import Utils, Survey, np
from SimPEG.Survey import BaseSurvey
from simpegEM.Utils import SrcUtils
from SimPEG.EM.Utils import SrcUtils
from BaseTDEM import FieldsTDEM