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

ProblemFDEM_x.getSource now returns two matrices, S_m and S_e, which cleans up getRHS and the Fields

Browse Source
This commit is contained in:
Lindsey
2015-04-17 10:53:25 -07:00
parent e980477031
commit 93778d13f1
4 changed files with 60 additions and 98 deletions
Download Patch File Download Diff File Expand all files Collapse all files
simpegEM/FDEM/FDEM.py
+23 -50
View File
@@ -109,13 +109,21 @@ class BaseFDEMProblem(BaseEMProblem):
:return: RHS
"""
Txs = self.survey.getTransmitters(freq)
j_m = range(len(Txs))
j_e = range(len(Txs))
for i, tx in enumerate(Txs):
j_m[i], j_e[i] = tx.getSource(self)
if self._eqLocs is 'FE':
S_m = 1j*np.zeros((self.mesh.nF,len(Txs)))
S_e = 1j*np.zeros((self.mesh.nE,len(Txs)))
elif self._eqLocs is 'EF':
S_m = 1j*np.zeros((self.mesh.nE,len(Txs)))
S_e = 1j*np.zeros((self.mesh.nF,len(Txs)))
return j_m, j_e
# return np.concatenate(rhs).reshape((-1, len(Txs)), order='F') #, np.concatenate(j_e).reshape((-1, len(Txs)), order='F')
for i, tx in enumerate(Txs):
smi, sei = tx.getSource(self)
if smi is not None:
S_m[:,i] = smi
if sei is not None:
S_e[:,i] = sei
return S_m, S_e
def getSourceDeriv(self,freq,adjoint=False):
return None, None
@@ -181,20 +189,11 @@ class ProblemFDEM_e(BaseFDEMProblem):
:return: RHS
"""
j_m, j_g = self.getSource(freq)
nTx_freq = self.survey.nTxByFreq[freq]
RHS = 1j*np.zeros([self.mesh.nE, nTx_freq])
S_m, S_e = self.getSource(freq)
C = self.mesh.edgeCurl
MfMui = self.MfMui
for ii in range(nTx_freq):
if j_m[ii] is not None:
RHS[:, ii] += C.T * (MfMui * j_m[ii])
if j_g[ii] is not None:
RHS[:, ii] += -1j*omega(freq)*j_g[ii]
RHS = C.T * (MfMui * S_m) -1j*omega(freq)*S_e
return RHS
@@ -256,20 +255,11 @@ class ProblemFDEM_b(BaseFDEMProblem):
:return: RHS
"""
j_m, j_g = self.getSource(freq)
nTx_freq = self.survey.nTxByFreq[freq]
RHS = 1j*np.zeros([self.mesh.nF, nTx_freq])
S_m, S_e = self.getSource(freq)
C = self.mesh.edgeCurl
MfSigmai = self.MfSigmai
for ii in range(nTx_freq):
if j_m[ii] is not None:
RHS[:,ii] += j_m[ii]
if j_g[ii] is not None:
RHS[:,ii] += C * ( MfSigmai * j_g[ii] )
MeSigmaI = self.MeSigmaI
RHS = S_m + C * ( MeSigmaI * S_e )
if self._makeASymmetric is True:
mui = self.MfMui
@@ -373,21 +363,12 @@ class ProblemFDEM_j(BaseFDEMProblem):
:return: RHS
"""
j_m, j_g = self.getSource(freq)
nTx_freq = self.survey.nTxByFreq[freq]
RHS = 1j*np.zeros([self.mesh.nF, nTx_freq])
S_m, S_e = self.getSource(freq)
C = self.mesh.edgeCurl
MeMuI = self.MeMuI
for ii in range(nTx_freq):
if j_m[ii] is not None:
RHS[:,ii] += C * (MeMuI * j_m[ii])
if j_g[ii] is not None:
RHS[:,ii] += -1j * omega(freq) * j_g[ii]
RHS = C * (MeMuI * S_m) - 1j * omega(freq) * S_e
if self._makeASymmetric is True:
MfSigi = self.MfSigmai
return MfSigi.T*RHS
@@ -465,19 +446,11 @@ class ProblemFDEM_h(BaseFDEMProblem):
:return: RHS
"""
j_m, j_g = self.getSource(freq)
nTx_freq = self.survey.nTxByFreq[freq]
RHS = 1j*np.zeros([self.mesh.nE, nTx_freq])
S_m, S_e = self.getSource(freq)
C = self.mesh.edgeCurl
MfSigmai = self.MfSigmai
for ii in range(nTx_freq):
if j_m[ii] is not None:
RHS[:,ii] += j_m[ii]
if j_g[ii] is not None:
RHS[:,ii] += C.T * ( MfSigmai * j_g[ii] )
RHS = S_m + C.T * ( MfSigmai * S_e )
return RHS
simpegEM/FDEM/FieldsFDEM.py
+34 -42
View File
@@ -4,9 +4,9 @@ from simpegEM.Utils.EMUtils import omega
class FieldsFDEM(Problem.Fields):
"""Fancy Field Storage for a FDEM survey."""
knownFields = None
knownFields = {}
dtype = complex
class FieldsFDEM_e(FieldsFDEM):
knownFields = {'e':'E'}
@@ -30,23 +30,21 @@ class FieldsFDEM_e(FieldsFDEM):
return None
def _b(self, e, tx): #adjoint=False
b = self._b_sec(e,tx)
j_m,_ = self.getSource(tx.freq)
if j_m[0] is not None:
b += 1./(1j*omega(tx.freq)) * np.array([j_m[0]]).T
return b
b_sec = self._b_sec(e,tx)
S_m,_ = self.getSource(tx.freq)
return b_sec + 1./(1j*omega(tx.freq)) * S_m
def _bDeriv(self, e, tx, adjoint=False):
j_mDeriv,_ = self.getSourceDeriv(tx.freq, adjoint)
S_mDeriv,_ = self.getSourceDeriv(tx.freq, adjoint)
b_secDeriv = self._b_secDeriv(e,tx.freq,adjoint)
if j_mDeriv is None & b_secDeriv is None:
if S_mDeriv is None & b_secDeriv is None:
return None
elif b_secDeriv is None:
return 1./(1j*omega(tx.freq)) * j_mDeriv
elif j_mDeriv is None:
return 1./(1j*omega(tx.freq)) * S_mDeriv
elif S_mDeriv is None:
return b_secDeriv
else:
return 1./(1j*omega(tx.freq)) * j_mDeriv + b_secDeriv
return 1./(1j*omega(tx.freq)) * S_mDeriv + b_secDeriv
class FieldsFDEM_b(FieldsFDEM):
@@ -73,24 +71,22 @@ class FieldsFDEM_b(FieldsFDEM):
return None
def _e(self, b, tx):
e = self._e_sec(b,tx)
_, j_g = self.getSource(tx.freq)
if j_g[0] is not None:
e += -np.array([j_g[0]]).T
return e
e_sec = self._e_sec(b,tx)
_, S_e = self.getSource(tx.freq)
return e_sec + S_e
def _eDeriv(self, b, tx, adjoint=False):
_,j_gDeriv = self.getSourceDeriv(tx.freq, adjoint)
_,S_eDeriv = self.getSourceDeriv(tx.freq, adjoint)
e_secDeriv = self._e_secDeriv(b,tx,adjoint)
if j_gDeriv is None & e_secDeriv is None:
if S_eDeriv is None & e_secDeriv is None:
return None
elif e_secDeriv is None:
return -j_gDeriv
elif j_gDeriv is None:
return -S_eDeriv
elif S_eDeriv is None:
return e_secDeriv
else:
return e_secDeriv - j_gDeriv
return e_secDeriv - S_eDeriv
class FieldsFDEM_j(FieldsFDEM):
@@ -129,23 +125,21 @@ class FieldsFDEM_j(FieldsFDEM):
raise NotImplementedError('fieldType "%s" is not implemented.' % fieldType)
def _h(self, j, tx): #adjoint=False
h = self._h_sec(j,tx)
j_m,_ = self.getSource(tx.freq)
if j_m[0] is not None:
h += 1./(1j*omega(tx.freq)) * self.MeMuI * np.array([j_m[0]]).T
return h
h_sec = self._h_sec(j,tx)
S_m,_ = self.getSource(tx.freq)
return h_sec + 1./(1j*omega(tx.freq)) * self.MeMuI * S_m
def _hDeriv(self, j, tx, adjoint=False):
j_mDeriv,_ = self.getSourceDeriv(tx.freq, adjoint)
S_mDeriv,_ = self.getSourceDeriv(tx.freq, adjoint)
h_secDeriv = self._h_secDeriv(j,tx.freq,adjoint)
if j_mDeriv is None & h_secDeriv is None:
if S_mDeriv is None & h_secDeriv is None:
return None
elif h_secDeriv is None:
return 1./(1j*omega(tx.freq)) * j_mDeriv
elif j_mDeriv is None:
return 1./(1j*omega(tx.freq)) * S_mDeriv
elif S_mDeriv is None:
return h_secDeriv
else:
return 1./(1j*omega(tx.freq)) * j_mDeriv + h_secDeriv
return 1./(1j*omega(tx.freq)) * S_mDeriv + h_secDeriv
class FieldsFDEM_h(FieldsFDEM):
knownFields = {'h':'E'}
@@ -171,23 +165,21 @@ class FieldsFDEM_h(FieldsFDEM):
return None
def _j(self, h, tx): #adjoint=False
j = self._j_sec(h,tx)
_,j_g = self.getSource(tx.freq)
if j_g[0] is not None:
j += -np.array([j_g[0]]).T
return j
j_sec = self._j_sec(h,tx)
_,S_e = self.getSource(tx.freq)
return j_sec - S_e
def _jDeriv(self, h, tx, adjoint=False):
_,j_gDeriv = self.getSourceDeriv(tx.freq, adjoint)
_,S_eDeriv = self.getSourceDeriv(tx.freq, adjoint)
j_secDeriv = self._j_secDeriv(j,tx.freq,adjoint)
if j_gDeriv is None & j_secDeriv is None:
if S_eDeriv is None & j_secDeriv is None:
return None
elif j_secDeriv is None:
return - j_gDeriv
elif j_gDeriv is None:
return - S_eDeriv
elif S_eDeriv is None:
return j_secDeriv
else:
return - j_gDeriv + j_secDeriv
return - S_eDeriv + j_secDeriv
# def calcFields(self, sol, freq, fieldType, adjoint=False):
simpegEM/FDEM/SurveyFDEM.py
-3
View File
@@ -182,10 +182,7 @@ class TxFDEM(Survey.BaseTx):
a = SRC
b_0 = C*a
# if solType == 'b' or solType == 'h':
return -1j*omega(freq)*b_0, None
# elif solType == 'e' or solType == 'j':
# return -1j*omega(freq)*C.T*mui*b_0, None
class SimpleTxFDEM_g(TxFDEM):
simpegEM/Tests/test_FDEM.py
+3 -3
View File
@@ -396,7 +396,7 @@ class FDEM_DerivTests(unittest.TestCase):
def test_EB_CrossCheck_byr_Eform(self):
self.assertTrue(crossCheckTest('e', 'byr'))
def test_EB_CrossCheck_bzr_Eform(self):
self.assertTrue(crossCheckTest('e', 'bzr')) # Doesn't make sense to test this for p-s approach
self.assertTrue(crossCheckTest('e', 'bzr'))
def test_EB_CrossCheck_bxi_Eform(self):
self.assertTrue(crossCheckTest('e', 'bxi'))
def test_EB_CrossCheck_byi_Eform(self):
@@ -422,8 +422,8 @@ class FDEM_DerivTests(unittest.TestCase):
self.assertTrue(crossCheckTest('j', 'hxr'))
def test_HJ_CrossCheck_hyr_Jform(self):
self.assertTrue(crossCheckTest('j', 'hyr'))
# def test_HJ_CrossCheck_hzr_Jform(self):
# self.assertTrue(crossCheckTest('j', 'hzr')) # Doesn't make sense to test this for p-s approach
def test_HJ_CrossCheck_hzr_Jform(self):
self.assertTrue(crossCheckTest('j', 'hzr'))
def test_HJ_CrossCheck_hxi_Jform(self):
self.assertTrue(crossCheckTest('j', 'hxi'))
def test_HJ_CrossCheck_hyi_Jform(self):