.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "content/tutorials/08-tdem/plot_fwd_1_em1dtm_dispersive.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_tutorials_08-tdem_plot_fwd_1_em1dtm_dispersive.py: 1D Forward Simulation with Chargeable and/or Magnetic Viscosity =============================================================== Here we use the module *simpeg.electromangetics.time_domain_1d* to compare predicted time domain data for a single sounding when the Earth is purely conductive, chargeable and/or magnetically viscous. In this tutorial, we focus on: - Defining receivers, sources, waveform and the survey - Defining physical properties when the Earth is chargeable and/or magnetically viscous - Setting physical property values as constant in the simulation Our survey geometry consists of a horizontal loop source with a radius of 10 m located 0.5 m above the Earth's surface. The receiver is located at the centre of the loop and measures the vertical component of the response. .. GENERATED FROM PYTHON SOURCE LINES 22-25 Import Modules -------------- .. GENERATED FROM PYTHON SOURCE LINES 25-35 .. code-block:: Python import numpy as np from matplotlib import pyplot as plt from simpeg import maps import simpeg.electromagnetics.time_domain as tdem from simpeg.electromagnetics.utils.em1d_utils import ColeCole, LogUniform # sphinx_gallery_thumbnail_number = 3 .. GENERATED FROM PYTHON SOURCE LINES 36-45 Create Survey ------------- Here we demonstrate a general way to define the receivers, sources, waveforms and survey. For this tutorial, the source is a horizontal loop whose current waveform is a unit step-off. Receivers are defined to measure the vertical component of the magnetic flux density and its time-derivative at the loop's center. .. GENERATED FROM PYTHON SOURCE LINES 45-86 .. code-block:: Python source_location = np.array([0.0, 0.0, 0.5]) source_orientation = "z" # "x", "y" or "z" current_amplitude = 1.0 # maximum amplitude of source current source_radius = 10.0 # loop radius receiver_location = np.array([0.0, 0.0, 0.5]) receiver_orientation = "z" # "x", "y" or "z" times = np.logspace(-6, -1, 51) # time channels (s) # Receiver list receiver_list = [] receiver_list.append( tdem.receivers.PointMagneticFluxDensity( receiver_location, times, orientation=receiver_orientation ) ) receiver_list.append( tdem.receivers.PointMagneticFluxTimeDerivative( receiver_location, times, orientation=receiver_orientation ) ) # Waveform waveform = tdem.sources.StepOffWaveform() # Sources source_list = [ tdem.sources.CircularLoop( receiver_list=receiver_list, location=source_location, waveform=waveform, current=current_amplitude, radius=source_radius, ) ] # Survey survey = tdem.Survey(source_list) .. GENERATED FROM PYTHON SOURCE LINES 87-100 Defining a 1D Layered Earth Model --------------------------------- Here, we define the layer thicknesses and physical properties for our 1D simulation. If we have N layers, parameters for the physical properties must be defined for each layer and we must provide N-1 layer thicknesses. The lowest layer is assumed to extend to infinity. For this tutorial, we predict the response for a halfspace model, however the script has been generalized to work for an arbitrary number of layers. If the Earth is a halfspace, the thicknesses could instead be defined by an empty array, and each physical property value by an array of length 1. .. GENERATED FROM PYTHON SOURCE LINES 100-175 .. code-block:: Python # Layer thicknesses thicknesses = np.array([40.0, 40.0]) n_layer = len(thicknesses) + 1 # In SimPEG, the Cole-Cole model is used to define a frequency-dependent # electrical conductivity when the Earth is chargeable. sigma = 1e-1 # infinite conductivity in S/m eta = 0.5 # intrinsice chargeability [0, 1] tau = 0.01 # central time-relaxation constant in seconds c = 0.75 # phase constant [0, 1] # In SimPEG, the a log-uniform distribution of time-relaxation constants is used # to define a frequency-dependent susceptibility when the Earth exhibits # magnetic viscosity chi = 0.001 # infinite susceptibility in SI dchi = 0.001 # amplitude of frequency-dependent susceptibility contribution tau1 = 1e-7 # lower limit for time relaxation constants in seconds tau2 = 1.0 # upper limit for time relaxation constants in seconds # For each physical property, the parameters must be defined for each layer. # In this case, we must define all parameters for the Cole-Cole conductivity # as well as the frequency-dependent magnetic susceptibility. sigma_model = sigma * np.ones(n_layer) eta_model = eta * np.ones(n_layer) tau_model = tau * np.ones(n_layer) c_model = c * np.ones(n_layer) chi_model = chi * np.ones(n_layer) dchi_model = dchi * np.ones(n_layer) tau1_model = tau1 * np.ones(n_layer) tau2_model = tau2 * np.ones(n_layer) # Here, we let the infinite conductivity be the model. As a result, we only # need to define the mapping for this parameter. All other parameters used # to define physical properties will be fixed when creating the simulation. model_mapping = maps.IdentityMap(nP=n_layer) # Compute and plot complex conductivity at all frequencies frequencies = np.logspace(-3, 6, 91) sigma_complex = ColeCole(frequencies, sigma, eta, tau, c) chi_complex = LogUniform(frequencies, chi, dchi, tau1, tau2) fig = plt.figure(figsize=(8, 6)) ax = fig.add_axes([0.15, 0.1, 0.8, 0.75]) ax.semilogx(frequencies, sigma * np.ones(len(frequencies)), "b", lw=3) ax.semilogx(frequencies, np.real(sigma_complex), "r", lw=3) ax.semilogx(frequencies, np.imag(sigma_complex), "r--", lw=3) ax.grid() ax.set_xlim(np.min(frequencies), np.max(frequencies)) ax.set_ylim(0.0, 1.1 * sigma) ax.set_xlabel("Frequency (Hz)") ax.set_ylabel("Conductivity") ax.set_title("Dispersive Electrical Conductivity") ax.legend( [r"$\sigma_{DC}$", r"$Re[\sigma (\omega)]$", r"$Im[\sigma (\omega)]$"], loc="center right", ) fig = plt.figure(figsize=(8, 6)) ax = fig.add_axes([0.15, 0.1, 0.8, 0.75]) ax.semilogx(frequencies, chi * np.ones(len(frequencies)), "b", lw=3) ax.semilogx(frequencies, np.real(chi_complex), "r", lw=3) ax.semilogx(frequencies, np.imag(chi_complex), "r--", lw=3) ax.grid() ax.set_xlim(np.min(frequencies), np.max(frequencies)) ax.set_ylim(-1.1 * chi, 1.1 * (chi + dchi)) ax.set_xlabel("Frequency (Hz)") ax.set_ylabel("Susceptibility") ax.set_title("Dispersive Magnetic Susceptibility") ax.legend( [r"$\chi_{DC}$", r"$Re[\chi (\omega)]$", r"$Im[\chi (\omega)]$"], loc="lower left" ) .. rst-class:: sphx-glr-horizontal * .. image-sg:: /content/tutorials/08-tdem/images/sphx_glr_plot_fwd_1_em1dtm_dispersive_001.png :alt: Dispersive Electrical Conductivity :srcset: /content/tutorials/08-tdem/images/sphx_glr_plot_fwd_1_em1dtm_dispersive_001.png :class: sphx-glr-multi-img * .. image-sg:: /content/tutorials/08-tdem/images/sphx_glr_plot_fwd_1_em1dtm_dispersive_002.png :alt: Dispersive Magnetic Susceptibility :srcset: /content/tutorials/08-tdem/images/sphx_glr_plot_fwd_1_em1dtm_dispersive_002.png :class: sphx-glr-multi-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 176-193 Define the Forward Simulation and Predict Data ---------------------------------------------- Here we predict the TDEM sounding for several halfspace models (conductive, chargeable, magnetically viscous). Since the physical properties defining the Earth are different, it requires a separate simulation object be created for each case. Each simulation requires the user define the survey, the layer thicknesses and a mapping. A universal mapping was created by letting sigma be the model. All other parameters used to define the physical properties are permanently set when defining the simulation. When using the *simpeg.electromagnetics.time_domain_1d* module, note that predicted data are organized by source, then by receiver, then by time channel. .. GENERATED FROM PYTHON SOURCE LINES 193-228 .. code-block:: Python # Simulate response for static conductivity simulation_conductive = tdem.Simulation1DLayered( survey=survey, thicknesses=thicknesses, sigmaMap=model_mapping ) dpred_conductive = simulation_conductive.dpred(sigma_model) # Simulate response for a chargeable Earth simulation_chargeable = tdem.Simulation1DLayered( survey=survey, thicknesses=thicknesses, sigmaMap=model_mapping, eta=eta, tau=tau, c=c, ) dpred_chargeable = simulation_chargeable.dpred(sigma_model) # Simulate response for viscous remanent magnetization mu0 = 4 * np.pi * 1e-7 mu = mu0 * (1 + chi) simulation_vrm = tdem.Simulation1DLayered( survey=survey, thicknesses=thicknesses, sigmaMap=model_mapping, mu=mu, dchi=dchi, tau1=tau1, tau2=tau2, ) dpred_vrm = simulation_vrm.dpred(sigma_model) .. rst-class:: sphx-glr-script-out .. code-block:: none /usr/share/miniconda/envs/simpeg-test/lib/python3.10/site-packages/numpy/_core/_asarray.py:126: ComplexWarning: Casting complex values to real discards the imaginary part .. GENERATED FROM PYTHON SOURCE LINES 229-232 Plotting Results ------------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 232-256 .. code-block:: Python fig = plt.figure(figsize=(12, 5)) ax1 = fig.add_axes([0.1, 0.1, 0.38, 0.85]) ax1.loglog(times, np.abs(dpred_conductive[0 : len(times)]), "k", lw=3) ax1.loglog(times, np.abs(dpred_chargeable[0 : len(times)]), "r", lw=3) ax1.loglog(times, np.abs(dpred_vrm[0 : len(times)]), "b", lw=3) ax1.set_xlim([times.min(), times.max()]) ax1.grid() ax1.legend(["Purely Inductive", "Chargeable", "Magnetically Viscous"]) ax1.set_xlabel("Times (s)") ax1.set_ylabel("|B| (T)") ax1.set_title("Magnetic Flux") ax2 = fig.add_axes([0.6, 0.1, 0.38, 0.85]) ax2.loglog(times, np.abs(dpred_conductive[len(times) :]), "k", lw=3) ax2.loglog(times, np.abs(dpred_chargeable[len(times) :]), "r", lw=3) ax2.loglog(times, np.abs(dpred_vrm[len(times) :]), "b", lw=3) ax2.set_xlim([times.min(), times.max()]) ax2.grid() ax2.legend(["Purely Inductive", "Chargeable", "Magnetically Viscous"]) ax2.set_xlabel("Times (s)") ax2.set_ylabel("|dB/dt| (T/s)") ax2.set_title("Time-Derivative of Magnetic Flux") .. image-sg:: /content/tutorials/08-tdem/images/sphx_glr_plot_fwd_1_em1dtm_dispersive_003.png :alt: Magnetic Flux, Time-Derivative of Magnetic Flux :srcset: /content/tutorials/08-tdem/images/sphx_glr_plot_fwd_1_em1dtm_dispersive_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Text(0.5, 1.0, 'Time-Derivative of Magnetic Flux') .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 4.008 seconds) **Estimated memory usage:** 288 MB .. _sphx_glr_download_content_tutorials_08-tdem_plot_fwd_1_em1dtm_dispersive.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_fwd_1_em1dtm_dispersive.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_fwd_1_em1dtm_dispersive.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_fwd_1_em1dtm_dispersive.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_