.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/tutorials/gradient.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_gallery_tutorials_gradient.py: 04. Gradient of the misfit function =================================== A basic example how to use the :attr:`emg3d.simulations.Simulation.gradient` routine to compute the adjoint-state gradient of the misfit function. Here we just show its usage. For this example we use the survey and data as obtained in the example :ref:`sphx_glr_gallery_tutorials_simulation.py`. .. GENERATED FROM PYTHON SOURCE LINES 13-25 .. code-block:: Python import os import pooch import emg3d import numpy as np import matplotlib.pyplot as plt from matplotlib.colors import LogNorm, SymLogNorm # Adjust this path to a folder of your choice. data_path = os.path.join('..', 'download', '') .. GENERATED FROM PYTHON SOURCE LINES 27-32 Load survey and data -------------------- First we load the survey and accompanying data as obtained in the example :ref:`sphx_glr_gallery_tutorials_simulation.py`. .. GENERATED FROM PYTHON SOURCE LINES 32-46 .. code-block:: Python fname = 'GemPy-II-survey-A.h5' pooch.retrieve( 'https://raw.github.com/emsig/data/2021-05-21/emg3d/surveys/'+fname, '5f2ed0b959a4f80f5378a071e6f729c6b7446898be7689ddc9bbd100a8f5bce7', fname=fname, path=data_path, ) survey = emg3d.load(data_path + fname)['survey'] # Let's have a look survey .. rst-class:: sphx-glr-script-out .. code-block:: none Data loaded from «/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II-survey-A.h5» [emg3d v0.17.1.dev22+ge23c468 (format 1.0) on 2021-04-14T22:24:03.229285]. .. raw:: html

Survey «GemPy-II Survey A»

<xarray.Dataset> Size: 10kB
    Dimensions:    (src: 3, rec: 45, freq: 2)
    Coordinates:
      * src        (src) <U6 72B 'TxED-1' 'TxED-2' 'TxED-3'
      * rec        (rec) <U7 1kB 'RxEP-01' 'RxEP-02' ... 'RxEP-44' 'RxEP-45'
      * freq       (freq) <U3 24B 'f-1' 'f-2'
    Data variables:
        observed   (src, rec, freq) complex128 4kB (-5.965512862515066e-14-5.5221...
        synthetic  (src, rec, freq) complex128 4kB (-1.821164983623268e-14-5.4279...
    Attributes:
        noise_floor:     1e-15
        relative_error:  0.05


.. GENERATED FROM PYTHON SOURCE LINES 47-56 We can see that the survey consists of three sources, 45 receivers, and two frequencies. Create an initial model ----------------------- To create an initial model we load the true model, but set all subsurface resistivities to 1 Ohm.m. So we are left with a homogeneous three-layer model air-seawater-subsurface, which includes the topography of the seafloor. .. GENERATED FROM PYTHON SOURCE LINES 57-70 .. code-block:: Python # Load true model fname = "GemPy-II.h5" pooch.retrieve( 'https://raw.github.com/emsig/data/2021-05-21/emg3d/models/'+fname, 'ea8c23be80522d3ca8f36742c93758370df89188816f50cb4e1b2a6a3012d659', fname=fname, path=data_path, ) model = emg3d.load(data_path + fname)['model'] grid = model.grid .. rst-class:: sphx-glr-script-out .. code-block:: none Data loaded from «/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II.h5» [emg3d v1.0.0rc3.dev5+g0cd9e09 (format 1.0) on 2021-05-21T18:40:16.721968]. .. GENERATED FROM PYTHON SOURCE LINES 71-96 .. code-block:: Python # Overwrite all subsurface resistivity values with 1.0 res = model.property_x subsurface = (res > 0.5) & (res < 1000) res[subsurface] = 1.0 model.property_x = res # QC the initial model and the survey. grid.plot_3d_slicer(model.property_x, xslice=12000, yslice=7500, pcolor_opts={'norm': LogNorm(vmin=0.3, vmax=200)}) # Plot survey in figure above fig = plt.gcf() fig.suptitle('Initial resistivity model (Ohm.m)') axs = fig.get_children() rec_coords = survey.receiver_coordinates() src_coords = survey.source_coordinates() axs[1].plot(rec_coords[0], rec_coords[1], 'bv') axs[2].plot(rec_coords[0], rec_coords[2], 'bv') axs[3].plot(rec_coords[2], rec_coords[1], 'bv') axs[1].plot(src_coords[0], src_coords[1], 'r*') axs[2].plot(src_coords[0], src_coords[2], 'r*') axs[3].plot(src_coords[2], src_coords[1], 'r*') .. image-sg:: /gallery/tutorials/images/sphx_glr_gradient_001.png :alt: Initial resistivity model (Ohm.m) :srcset: /gallery/tutorials/images/sphx_glr_gradient_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 97-100 Options for automatic gridding ------------------------------ .. GENERATED FROM PYTHON SOURCE LINES 100-115 .. code-block:: Python gridding_opts = { 'center': (src_coords[0][1], src_coords[1][1], -2200), 'properties': [0.3, 10, 1, 0.3], 'domain': ( [rec_coords[0].min()-100, rec_coords[0].max()+100], [rec_coords[1].min()-100, rec_coords[1].max()+100], [-5500, -2000] ), 'min_width_limits': (100, 100, 50), 'stretching': (None, None, [1.05, 1.5]), 'center_on_edge': False, } .. GENERATED FROM PYTHON SOURCE LINES 116-118 Create the Simulation --------------------- .. GENERATED FROM PYTHON SOURCE LINES 118-134 .. code-block:: Python simulation = emg3d.simulations.Simulation( name="Initial Model", # A name for this simulation survey=survey, # Our survey instance model=model, # The model gridding='both', # Src- and freq-dependent grids max_workers=4, # How many parallel jobs # solver_opts=..., # Any parameter to pass to emg3d.solve gridding_opts=gridding_opts, receiver_interpolation='linear', # For proper adjoint-state gradient ) # Let's QC our Simulation instance simulation .. raw:: html

Simulation «Initial Model»



.. GENERATED FROM PYTHON SOURCE LINES 135-137 Compute Gradient ---------------- .. GENERATED FROM PYTHON SOURCE LINES 137-141 .. code-block:: Python grad = simulation.gradient .. rst-class:: sphx-glr-script-out .. code-block:: none Compute efields 0/6 [00:00] Compute efields █▋ 1/6 [00:21] Compute efields █████ 3/6 [00:21] Compute efields ████████▎ 5/6 [00:34] Compute efields ██████████ 6/6 [00:34] Back-propagate 0/6 [00:00] Back-propagate █▋ 1/6 [00:21] Back-propagate ████████▎ 5/6 [00:34] Back-propagate ██████████ 6/6 [00:34] .. GENERATED FROM PYTHON SOURCE LINES 142-144 QC Gradient ''''''''''' .. GENERATED FROM PYTHON SOURCE LINES 144-170 .. code-block:: Python # Set the gradient of air and water to NaN. # This will eventually move directly into emgd3 (active and inactive cells). grad[~subsurface] = np.nan # Plot the gradient grid.plot_3d_slicer( grad.ravel('F'), xslice=12000, yslice=7500, zslice=-4000, pcolor_opts={'cmap': 'RdBu_r', 'norm': SymLogNorm( linthresh=1e-2, base=10, vmin=-1e1, vmax=1e1)} ) # Add survey fig = plt.gcf() fig.suptitle('Gradient of the misfit function') axs = fig.get_children() axs[1].plot(rec_coords[0], rec_coords[1], 'bv') axs[2].plot(rec_coords[0], rec_coords[2], 'bv') axs[3].plot(rec_coords[2], rec_coords[1], 'bv') axs[1].plot(src_coords[0], src_coords[1], 'r*') axs[2].plot(src_coords[0], src_coords[2], 'r*') axs[3].plot(src_coords[2], src_coords[1], 'r*') .. image-sg:: /gallery/tutorials/images/sphx_glr_gradient_002.png :alt: Gradient of the misfit function :srcset: /gallery/tutorials/images/sphx_glr_gradient_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 171-173 .. code-block:: Python emg3d.Report() .. raw:: html
Thu Jul 04 15:35:14 2024 CEST
OS Linux (Ubuntu 22.04) CPU(s) 16 Machine x86_64
Architecture 64bit RAM 31.0 GiB Environment Python
File system ext4
Python 3.12.4 | packaged by conda-forge | (main, Jun 17 2024, 10:23:07) [GCC 12.3.0]
numpy 1.26.4 scipy 1.14.0 numba 0.60.0
emg3d 1.8.3 empymod 2.3.1 xarray 2024.6.1.dev31+g6c2d8c33
discretize 0.10.0 h5py 3.11.0 matplotlib 3.8.4
tqdm 4.66.4
Intel(R) oneAPI Math Kernel Library Version 2023.2-Product Build 20230613 for Intel(R) 64 architecture applications


.. rst-class:: sphx-glr-timing **Total running time of the script:** (1 minutes 12.759 seconds) **Estimated memory usage:** 493 MB .. _sphx_glr_download_gallery_tutorials_gradient.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: gradient.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: gradient.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_