wassname/simpeg
1
0
Fork 0
mirror of https://github.com/wassname/simpeg.git synced 2026-07-29 11:27:23 +08:00
Code Issues Packages Projects Releases Wiki Activity

made loc a kwarg

Browse Source
This commit is contained in:
Lindsey
2015-05-06 16:58:46 -07:00
parent 0869777284
commit 569d22bf6a
4 changed files with 115 additions and 67 deletions
Download Patch File Download Diff File Expand all files Collapse all files
simpegEM/FDEM/FDEM.py
+10 -10
View File
@@ -28,7 +28,7 @@ class BaseFDEMProblem(BaseEMProblem):
rhs = RHS(freq)
Ainv = self.Solver(A, **self.solverOpts)
sol = Ainv * rhs
Srcs = self.survey.getSource(freq)
Srcs = self.survey.getSrcByFreq(freq)
F[Srcs, self._fieldType] = sol
return F
@@ -102,13 +102,13 @@ class BaseFDEMProblem(BaseEMProblem):
return Jtv
def getSource(self, freq):
def getSourceTerm(self, freq):
"""
:param float freq: Frequency
:rtype: numpy.ndarray (nE or nF, nSrc)
:return: RHS
"""
Srcs = self.survey.getSource(freq)
Srcs = self.survey.getSrcByFreq(freq)
if self._eqLocs is 'FE':
S_m = np.zeros((self.mesh.nF,len(Srcs)), dtype=complex)
S_e = np.zeros((self.mesh.nE,len(Srcs)), dtype=complex)
@@ -117,7 +117,7 @@ class BaseFDEMProblem(BaseEMProblem):
S_e = np.zeros((self.mesh.nF,len(Srcs)), dtype=complex)
for i, src in enumerate(Srcs):
smi, sei = src.getSource(self)
smi, sei = src.eval(self)
if smi is not None:
S_m[:,i] = smi
if sei is not None:
@@ -125,8 +125,8 @@ class BaseFDEMProblem(BaseEMProblem):
return S_m, S_e
def getSourceDeriv(self,freq,m,v,u=None,adjoint=False):
raise NotImplementedError('getSourceDeriv not implemented yet')
def getSourceTermDeriv(self,freq,m,v,u=None,adjoint=False):
raise NotImplementedError('getSourceTermDeriv not implemented yet')
return None, None
@@ -190,7 +190,7 @@ class ProblemFDEM_e(BaseFDEMProblem):
:return: RHS
"""
S_m, S_e = self.getSource(freq)
S_m, S_e = self.getSourceTerm(freq)
C = self.mesh.edgeCurl
MfMui = self.MfMui
@@ -259,7 +259,7 @@ class ProblemFDEM_b(BaseFDEMProblem):
:return: RHS
"""
S_m, S_e = self.getSource(freq)
S_m, S_e = self.getSourceTerm(freq)
C = self.mesh.edgeCurl
MeSigmaI = self.MeSigmaI
@@ -371,7 +371,7 @@ class ProblemFDEM_j(BaseFDEMProblem):
:return: RHS
"""
S_m, S_e = self.getSource(freq)
S_m, S_e = self.getSourceTerm(freq)
C = self.mesh.edgeCurl
MeMuI = self.MeMuI
@@ -457,7 +457,7 @@ class ProblemFDEM_h(BaseFDEMProblem):
:return: RHS
"""
S_m, S_e = self.getSource(freq)
S_m, S_e = self.getSourceTerm(freq)
C = self.mesh.edgeCurl
MfSigmai = self.MfSigmai
simpegEM/FDEM/FieldsFDEM.py
+33 -24
View File
@@ -20,8 +20,8 @@ class FieldsFDEM_e(FieldsFDEM):
def startup(self):
self._edgeCurl = self.survey.prob.mesh.edgeCurl
self._getSource = self.survey.prob.getSource
self._getSourceDeriv = self.survey.prob.getSourceDeriv
# self._getSource = self.survey.prob.getSource
# self._getSourceDeriv = self.survey.prob.getSourceDeriv
def _b_sec(self, e, src): #adjoint=False
return - 1./(1j*omega(src.freq)) * (self._edgeCurl * e)
@@ -30,12 +30,15 @@ class FieldsFDEM_e(FieldsFDEM):
return None
def _b(self, e, src): #adjoint=False
b_sec = self._b_sec(e,src)
S_m,_ = self._getSource(src.freq)
return b_sec + 1./(1j*omega(src.freq)) * S_m
b = self._b_sec(e,src)
S_m = src._getS_m(self.survey.prob)
print S_m.shape
if S_m is not None:
b += 1./(1j*omega(src.freq)) * S_m
return b
def _bDeriv(self, e, src, v, adjoint=False):
S_mDeriv,_ = self._getSourceDeriv(src.freq, v, adjoint)
S_mDeriv,_ = src.getSourceDeriv(self.survey.prob, v, adjoint)
b_secDeriv = self._b_secDeriv(e, src.freq, v, adjoint)
if S_mDeriv is None & b_secDeriv is None:
return None
@@ -61,8 +64,8 @@ class FieldsFDEM_b(FieldsFDEM):
self._edgeCurl = self.survey.prob.mesh.edgeCurl
self._MeSigmaI = self.survey.prob.MeSigmaI
self._MfMui = self.survey.prob.MfMui
self._getSource = self.survey.prob.getSource
self._getSourceDeriv = self.survey.prob.getSourceDeriv
# self._getSource = self.survey.prob.getSource
# self._getSourceDeriv = self.survey.prob.getSourceDeriv
def _e_sec(self, b, src):
return self._MeSigmaI * ( self._edgeCurl.T * ( self._MfMui * b) )
@@ -71,12 +74,14 @@ class FieldsFDEM_b(FieldsFDEM):
return None
def _e(self, b, src):
e_sec = self._e_sec(b,src)
_, S_e = self._getSource(src.freq)
return e_sec + S_e
e = self._e_sec(b,src)
S_e = src._getS_e(self.survey.prob)
if S_e is not None:
e += S_e
return e
def _eDeriv(self, b, src, v, adjoint=False):
_,S_eDeriv = self._getSourceDeriv(src.freq, v, adjoint)
_,S_eDeriv = src.getSourceDeriv(self.survey.prob, v, adjoint)
e_secDeriv = self._e_secDeriv(b, src, v, adjoint)
if S_eDeriv is None & e_secDeriv is None:
@@ -103,8 +108,8 @@ class FieldsFDEM_j(FieldsFDEM):
self._edgeCurl = self.survey.prob.mesh.edgeCurl
self._MeMuI = self.survey.prob.MeMuI
self._MfSigmai = self.survey.prob.MfSigmai
self._getSource = self.survey.prob.getSource
self._getSourceDeriv = self.survey.prob.getSourceDeriv
# self._getSource = self.survey.prob.getSource
# self._getSourceDeriv = self.survey.prob.getSourceDeriv
self._curModel = self.survey.prob.curModel
def _h_sec(self, j, src): #v, adjoint=False
@@ -125,12 +130,14 @@ class FieldsFDEM_j(FieldsFDEM):
return -(1./(1j*omega(freq))) * dsig_dm.T * ( dsigi_dsig.T * ( dMf_dsigi.T * ( C * ( MeMuI.T * v ) ) ) )
def _h(self, j, src): #v, adjoint=False
h_sec = self._h_sec(j,src)
S_m,_ = self._getSource(src.freq)
return h_sec + 1./(1j*omega(src.freq)) * self._MeMuI * S_m
h = self._h_sec(j,src)
S_m = src._getS_m(self.survey.prob)
if S_m is not None:
h += 1./(1j*omega(src.freq)) * self._MeMuI * S_m
return h
def _hDeriv(self, j, src, v, adjoint=False):
S_mDeriv,_ = self._getSourceDeriv(src.freq, v, adjoint)
S_mDeriv,_ = src.getSourceDeriv(self.survey.prob, v, adjoint)
h_secDeriv = self._h_secDeriv(j,src.freq, v, adjoint)
if S_mDeriv is None & h_secDeriv is None:
return None
@@ -155,8 +162,8 @@ class FieldsFDEM_h(FieldsFDEM):
self._edgeCurl = self.survey.prob.mesh.edgeCurl
self._MeMuI = self.survey.prob.MeMuI
self._MfSigmai = self.survey.prob.MfSigmai
self._getSource = self.survey.prob.getSource
self._getSourceDeriv = self.survey.prob.getSourceDeriv
# self._getSource = self.survey.prob.getSource
# self._getSourceDeriv = self.survey.prob.getSourceDeriv
def _j_sec(self, h, src): # adjoint=False
return self._edgeCurl*h
@@ -165,12 +172,14 @@ class FieldsFDEM_h(FieldsFDEM):
return None
def _j(self, h, src): # adjoint=False
j_sec = self._j_sec(h,src)
_,S_e = self._getSource(src.freq)
return j_sec - S_e
j = self._j_sec(h,src)
S_e = src._getS_e(self.survey.prob)
if S_e is not None:
j += -S_e
return j
def _jDeriv(self, h, src, v, adjoint=False):
_,S_eDeriv = self._getSourceDeriv(src.freq, v, adjoint)
_,S_eDeriv = src.getSourceDeriv(self.survey.prob, v, adjoint)
j_secDeriv = self._j_secDeriv(j,src.freq, v, adjoint)
if S_eDeriv is None & j_secDeriv is None:
return None
simpegEM/FDEM/SurveyFDEM.py
+70 -31
View File
@@ -95,6 +95,12 @@ class SrcFDEM(Survey.BaseSrc):
freq = None
rxPair = RxFDEM
def eval(self, prob):
return self._getS_m(prob), self._getS_e(prob)
def evalDeriv(self, prob, v, adjoint=None):
return self._getS_mDeriv(prob,v,adjoint), self._getS_eDeriv(prob,v,adjoint)
class SrcFDEM_RawVec_e(SrcFDEM):
"""
@@ -105,16 +111,22 @@ class SrcFDEM_RawVec_e(SrcFDEM):
:param rxList: receiver list
"""
def __init__(self, S_e, freq, rxList):
def __init__(self, rxList, freq, S_e):
self.S_e = np.array(S_e,dtype=float)
self.freq = float(freq)
SrcFDEM.__init__(self, None, 'RawVec', rxList)
SrcFDEM.__init__(self, rxList)
def getSource(self, prob):
return None, self.S_e
def _getS_m(self, prob):
return None
def getSourceDeriv(self, prob, v, adjoint = False):
return None, None
def _getS_e(self, prob):
return self.S_e
def _getS_mDeriv(self, prob, v, adjoint = False):
return None
def _getS_eDeriv(self, prob, v, adjoint = False):
return None
class SrcFDEM_RawVec_m(SrcFDEM):
@@ -126,16 +138,22 @@ class SrcFDEM_RawVec_m(SrcFDEM):
:param rxList: receiver list
"""
def __init__(self, S_m, freq, rxList):
def __init__(self, rxList, freq, S_m):
self.S_m = np.array(S_m,dtype=float)
self.freq = float(freq)
SrcFDEM.__init__(self, None, 'RawVec', rxList)
SrcFDEM.__init__(self, rxList)
def getSource(self, prob):
return self.S_m, None
def _getS_m(self, prob):
return self.S_m
def getSourceDeriv(self, prob, v, adjoint = False):
return None, None
def _getS_e(self, prob):
return None
def _getS_mDeriv(self, prob, v, adjoint = False):
return None
def _getS_eDeriv(self, prob, v, adjoint = False):
return None
class SrcFDEM_RawVec(SrcFDEM):
@@ -147,30 +165,35 @@ class SrcFDEM_RawVec(SrcFDEM):
:param float freq: frequency
:param rxList: receiver list
"""
def __init__(self, S_m, S_e, freq, rxList):
def __init__(self, rxList, freq, S_m, S_e):
self.S_m = np.array(S_m,dtype=float)
self.S_e = np.array(S_e,dtype=float)
self.freq = float(freq)
SrcFDEM.__init__(self, None, 'RawVec', rxList)
SrcFDEM.__init__(self, rxList)
def getSource(self, prob):
return self.S_m, self.S_e
def _getS_m(self,prob):
return self.S_m
def getSourceDeriv(self, prob, v, adjoint=None):
return None, None
def _getS_e(self,prob):
return self.S_e
def _getS_mDeriv(self, prob, v, adjoint = False):
return None
def _getS_eDeriv(self, prob, v, adjoint = False):
return None
class SrcFDEM_MagDipole(SrcFDEM):
#TODO: right now, orientation doesn't actually do anything! The methods in SrcUtils should take care of that
def __init__(self, loc, freq, rxList, orientation='Z', moment=1.):
def __init__(self, rxList, freq, loc, orientation='Z', moment=1.):
self.freq = float(freq)
self.loc = loc
self.orientation = orientation
self.moment = moment
SrcFDEM.__init__(self, loc, 'MagDipole', rxList)
SrcFDEM.__init__(self, rxList)
def getSource(self, prob):
def _getS_m(self,prob):
eqLocs = prob._eqLocs
if eqLocs is 'FE':
@@ -201,8 +224,10 @@ class SrcFDEM_MagDipole(SrcFDEM):
S_m = -1j*omega(self.freq)*C*a
return S_m, None
return S_m
def _getS_e(self,prob):
return None
def getSourceDeriv(self, prob, v, adjoint=None):
return None, None
@@ -211,12 +236,13 @@ class SrcFDEM_MagDipole(SrcFDEM):
class SrcFDEM_MagDipole_Bfield(SrcFDEM):
#TODO: right now, orientation doesn't actually do anything! The methods in SrcUtils should take care of that
def __init__(self, loc, freq, rxList, orientation='Z'):
#TODO: neither does moment
def __init__(self, rxList, freq, loc, orientation='Z', moment=1.):
self.freq = float(freq)
self.orientation = orientation
SrcFDEM.__init__(self, loc, 'MagDipole', rxList)
SrcFDEM.__init__(self, rxList)
def getSource(self, prob):
def _getS_m(self,prob):
eqLocs = prob._eqLocs
if eqLocs is 'FE':
@@ -245,20 +271,33 @@ class SrcFDEM_MagDipole_Bfield(SrcFDEM):
bz = srcfct(self.loc, gridZ, 'z')
b = np.concatenate((bx,by,bz))
return -1j*omega(self.freq)*b, None
return -1j*omega(self.freq)*b
def getSourceDeriv(self, prob, v, adjoint=None):
return None, None
def _getS_e(self,prob):
return None
def _getS_mDeriv(self, prob, v, adjoint = False):
return None
def _getS_eDeriv(self, prob, v, adjoint = False):
return None
class SrcFDEM_CircularLoop(SrcFDEM):
#TODO: right now, orientation doesn't actually do anything! The methods in SrcUtils should take care of that
def __init__(self, loc, freq, rxList, orientation='Z', radius = 1.):
def __init__(self, rxList, freq, loc, orientation='Z', radius = 1.):
self.freq = float(freq)
self.orientation = orientation
self.radius = radius
SrcFDEM.__init__(self, loc, 'MagDipole', rxList)
SrcFDEM.__init__(self, rxList)
def _getS_mDeriv(self, prob, v, adjoint = False):
return None
def _getS_eDeriv(self, prob, v, adjoint = False):
return None
def getSource(self, prob):
eqLocs = prob._eqLocs
@@ -334,7 +373,7 @@ class SurveyFDEM(Survey.BaseSurvey):
self._nSrcByFreq[freq] = len(self.getSource(freq))
return self._nSrcByFreq
def getSource(self, freq):
def getSrcByFreq(self, freq):
"""Returns the sources associated with a specific frequency."""
assert freq in self._freqDict, "The requested frequency is not in this survey."
return self._freqDict[freq]
simpegEM/Tests/test_FDEM.py
+2 -2
View File
@@ -11,7 +11,7 @@ testAdjoint = False
testEB = True
testHJ = True
verbose = False
verbose = True
TOL = 1e-4
FLR = 1e-20 # "zero", so if residual below this --> pass regardless of order
@@ -35,7 +35,7 @@ def getProblem(fdemType, comp):
x = np.array([np.linspace(-30,-15,3),np.linspace(15,30,3)]) #don't sample right by the source
XYZ = Utils.ndgrid(x,x,np.r_[0.])
Rx0 = EM.FDEM.RxFDEM(XYZ, comp)
Src0 = EM.FDEM.SrcFDEM_MagDipole(np.r_[0.,0.,0.], freq, [Rx0])
Src0 = EM.FDEM.SrcFDEM_MagDipole([Rx0],freq=freq, loc=np.r_[0.,0.,0.])
survey = EM.FDEM.SurveyFDEM([Src0])