<commit_before>"""
<NAME> al., 2017 1D FDEM and TDEM inversions
==============================================

Here, we perform a 1D inversion using both the frequency and time domain
codes. The forward simulations are conducted on a cylindrically symmetric
mesh. The source is a point magnetic dipole source.

This example is used in the paper

    <NAME>, <NAME>, <NAME>, <NAME>, <NAME>, A framework for simulation and inversion in electromagnetics,
    Computers & Geosciences, Volume 107, 2017, Pages 1-19, ISSN 0098-3004,
    http://dx.doi.org/10.1016/j.cageo.2017.06.018.

This example is on figshare:
https://doi.org/10.6084/m9.figshare.5035175

This example was updated for SimPEG 0.14.0 on January 31st, 2020 by <NAME>
"""
import discretize
from SimPEG import (
    maps, utils, data_misfit, regularization,
    optimization, inversion, inverse_problem, directives
)
import numpy as np
from SimPEG.electromagnetics import frequency_domain as FDEM, time_domain as TDEM, mu_0
import matplotlib.pyplot as plt
import matplotlib
try:
    from pymatsolver import Pardiso as Solver
except ImportError:
    from SimPEG import SolverLU as Solver


def run(plotIt=True, saveFig=False):

    # Set up cylindrically symmeric mesh
    cs, ncx, ncz, npad = 10., 15, 25, 13  # padded cyl mesh
    hx = [(cs, ncx), (cs, npad, 1.3)]
    hz = [(cs, npad, -1.3), (cs, ncz), (cs, npad, 1.3)]
    mesh = discretize.CylMesh([hx, 1, hz], '00C')

    # Conductivity model
    layerz = np.r_[-200., -100.]
    layer = (mesh.vectorCCz >= layerz[0]) & (mesh.vectorCCz <= layerz[1])
    active = mesh.vectorCCz < 0.
    sig_half = 1e-2  # Half-space conductivity
    sig_air = 1e-8  # Air conductivity
    sig_layer = 5e-2  # Layer conductivity
    sigma = np.ones(mesh.nCz)*sig_air
    sigma[active] = sig_half
    sigma[layer] = sig_layer

    # Mapping
    actMap = maps.InjectActiveCells(mesh, active, np.log(1e-8), nC=mesh.nCz)
    mapping = maps.ExpMap(mesh) * maps.SurjectVertical1D(mesh) * actMap
    mtrue = np.log(sigma[active])

    # ----- FDEM problem & survey ----- #
    rxlocs = utils.ndgrid([np.r_[50.], np.r_[0], np.r_[0.]])
    bzr = FDEM.Rx.PointMagneticFluxDensitySecondary(rxlocs, 'z', 'real')
    bzi = FDEM.Rx.PointMagneticFluxDensitySecondary(rxlocs, 'z', 'imag')

    freqs = np.logspace(2, 3, 5)
    srcLoc = np.array([0., 0., 0.])

    print('min skin depth = ', 500./np.sqrt(freqs.max() * sig_half),
          'max skin depth = ', 500./np.sqrt(freqs.min() * sig_half))
    print('max x ', mesh.vectorCCx.max(), 'min z ', mesh.vectorCCz.min(),
          'max z ', mesh.vectorCCz.max())

    srcList = [
        FDEM.Src.MagDipole([bzr, bzi], freq, srcLoc, orientation='Z')
        for freq in freqs
    ]

    surveyFD = FDEM.Survey(srcList)
    prbFD = FDEM.Simulation3DMagneticFluxDensity(
        mesh, survey=surveyFD, sigmaMap=mapping, solver=Solver)
    std = 0.03
    dataFD = prbFD.make_synthetic_data(mtrue, standard_deviation=std, add_noise=True)
    dataFD.noise_floor = np.linalg.norm(dataFD.dclean)*1e-5

    # FDEM inversion
    np.random.seed(1)
    dmisfit = data_misfit.l2_DataMisfit(simulation=prbFD, data=dataFD)
    regMesh = discretize.TensorMesh([mesh.hz[mapping.maps[-1].indActive]])
    reg = regularization.Simple(regMesh)
    opt = optimization.InexactGaussNewton(maxIterCG=10)
    invProb = inverse_problem.BaseInvProblem(dmisfit, reg, opt)<commit_msg>l2_DataMisft to new style L2DataMisft<commit_after>"""
<NAME> al., 2017 1D FDEM and TDEM inversions
==============================================

Here, we perform a 1D inversion using both the frequency and time domain
codes. The forward simulations are conducted on a cylindrically symmetric
mesh. The source is a point magnetic dipole source.

This example is used in the paper

    <NAME>, <NAME>, <NAME>, <NAME>, <NAME>, A framework for simulation and inversion in electromagnetics,
    Computers & Geosciences, Volume 107, 2017, Pages 1-19, ISSN 0098-3004,
    http://dx.doi.org/10.1016/j.cageo.2017.06.018.

This example is on figshare:
https://doi.org/10.6084/m9.figshare.5035175

This example was updated for SimPEG 0.14.0 on January 31st, 2020 by <NAME>riotti
"""
import discretize
from SimPEG import (
    maps, utils, data_misfit, regularization,
    optimization, inversion, inverse_problem, directives
)
import numpy as np
from SimPEG.electromagnetics import frequency_domain as FDEM, time_domain as TDEM, mu_0
import matplotlib.pyplot as plt
import matplotlib
try:
    from pymatsolver import Pardiso as Solver
except ImportError:
    from SimPEG import SolverLU as Solver


def run(plotIt=True, saveFig=False):

    # Set up cylindrically symmeric mesh
    cs, ncx, ncz, npad = 10., 15, 25, 13  # padded cyl mesh
    hx = [(cs, ncx), (cs, npad, 1.3)]
    hz = [(cs, npad, -1.3), (cs, ncz), (cs, npad, 1.3)]
    mesh = discretize.CylMesh([hx, 1, hz], '00C')

    # Conductivity model
    layerz = np.r_[-200., -100.]
    layer = (mesh.vectorCCz >= layerz[0]) & (mesh.vectorCCz <= layerz[1])
    active = mesh.vectorCCz < 0.
    sig_half = 1e-2  # Half-space conductivity
    sig_air = 1e-8  # Air conductivity
    sig_layer = 5e-2  # Layer conductivity
    sigma = np.ones(mesh.nCz)*sig_air
    sigma[active] = sig_half
    sigma[layer] = sig_layer

    # Mapping
    actMap = maps.InjectActiveCells(mesh, active, np.log(1e-8), nC=mesh.nCz)
    mapping = maps.ExpMap(mesh) * maps.SurjectVertical1D(mesh) * actMap
    mtrue = np.log(sigma[active])

    # ----- FDEM problem & survey ----- #
    rxlocs = utils.ndgrid([np.r_[50.], np.r_[0], np.r_[0.]])
    bzr = FDEM.Rx.PointMagneticFluxDensitySecondary(rxlocs, 'z', 'real')
    bzi = FDEM.Rx.PointMagneticFluxDensitySecondary(rxlocs, 'z', 'imag')

    freqs = np.logspace(2, 3, 5)
    srcLoc = np.array([0., 0., 0.])

    print('min skin depth = ', 500./np.sqrt(freqs.max() * sig_half),
          'max skin depth = ', 500./np.sqrt(freqs.min() * sig_half))
    print('max x ', mesh.vectorCCx.max(), 'min z ', mesh.vectorCCz.min(),
          'max z ', mesh.vectorCCz.max())

    srcList = [
        FDEM.Src.MagDipole([bzr, bzi], freq, srcLoc, orientation='Z')
        for freq in freqs
    ]

    surveyFD = FDEM.Survey(srcList)
    prbFD = FDEM.Simulation3DMagneticFluxDensity(
        mesh, survey=surveyFD, sigmaMap=mapping, solver=Solver)
    std = 0.03
    dataFD = prbFD.make_synthetic_data(mtrue, standard_deviation=std, add_noise=True)
    dataFD.noise_floor = np.linalg.norm(dataFD.dclean)*1e-5

    # FDEM inversion
    np.random.seed(1)
    dmisfit = data_misfit.L2DataMisfit(simulation=prbFD, data=dataFD)
    regMesh = discretize.TensorMesh([mesh.hz[mapping.maps[-1].indActive]])
    reg = regularization.Simple(regMesh)
    opt = optimization.InexactGaussNewton(maxIterCG=10)
    invProb = inverse_problem.BaseInvProblem(dmisfit, reg, opt)
{'fraction_non_alphanumeric': 0.14220655207874247, 'fraction_numerical': 0.037314529161157635, 'mean_word_length': 3.5538461538461537, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 19, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '3137009', 'n_tokens_mistral': 2701, 'n_tokens_neox': 2433, 'n_words': 611}