mirror of
https://github.com/wassname/simpeg.git
synced 2026-07-08 19:12:57 +08:00
b from h and h from b
This commit is contained in:
@@ -34,6 +34,7 @@ class Fields_e(Fields):
|
||||
self._aveF2CCV = self.survey.prob.mesh.aveF2CCV
|
||||
self._sigma = self.survey.prob.curModel.sigma
|
||||
self._sigmaDeriv = self.survey.prob.curModel.sigmaDeriv
|
||||
self._mui = self.survey.prob.curModel.mui
|
||||
self._nC = self.survey.prob.mesh.nC
|
||||
|
||||
def _GLoc(self,fieldType):
|
||||
@@ -108,12 +109,25 @@ class Fields_e(Fields):
|
||||
sigma = self._sigma
|
||||
aveE2CCV = self._aveE2CCV
|
||||
n = int(aveE2CCV.shape[0] / self._nC) #TODO: This is a bit sloppy
|
||||
# Sigma = sdiag(np.kron(np.ones(n), sigma))
|
||||
Sigma = self.prob.MeSigma
|
||||
VI = sdiag(1./np.kron(np.ones(n), self.prob.mesh.vol))
|
||||
Sigma = sdiag(np.kron(np.ones(n), sigma))
|
||||
|
||||
e = self._e(eSolution, srcList)
|
||||
|
||||
return Sigma * (aveE2CCV * e)
|
||||
return VI * (aveE2CCV * (Sigma *e) )
|
||||
|
||||
def _h(self, eolution, srcList):
|
||||
b = self._b(eSolution, srcList)
|
||||
Mui = self.survey.prob.MfMui
|
||||
aveF2CCV = self._aveF2CCV
|
||||
n = int(aveF2CCV.shape[0] / self._nC) #TODO: This is a bit sloppy
|
||||
# Mui = sdiag(sp.kron(np.ones(n), mui))
|
||||
VI = sdiag(1./np.kron(np.ones(n), self.prob.mesh.vol))
|
||||
|
||||
|
||||
return VI * (aveF2CCV * (Mui * b))
|
||||
|
||||
|
||||
def _jDeriv_u(self, src, v, adjoint=False):
|
||||
raise NotImplementedError
|
||||
@@ -149,6 +163,8 @@ class Fields_b(Fields):
|
||||
'e' : ['bSolution','E','_e'],
|
||||
'ePrimary' : ['bSolution','E','_ePrimary'],
|
||||
'eSecondary' : ['bSolution','E','_eSecondary'],
|
||||
'j' : ['bSolution','C','_j'],
|
||||
'h' : ['bSolution','C','_h'],
|
||||
}
|
||||
|
||||
def __init__(self,mesh,survey,**kwargs):
|
||||
@@ -161,6 +177,13 @@ class Fields_b(Fields):
|
||||
self._MfMui = self.survey.prob.MfMui
|
||||
self._MeSigmaIDeriv = self.survey.prob.MeSigmaIDeriv
|
||||
self._Me = self.survey.prob.Me
|
||||
self._aveF2CCV = self.survey.prob.mesh.aveF2CCV
|
||||
self._aveE2CCV = self.survey.prob.mesh.aveE2CCV
|
||||
self._sigma = self.survey.prob.curModel.sigma
|
||||
self._mui = self.survey.prob.curModel.mui
|
||||
self._nC = self.survey.prob.mesh.nC
|
||||
|
||||
|
||||
|
||||
def _GLoc(self,fieldType):
|
||||
if fieldType == 'e':
|
||||
@@ -245,6 +268,29 @@ class Fields_b(Fields):
|
||||
# assuming primary doesn't depend on model
|
||||
return self._eSecondaryDeriv_m(src, v, adjoint)
|
||||
|
||||
def _j(self, bSolution, srcList):
|
||||
sigma = self._sigma
|
||||
aveE2CCV = self._aveE2CCV
|
||||
n = int(aveE2CCV.shape[0] / self._nC) #TODO: This is a bit sloppy
|
||||
# Sigma = sdiag(np.kron(np.ones(n), sigma))
|
||||
Sigma = self.prob.MeSigma
|
||||
VI = sdiag(1./np.kron(np.ones(n), self.prob.mesh.vol))
|
||||
|
||||
e = self._e(bSolution, srcList)
|
||||
|
||||
return VI * (aveE2CCV * (Sigma *e) )
|
||||
|
||||
def _h(self, bSolution, srcList):
|
||||
b = self._b(bSolution, srcList)
|
||||
Mui = self.survey.prob.MfMui
|
||||
aveF2CCV = self._aveF2CCV
|
||||
n = int(aveF2CCV.shape[0] / self._nC) #TODO: This is a bit sloppy
|
||||
# Mui = sdiag(sp.kron(np.ones(n), mui))
|
||||
VI = sdiag(1./np.kron(np.ones(n), self.prob.mesh.vol))
|
||||
|
||||
|
||||
return VI * (aveF2CCV * (Mui * b))
|
||||
|
||||
|
||||
class Fields_j(Fields):
|
||||
knownFields = {'jSolution':'F'}
|
||||
@@ -256,6 +302,7 @@ class Fields_j(Fields):
|
||||
'hPrimary' : ['jSolution','E','_hPrimary'],
|
||||
'hSecondary' : ['jSolution','E','_hSecondary'],
|
||||
'e' : ['jSolution','C','_e'],
|
||||
'b' : ['jSolution','C','_b'],
|
||||
}
|
||||
|
||||
def __init__(self,mesh,survey,**kwargs):
|
||||
@@ -269,7 +316,9 @@ class Fields_j(Fields):
|
||||
self._MfRhoDeriv = self.survey.prob.MfRhoDeriv
|
||||
self._Me = self.survey.prob.Me
|
||||
self._rho = self.survey.prob.curModel.rho
|
||||
self._mu = self.survey.prob.curModel.mui
|
||||
self._aveF2CCV = self.survey.prob.mesh.aveF2CCV
|
||||
self._aveE2CCV = self.survey.prob.mesh.aveE2CCV
|
||||
self._nC = self.survey.prob.mesh.nC
|
||||
|
||||
def _GLoc(self,fieldType):
|
||||
@@ -361,11 +410,14 @@ class Fields_j(Fields):
|
||||
rho = self._rho
|
||||
aveF2CCV = self._aveF2CCV
|
||||
n = int(aveF2CCV.shape[0] / self._nC) #TODO: This is a bit sloppy
|
||||
Rho = sdiag(np.kron(np.ones(n), rho))
|
||||
|
||||
Rho = self.prob.MfRho
|
||||
VI = sdiag(1./np.kron(np.ones(n), self.prob.mesh.vol))
|
||||
|
||||
|
||||
j = self._j(jSolution, srcList)
|
||||
|
||||
return 1./self.survey.prob.mesh.dim * Rho * (aveF2CCV * j)
|
||||
return VI * (aveF2CCV * (Rho * j))
|
||||
|
||||
def _eDeriv_u(self, src, v, adjoint=False):
|
||||
raise NotImplementedError
|
||||
@@ -373,6 +425,16 @@ class Fields_j(Fields):
|
||||
def _eDeriv_m(self, src, v, adjoint=False):
|
||||
raise NotImplementedError
|
||||
|
||||
def _b(self, jSolution, srcList):
|
||||
h = self._h(jSolution, srcList)
|
||||
Mu = self.prob.MeMu
|
||||
aveE2CCV = self._aveE2CCV
|
||||
n = int(aveE2CCV.shape[0] / self._nC) #TODO: This is a bit sloppy
|
||||
# Mu = sdiag(sp.kron(np.ones(n), mu))
|
||||
VI = sdiag(1./np.kron(np.ones(n), self.prob.mesh.vol))
|
||||
|
||||
return VI * (aveE2CCV * (Mu * h))
|
||||
|
||||
|
||||
class Fields_h(Fields):
|
||||
knownFields = {'hSolution':'E'}
|
||||
@@ -382,7 +444,9 @@ class Fields_h(Fields):
|
||||
'hSecondary' : ['hSolution','E','_hSecondary'],
|
||||
'j' : ['hSolution','F','_j'],
|
||||
'jPrimary' : ['hSolution','F','_jPrimary'],
|
||||
'jSecondary' : ['hSolution','F','_jSecondary']
|
||||
'jSecondary' : ['hSolution','F','_jSecondary'],
|
||||
'e' : ['hSolution','C','_e'],
|
||||
'b' : ['hSolution','C','_b'],
|
||||
}
|
||||
|
||||
def __init__(self,mesh,survey,**kwargs):
|
||||
@@ -393,6 +457,11 @@ class Fields_h(Fields):
|
||||
self._edgeCurl = self.survey.prob.mesh.edgeCurl
|
||||
self._MeMuI = self.survey.prob.MeMuI
|
||||
self._MfRho = self.survey.prob.MfRho
|
||||
self._rho = self.survey.prob.curModel.rho
|
||||
self._mu = self.survey.prob.curModel.mui
|
||||
self._aveF2CCV = self.survey.prob.mesh.aveF2CCV
|
||||
self._aveE2CCV = self.survey.prob.mesh.aveE2CCV
|
||||
self._nC = self.survey.prob.mesh.nC
|
||||
|
||||
def _GLoc(self,fieldType):
|
||||
if fieldType == 'h':
|
||||
@@ -458,3 +527,32 @@ class Fields_h(Fields):
|
||||
def _jDeriv_m(self, src, v, adjoint=False):
|
||||
# assuming the primary does not depend on the model
|
||||
return self._jSecondaryDeriv_m(src,v,adjoint)
|
||||
|
||||
def _e(self, hSolution, srcList):
|
||||
rho = self._rho
|
||||
aveF2CCV = self._aveF2CCV
|
||||
n = int(aveF2CCV.shape[0] / self._nC) #TODO: This is a bit sloppy
|
||||
|
||||
Rho = self.prob.MfRho
|
||||
VI = sdiag(1./np.kron(np.ones(n), self.prob.mesh.vol))
|
||||
|
||||
|
||||
j = self._j(hSolution, srcList)
|
||||
|
||||
return VI * (aveF2CCV * (Rho * j))
|
||||
|
||||
def _eDeriv_u(self, src, v, adjoint=False):
|
||||
raise NotImplementedError
|
||||
|
||||
def _eDeriv_m(self, src, v, adjoint=False):
|
||||
raise NotImplementedError
|
||||
|
||||
def _b(self, hSolution, srcList):
|
||||
h = self._h(hSolution, srcList)
|
||||
Mu = self.prob.MeMu
|
||||
aveE2CCV = self._aveE2CCV
|
||||
n = int(aveE2CCV.shape[0] / self._nC) #TODO: This is a bit sloppy
|
||||
# Mu = sdiag(sp.kron(np.ones(n), mu))
|
||||
VI = sdiag(1./np.kron(np.ones(n), self.prob.mesh.vol))
|
||||
|
||||
return VI * (aveE2CCV * (Mu * h))
|
||||
|
||||
@@ -32,15 +32,15 @@ def getFDEMProblem(fdemType, comp, SrcList, freq, verbose=False):
|
||||
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])] = 1.
|
||||
S_e[Utils.closestPoints(mesh,[0.,0.,0.],'Ez') + np.sum(mesh.vnE[:1])] = 1.
|
||||
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))
|
||||
|
||||
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])] = 1.
|
||||
S_e[Utils.closestPoints(mesh,[0.,0.,0.],'Fz') + np.sum(mesh.vnF[:1])] = 1.
|
||||
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))
|
||||
|
||||
if verbose:
|
||||
|
||||
@@ -8,6 +8,7 @@ from SimPEG.EM.Utils.testingUtils import getFDEMProblem
|
||||
testEB = True
|
||||
testHJ = True
|
||||
testEJ = True
|
||||
testBH = True
|
||||
verbose = False
|
||||
|
||||
TOLEBHJ = 1e-5
|
||||
@@ -145,5 +146,84 @@ class FDEM_CrossCheck(unittest.TestCase):
|
||||
def test_EJ_CrossCheck_ezi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'ezi', TOL=TOLEJHB))
|
||||
|
||||
def test_EJ_CrossCheck_bxr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'bxr', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_byr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'byr', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_bzr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'bzr', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_bxi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'bxi', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_byi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'byi', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_bzi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'bzi', TOL=TOLEJHB))
|
||||
|
||||
def test_EJ_CrossCheck_hxr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'hxr', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_hyr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'hyr', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_hzr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'hzr', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_hxi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'hxi', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_hyi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'hyi', TOL=TOLEJHB))
|
||||
def test_EJ_CrossCheck_hzi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('e', 'j', 'hzi', TOL=TOLEJHB))
|
||||
|
||||
if testBH:
|
||||
def test_BH_CrossCheck_jxr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'jxr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_jyr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'jyr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_jzr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'jzr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_jxi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'jxi', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_jyi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'jyi', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_jzi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'jzi', TOL=TOLEJHB))
|
||||
|
||||
def test_BH_CrossCheck_exr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'exr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_eyr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'eyr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_ezr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'ezr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_exi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'exi', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_eyi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'eyi', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_ezi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'ezi', TOL=TOLEJHB))
|
||||
|
||||
def test_BH_CrossCheck_bxr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'bxr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_byr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'byr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_bzr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'bzr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_bxi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'bxi', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_byi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'byi', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_bzi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'bzi', TOL=TOLEJHB))
|
||||
|
||||
def test_BH_CrossCheck_hxr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'hxr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_hyr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'hyr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_hzr_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'hzr', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_hxi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'hxi', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_hyi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'hyi', TOL=TOLEJHB))
|
||||
def test_BH_CrossCheck_hzi_Jform(self):
|
||||
self.assertTrue(crossCheckTest('b', 'h', 'hzi', TOL=TOLEJHB))
|
||||
|
||||
if __name__ == '__main__':
|
||||
unittest.main()
|
||||
Reference in New Issue
Block a user