.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "content/examples/06-tdem/plot_inv_tdem_1D.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_content_examples_06-tdem_plot_inv_tdem_1D.py: EM: TDEM: 1D: Inversion ======================= Here we will create and run a TDEM 1D inversion. .. GENERATED FROM PYTHON SOURCE LINES 7-97 .. image-sg:: /content/examples/06-tdem/images/sphx_glr_plot_inv_tdem_1D_001.png :alt: plot inv tdem 1D :srcset: /content/examples/06-tdem/images/sphx_glr_plot_inv_tdem_1D_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none SimPEG.InvProblem will set Regularization.reference_model to m0. SimPEG.InvProblem will set Regularization.reference_model to m0. SimPEG.InvProblem will set Regularization.reference_model to m0. SimPEG.InvProblem is setting bfgsH0 to the inverse of the eval2Deriv. ***Done using same Solver, and solver_opts as the Simulation3DElectricField problem*** /home/vsts/conda/envs/simpeg-test/lib/python3.8/site-packages/pymatsolver/direct.py:23: PardisoTypeConversionWarning: Converting csc_matrix matrix to CSR format. /home/vsts/conda/envs/simpeg-test/lib/python3.8/site-packages/pymatsolver/direct.py:73: PardisoTypeConversionWarning: Converting csc_matrix matrix to CSR format. model has any nan: 0 ============================ Inexact Gauss Newton ============================ # beta phi_d phi_m f |proj(x-g)-x| LS Comment ----------------------------------------------------------------------------- x0 has any nan: 0 0 9.59e+02 5.26e+03 0.00e+00 5.26e+03 6.56e+03 0 1 9.59e+02 3.94e+02 1.65e-01 5.52e+02 1.01e+03 0 2 1.92e+02 4.51e+01 3.00e-01 1.03e+02 2.71e+02 0 Skip BFGS 3 1.92e+02 3.92e+00 3.72e-01 7.52e+01 3.29e+01 0 Skip BFGS 4 3.83e+01 3.65e+00 3.66e-01 1.77e+01 6.33e+01 0 Skip BFGS 5 3.83e+01 4.31e-01 3.91e-01 1.54e+01 3.03e+00 0 ------------------------- STOP! ------------------------- 1 : |fc-fOld| = 2.2481e+00 <= tolF*(1+|f0|) = 5.2595e+02 1 : |xc-x_last| = 8.9019e-02 <= tolX*(1+|x0|) = 3.0149e+00 0 : |proj(x-g)-x| = 3.0311e+00 <= tolG = 1.0000e-01 0 : |proj(x-g)-x| = 3.0311e+00 <= 1e3*eps = 1.0000e-02 1 : maxIter = 5 <= iter = 5 ------------------------- DONE! ------------------------- | .. code-block:: Python import numpy as np from SimPEG.electromagnetics import time_domain from SimPEG import ( optimization, discretize, maps, data_misfit, regularization, inverse_problem, inversion, directives, utils, ) import matplotlib.pyplot as plt def run(plotIt=True): cs, ncx, ncz, npad = 5.0, 25, 15, 15 hx = [(cs, ncx), (cs, npad, 1.3)] hz = [(cs, npad, -1.3), (cs, ncz), (cs, npad, 1.3)] mesh = discretize.CylindricalMesh([hx, 1, hz], "00C") active = mesh.cell_centers_z < 0.0 layer = (mesh.cell_centers_z < 0.0) & (mesh.cell_centers_z >= -100.0) actMap = maps.InjectActiveCells(mesh, active, np.log(1e-8), nC=mesh.shape_cells[2]) mapping = maps.ExpMap(mesh) * maps.SurjectVertical1D(mesh) * actMap sig_half = 2e-3 sig_air = 1e-8 sig_layer = 1e-3 sigma = np.ones(mesh.shape_cells[2]) * sig_air sigma[active] = sig_half sigma[layer] = sig_layer mtrue = np.log(sigma[active]) rxOffset = 1e-3 rx = time_domain.Rx.PointMagneticFluxTimeDerivative( np.array([[rxOffset, 0.0, 30]]), np.logspace(-5, -3, 31), "z" ) src = time_domain.Src.MagDipole([rx], location=np.array([0.0, 0.0, 80])) survey = time_domain.Survey([src]) time_steps = [(1e-06, 20), (1e-05, 20), (0.0001, 20)] simulation = time_domain.Simulation3DElectricField( mesh, sigmaMap=mapping, survey=survey, time_steps=time_steps ) # d_true = simulation.dpred(mtrue) # create observed data rel_err = 0.05 data = simulation.make_synthetic_data(mtrue, relative_error=rel_err) dmisfit = data_misfit.L2DataMisfit(simulation=simulation, data=data) regMesh = discretize.TensorMesh([mesh.h[2][mapping.maps[-1].indActive]]) reg = regularization.WeightedLeastSquares(regMesh, alpha_s=1e-2, alpha_x=1.0) opt = optimization.InexactGaussNewton(maxIter=5, LSshorten=0.5) invProb = inverse_problem.BaseInvProblem(dmisfit, reg, opt) # Create an inversion object beta = directives.BetaSchedule(coolingFactor=5, coolingRate=2) betaest = directives.BetaEstimate_ByEig(beta0_ratio=1e0) inv = inversion.BaseInversion(invProb, directiveList=[beta, betaest]) m0 = np.log(np.ones(mtrue.size) * sig_half) simulation.counter = opt.counter = utils.Counter() opt.remember("xc") mopt = inv.run(m0) if plotIt: fig, ax = plt.subplots(1, 2, figsize=(10, 6)) ax[0].loglog(rx.times, -invProb.dpred, "b.-") ax[0].loglog(rx.times, -data.dobs, "r.-") ax[0].legend(("Noisefree", "$d^{obs}$"), fontsize=16) ax[0].set_xlabel("Time (s)", fontsize=14) ax[0].set_ylabel("$B_z$ (T)", fontsize=16) ax[0].set_xlabel("Time (s)", fontsize=14) ax[0].grid(color="k", alpha=0.5, linestyle="dashed", linewidth=0.5) plt.semilogx(sigma[active], mesh.cell_centers_z[active]) plt.semilogx(np.exp(mopt), mesh.cell_centers_z[active]) ax[1].set_ylim(-600, 0) ax[1].set_xlim(1e-4, 1e-2) ax[1].set_xlabel("Conductivity (S/m)", fontsize=14) ax[1].set_ylabel("Depth (m)", fontsize=14) ax[1].grid(color="k", alpha=0.5, linestyle="dashed", linewidth=0.5) plt.legend([r"$\sigma_{true}$", r"$\sigma_{pred}$"]) if __name__ == "__main__": run() plt.show() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 27.201 seconds) **Estimated memory usage:** 9 MB .. _sphx_glr_download_content_examples_06-tdem_plot_inv_tdem_1D.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_inv_tdem_1D.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_inv_tdem_1D.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_