From cb9a70bacc41900107d2fa7ea7f80e14f260f6ce Mon Sep 17 00:00:00 2001 From: Lindsey Heagy Date: Sun, 7 Feb 2016 12:53:57 -0800 Subject: [PATCH] docs for source --- SimPEG/EM/FDEM/SrcFDEM.py | 293 +++++++++++++++++++++++++++++++++++--- 1 file changed, 274 insertions(+), 19 deletions(-) diff --git a/SimPEG/EM/FDEM/SrcFDEM.py b/SimPEG/EM/FDEM/SrcFDEM.py index b29768ac..ba57d39b 100644 --- a/SimPEG/EM/FDEM/SrcFDEM.py +++ b/SimPEG/EM/FDEM/SrcFDEM.py @@ -1,55 +1,141 @@ from SimPEG import Survey, Problem, Utils, np, sp from scipy.constants import mu_0 from SimPEG.EM.Utils import * -from SimPEG.Utils import Zero -# from SurveyFDEM import Rx - +from SimPEG.Utils import Zero class BaseSrc(Survey.BaseSrc): + """ + Base source class for FDEM Survey + """ + freq = None - # rxPair = Rx + # rxPair = RxFDEM integrate = True def eval(self, prob): + """ + Evaluate the source terms. + - :math:`S_m` : magnetic source term + - :math:`S_e` : electric source term + + :param Problem prob: FDEM Problem + :rtype: (numpy.ndarray, numpy.ndarray) + :return: tuple with magnetic source term and electric source term + """ S_m = self.S_m(prob) S_e = self.S_e(prob) return S_m, S_e - def evalDeriv(self, prob, v, adjoint=False): - return lambda v: self.S_mDeriv(prob,v,adjoint), lambda v: self.S_eDeriv(prob,v,adjoint) + def evalDeriv(self, prob, v=None, adjoint=False): + """ + Derivatives of the source terms with respect to the inversion model + - :code:`S_mDeriv` : derivative of the magnetic source term + - :code:`S_eDeriv` : derivative of the electric source term + + :param Problem prob: FDEM Problem + :param numpy.ndarray v: vector to take product with + :param bool adjoint: adjoint? + :rtype: (numpy.ndarray, numpy.ndarray) + :return: tuple with magnetic source term and electric source term derivatives times a vector + """ + if v is not None: + return self.S_mDeriv(prob,v,adjoint), self.S_eDeriv(prob,v,adjoint) + else: + return lambda v: self.S_mDeriv(prob,v,adjoint), lambda v: self.S_eDeriv(prob,v,adjoint) def bPrimary(self, prob): + """ + Primary magnetic flux density + + :param Problem prob: FDEM Problem + :rtype: numpy.ndarray + :return: primary magnetic flux density + """ return Zero() def hPrimary(self, prob): + """ + Primary magnetic field + + :param Problem prob: FDEM Problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ return Zero() def ePrimary(self, prob): + """ + Primary electric field + + :param Problem prob: FDEM Problem + :rtype: numpy.ndarray + :return: primary electric field + """ return Zero() def jPrimary(self, prob): + """ + Primary current density + + :param Problem prob: FDEM Problem + :rtype: numpy.ndarray + :return: primary current density + """ return Zero() def S_m(self, prob): + """ + Magnetic source term + + :param Problem prob: FDEM Problem + :rtype: numpy.ndarray + :return: magnetic source term on mesh + """ return Zero() def S_e(self, prob): + """ + Electric source term + + :param Problem prob: FDEM Problem + :rtype: numpy.ndarray + :return: electric source term on mesh + """ return Zero() def S_mDeriv(self, prob, v, adjoint = False): + """ + Derivative of magnetic source term with respect to the inversion model + + :param Problem prob: FDEM Problem + :param numpy.ndarray v: vector to take product with + :param bool adjoint: adjoint? + :rtype: numpy.ndarray + :return: product of magnetic source term derivative with a vector + """ + return Zero() def S_eDeriv(self, prob, v, adjoint = False): + """ + Derivative of electric source term with respect to the inversion model + + :param Problem prob: FDEM Problem + :param numpy.ndarray v: vector to take product with + :param bool adjoint: adjoint? + :rtype: numpy.ndarray + :return: product of electric source term derivative with a vector + """ return Zero() class RawVec_e(BaseSrc): """ - RawVec electric source. It is defined by the user provided vector S_e + RawVec electric source. It is defined by the user provided vector S_e - :param numpy.array S_e: electric source term - :param float freq: frequency - :param rxList: receiver list + :param list rxList: receiver list + :param float freq: frequency + :param numpy.array S_e: electric source term """ def __init__(self, rxList, freq, S_e): #, ePrimary=None, bPrimary=None, hPrimary=None, jPrimary=None): @@ -58,16 +144,17 @@ class RawVec_e(BaseSrc): BaseSrc.__init__(self, rxList) def S_e(self, prob): + return self._S_e class RawVec_m(BaseSrc): """ - RawVec magnetic source. It is defined by the user provided vector S_m + RawVec magnetic source. It is defined by the user provided vector S_m - :param numpy.array S_m: magnetic source term - :param float freq: frequency - :param rxList: receiver list + :param float freq: frequency + :param rxList: receiver list + :param numpy.array S_m: magnetic source term """ def __init__(self, rxList, freq, S_m, integrate = True): #ePrimary=Zero(), bPrimary=Zero(), hPrimary=Zero(), jPrimary=Zero()): @@ -78,17 +165,24 @@ class RawVec_m(BaseSrc): BaseSrc.__init__(self, rxList) def S_m(self, prob): + """ + Magnetic source term + + :param Problem prob: FDEM Problem + :rtype: numpy.ndarray + :return: magnetic source term on mesh + """ return self._S_m class RawVec(BaseSrc): """ - RawVec source. It is defined by the user provided vectors S_m, S_e + RawVec source. It is defined by the user provided vectors S_m, S_e - :param numpy.array S_m: magnetic source term - :param numpy.array S_e: electric source term - :param float freq: frequency - :param rxList: receiver list + :param rxList: receiver list + :param float freq: frequency + :param numpy.array S_m: magnetic source term + :param numpy.array S_e: electric source term """ def __init__(self, rxList, freq, S_m, S_e, integrate = True): self._S_m = np.array(S_m,dtype=complex) @@ -109,6 +203,51 @@ class RawVec(BaseSrc): class MagDipole(BaseSrc): + """ + Point magnetic dipole source calculated by taking the curl of a magnetic + vector potential. By taking the discrete curl, we ensure that the magnetic + flux density is divergence free (no magnetic monopoles!). + + This approach uses a primary-secondary in frequency. Here we show the + derivation for E-B formulation noting that similar steps are followed for + the H-J formulation. + + .. math:: + \mathbf{C} \mathbf{e} + i \omega \mathbf{b} = \mathbf{s_m} \\\\ + {\mathbf{C}^T \mathbf{M_{\mu^{-1}}^f} \mathbf{b} - \mathbf{M_{\sigma}^e} \mathbf{e} = \mathbf{s_e}} + + We split up the fields and :math:`\mu^{-1}` into primary (:math:`\mathbf{P}`) and secondary (:math:`\mathbf{S}`) components + + - :math:`\mathbf{e} = \mathbf{e^P} + \mathbf{e^S}` + - :math:`\mathbf{b} = \mathbf{b^P} + \mathbf{b^S}` + - :math:`\\boldsymbol{\mu}^{\mathbf{-1}} = \\boldsymbol{\mu}^{\mathbf{-1}^\mathbf{P}} + \\boldsymbol{\mu}^{\mathbf{-1}^\mathbf{S}}` + + and define a zero-frequency primary problem, noting that the source is + generated by a divergence free electric current + + .. math:: + \mathbf{C} \mathbf{e^P} = \mathbf{s_m^P} = 0 \\\\ + {\mathbf{C}^T \mathbf{{M_{\mu^{-1}}^f}^P} \mathbf{b^P} - \mathbf{M_{\sigma}^e} \mathbf{e^P} = \mathbf{M^e} \mathbf{s_e^P}} + + Since :math:`\mathbf{e^P}` is curl-free, divergence-free, we assume that there is no constant field background, the :math:`\mathbf{e^P} = 0`, so our primary problem is + + .. math:: + \mathbf{e^P} = 0 \\\\ + {\mathbf{C}^T \mathbf{M_{\mu^{-1}}^f^P} \mathbf{b^P} = \mathbf{s_e^P}} + + Our secondary problem is then + + .. math:: + \mathbf{C} \mathbf{e^S} + i \omega \mathbf{b^S} = - i \omega \mathbf{b^P} \\\\ + {\mathbf{C}^T \mathbf{M_{\mu^{-1}}^f} \mathbf{b^S} - \mathbf{M_{\sigma}^e} \mathbf{e^S} = -\mathbf{C}^T \mathbf{M_{\mu^{-1}^S}^f} \mathbf{b^P}} + + :param list rxList: receiver list + :param float freq: frequency + :param numpy.ndarray loc: source location (ie: :code:`np.r_[xloc,yloc,zloc]`) + :param string orientation: 'X', 'Y', 'Z' + :param float moment: magnetic dipole moment + :param float mu: background magnetic permeability + """ #TODO: right now, orientation doesn't actually do anything! The methods in SrcUtils should take care of that def __init__(self, rxList, freq, loc, orientation='Z', moment=1., mu = mu_0): @@ -121,6 +260,13 @@ class MagDipole(BaseSrc): BaseSrc.__init__(self, rxList) def bPrimary(self, prob): + """ + The primary magnetic flux density from a magnetic vector potential + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ eqLocs = prob._eqLocs if eqLocs is 'FE': @@ -152,14 +298,37 @@ class MagDipole(BaseSrc): return C*a def hPrimary(self, prob): + """ + The primary magnetic field from a magnetic vector potential + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ b = self.bPrimary(prob) return h_from_b(prob,b) def S_m(self, prob): + """ + The magnetic source term + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ + b_p = self.bPrimary(prob) return -1j*omega(self.freq)*b_p def S_e(self, prob): + """ + The electric source term + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ + if all(np.r_[self.mu] == np.r_[prob.curModel.mu]): return Zero() else: @@ -179,6 +348,21 @@ class MagDipole(BaseSrc): class MagDipole_Bfield(BaseSrc): + """ + Point magnetic dipole source calculated with the analytic solution for the + fields from a magnetic dipole. No discrete curl is taken, so the magnetic + flux density may not be strictly divergence free. + + This approach uses a primary-secondary in frequency in the same fashion as the MagDipole. + + :param list rxList: receiver list + :param float freq: frequency + :param numpy.ndarray loc: source location (ie: :code:`np.r_[xloc,yloc,zloc]`) + :param string orientation: 'X', 'Y', 'Z' + :param float moment: magnetic dipole moment + :param float mu: background magnetic permeability + """ + #TODO: right now, orientation doesn't actually do anything! The methods in SrcUtils should take care of that #TODO: neither does moment def __init__(self, rxList, freq, loc, orientation='Z', moment=1., mu = mu_0): @@ -190,6 +374,14 @@ class MagDipole_Bfield(BaseSrc): BaseSrc.__init__(self, rxList) def bPrimary(self, prob): + """ + The primary magnetic flux density from the analytic solution for magnetic fields from a dipole + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ + eqLocs = prob._eqLocs if eqLocs is 'FE': @@ -221,14 +413,35 @@ class MagDipole_Bfield(BaseSrc): return b def hPrimary(self, prob): + """ + The primary magnetic field from a magnetic vector potential + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ b = self.bPrimary(prob) return h_from_b(prob, b) def S_m(self, prob): + """ + The magnetic source term + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ b = self.bPrimary(prob) return -1j*omega(self.freq)*b def S_e(self, prob): + """ + The electric source term + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ if all(np.r_[self.mu] == np.r_[prob.curModel.mu]): return Zero() else: @@ -247,6 +460,20 @@ class MagDipole_Bfield(BaseSrc): class CircularLoop(BaseSrc): + """ + Circular loop magnetic source calculated by taking the curl of a magnetic + vector potential. By taking the discrete curl, we ensure that the magnetic + flux density is divergence free (no magnetic monopoles!). + + This approach uses a primary-secondary in frequency in the same fashion as the MagDipole. + + :param list rxList: receiver list + :param float freq: frequency + :param numpy.ndarray loc: source location (ie: :code:`np.r_[xloc,yloc,zloc]`) + :param string orientation: 'X', 'Y', 'Z' + :param float moment: magnetic dipole moment + :param float mu: background magnetic permeability + """ #TODO: right now, orientation doesn't actually do anything! The methods in SrcUtils should take care of that def __init__(self, rxList, freq, loc, orientation='Z', radius = 1., mu=mu_0): @@ -259,6 +486,13 @@ class CircularLoop(BaseSrc): BaseSrc.__init__(self, rxList) def bPrimary(self, prob): + """ + The primary magnetic flux density from a magnetic vector potential + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ eqLocs = prob._eqLocs if eqLocs is 'FE': @@ -289,14 +523,35 @@ class CircularLoop(BaseSrc): return C*a def hPrimary(self, prob): + """ + The primary magnetic field from a magnetic vector potential + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ b = self.bPrimary(prob) return 1./self.mu*b def S_m(self, prob): + """ + The magnetic source term + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ b = self.bPrimary(prob) return -1j*omega(self.freq)*b def S_e(self, prob): + """ + The electric source term + + :param Problem prob: FDEM problem + :rtype: numpy.ndarray + :return: primary magnetic field + """ if all(np.r_[self.mu] == np.r_[prob.curModel.mu]): return Zero() else: