From 569d22bf6aab50fd16c08d57162231c78996153b Mon Sep 17 00:00:00 2001 From: Lindsey Date: Wed, 6 May 2015 16:58:46 -0700 Subject: [PATCH] made loc a kwarg --- simpegEM/FDEM/FDEM.py | 20 +++---- simpegEM/FDEM/FieldsFDEM.py | 57 +++++++++++--------- simpegEM/FDEM/SurveyFDEM.py | 101 +++++++++++++++++++++++++----------- simpegEM/Tests/test_FDEM.py | 4 +- 4 files changed, 115 insertions(+), 67 deletions(-) diff --git a/simpegEM/FDEM/FDEM.py b/simpegEM/FDEM/FDEM.py index 72dca5eb..0c4ddd84 100644 --- a/simpegEM/FDEM/FDEM.py +++ b/simpegEM/FDEM/FDEM.py @@ -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 diff --git a/simpegEM/FDEM/FieldsFDEM.py b/simpegEM/FDEM/FieldsFDEM.py index 69160ebe..ced02d08 100644 --- a/simpegEM/FDEM/FieldsFDEM.py +++ b/simpegEM/FDEM/FieldsFDEM.py @@ -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 diff --git a/simpegEM/FDEM/SurveyFDEM.py b/simpegEM/FDEM/SurveyFDEM.py index 7dbe8df0..5d60ef94 100644 --- a/simpegEM/FDEM/SurveyFDEM.py +++ b/simpegEM/FDEM/SurveyFDEM.py @@ -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] diff --git a/simpegEM/Tests/test_FDEM.py b/simpegEM/Tests/test_FDEM.py index 1c5f4230..763a7bfe 100644 --- a/simpegEM/Tests/test_FDEM.py +++ b/simpegEM/Tests/test_FDEM.py @@ -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])