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

Adjoint test and Jtvec for multiple transmitters

Browse Source
This commit is contained in:
rowanc1
2014-03-17 10:28:43 -07:00
parent 5c4ec4a1cd
commit 95e5a09523
2 changed files with 46 additions and 18 deletions
Download Patch File Download Diff File Expand all files Collapse all files
simpegEM/FDEM/FDEM.py
+28 -9
View File
@@ -97,22 +97,20 @@ class ProblemFDEM_e(Problem.BaseProblem):
u = self.fields(m)
Jv = self.dataPair(self.survey)
dsig_dm = self.model.transformDeriv(m)
for i, freq in enumerate(self.survey.freqs):
e = u[freq, 'e']
A = self.getA(freq)
solver = Solver(A, options=self.solveOpts)
for tx in self.survey.getTransmitters(freq):
dMe_dsig = self.mesh.getEdgeInnerProductDeriv(m, v=e)
dsig_dm = self.model.transformDeriv(m)
for txi, tx in enumerate(self.survey.getTransmitters(freq)):
dMe_dsig = self.mesh.getEdgeInnerProductDeriv(m, v=e[:,txi])
P = tx.projectFieldsDeriv(self.mesh, u)
b = 1j*omega(freq) * ( dMe_dsig * ( dsig_dm * v ) )
Ab = solver.solve(b)
P = tx.projectFieldsDeriv(self.mesh, u)
Jv[tx] = -P*Ab
Ainvb = solver.solve(b)
Jv[tx] = -P*Ainvb
return Utils.mkvc(Jv)
@@ -120,7 +118,28 @@ class ProblemFDEM_e(Problem.BaseProblem):
def Jtvec(self, m, v, u=None):
if u is None:
u = self.fields(m)
raise NotImplementedError('Jtvec todo!')
# Ensure v is a data object.
if not isinstance(v, self.dataPair):
v = self.dataPair(self.survey, v)
Jtv = np.zeros(self.model.nP, dtype=complex)
dsig_dm = self.model.transformDeriv(m)
for i, freq in enumerate(self.survey.freqs):
e = u[freq, 'e']
AT = self.getA(freq).T
solver = Solver(AT, options=self.solveOpts)
for txi, tx in enumerate(self.survey.getTransmitters(freq)):
dMe_dsig = self.mesh.getEdgeInnerProductDeriv(m, v=e[:,txi])
P = tx.projectFieldsDeriv(self.mesh, u)
w = solver.solve(P.T * v[tx])
Jtv += - 1j*omega(freq) * ( dsig_dm.T * ( dMe_dsig.T * w ) )
return Jtv
simpegEM/Tests/test_FDEM.py
+18 -9
View File
@@ -2,6 +2,8 @@ import unittest
from SimPEG import *
import simpegEM as EM
TOL = 1e-10
class FDEM_bDerivTests(unittest.TestCase):
def setUp(self):
@@ -26,7 +28,12 @@ class FDEM_bDerivTests(unittest.TestCase):
rxList = EM.FDEM.RxListFDEM(XYZ, 'Ex')
Tx0 = EM.FDEM.TxFDEM(None, 'VMD', 3, rxList)
survey = EM.FDEM.SurveyFDEM([Tx0])
x = np.linspace(5,10,3)
XYZ = Utils.ndgrid(x,x,np.r_[0])
rxList = EM.FDEM.RxListFDEM(XYZ, 'Ey')
Tx1 = EM.FDEM.TxFDEM(None, 'VMD', 3, rxList)
survey = EM.FDEM.SurveyFDEM([Tx0, Tx1])
prb = EM.FDEM.ProblemFDEM_e(model)
prb.pair(survey)
@@ -37,7 +44,7 @@ class FDEM_bDerivTests(unittest.TestCase):
j_sx[4,4,4] = 1
j_s = np.r_[Utils.mkvc(j_sx),np.zeros(mesh.nEy+mesh.nEz)]
prb.j_s = j_s
prb.j_s = np.c_[j_s, j_s]
f = prb.fields(sigma)
self.sigma = sigma
@@ -51,13 +58,15 @@ class FDEM_bDerivTests(unittest.TestCase):
passed = Tests.checkDerivative(fun, x0, num=3, plotIt=False)
self.assertTrue(passed)
# def test_Jtvec(self):
# v = self.survey.nD
# x0 = self.sigma
# def fun(x):
# return self.survey.dpred(x), lambda x: self.prb.Jvec(x0, x)
# passed = Tests.checkDerivative(fun, x0, num=3, plotIt=False)
# self.assertTrue(passed)
def test_Jtvec_adjointTest(self):
v = np.random.rand(self.survey.nD)
w = np.random.rand(self.prb.model.nP)
m = self.sigma
u = self.prb.fields(m)
vJw = v.dot(self.prb.Jvec(m, w, u=u))
wJtv = w.dot(self.prb.Jtvec(m, v, u=u))
self.assertTrue(vJw - wJtv < TOL)
if __name__ == '__main__':