Thin plate inversion

Note

Go to the end to download the full example code.

Thin plate inversion#

Open In Colab

If you are running this notebook locally, make sure you’ve followed steps here to set up the environment. (This environment.yml file specifies a list of packages required to run the notebooks)

# -------------------------------------------------------- #
#                                                          #
#     Uncomment below to set up environment on "colab"     #
#                                                          #
# -------------------------------------------------------- #

# !pip install -U cofi
# !pip install git+https://github.com/JuergHauser/PyP223.git
# !git clone https://github.com/inlab-geo/cofi-examples.git
# %cd cofi-examples/examples/vtem_max

# If this notebook is run locally PyP223 needs to be installed separately by uncommenting the following line,
# that is by removing the # and the white space between it and the exclamation mark.
# !pip install git+https://github.com/JuergHauser/PyP223.git

import numpy
import matplotlib.pyplot as plt
from matplotlib.lines import Line2D
import arviz
import cofi

from vtem_max_forward_lib import (
    problem_setup,
    system_spec,
    survey_setup,
    true_model,
    ForwardWrapper,
    plot_transient,
    plot_plate_faces,
    plot_plate_faces_single
)

numpy.random.seed(42)

Background#

When modelling the electromagnetic response of subvertical bodies such as a VMS deposit they can be approximated using a thin plate in the halfspace of a layered earth, that is the forward solver computes the 3D response of a thin plate (e.g. Prikhodko et al. 2019). Here we develop a thin plate inversion method using CoFI to solve the inverse problem and P223 (Raiche et. al., 2007) to solve the forward problem.

Model parametrisation#

In the following we look at a thin plate with a conductance of \(2 \mathrm{S}\) located in a halfspace with a resistivity of \(1000 \mathrm{\Omega m}\) with a \(20 \mathrm{m}\) thick regolith that has a resistivity of \(300 \mathrm{\Omega m}\).

wpar8wl.png

wpar8wl.png#

The problem setup is imported from vtem_max_forward_lib.py but can be adjusted for other applications. The wrapper is created so that we can declare model parameters which are a subset of all the model parameters required by the forward solver. This allows to, for example, invert only for dip of the thin plate with all the other mode paremters assumed to be known. The thin plate is parameterised using the parametrisation introduced in (Hauser et. al. 2016). Compared to the commonly employed parametrisation with a plate reference point on the edge of the plate this parametrisation allows for a thin plate to grow and shrink around a plate refrerence point, without the need to move the reference point. This can be advantageous when there is for example a borehole intersecting a thin plate.

wpar7wl.png

wpar7wl.png#

Forward solver#

The forward solver is LeroiAir (Raiche et. al, 2007) and the code has been reorganised so that the response measured by an AEM system is give by a function that can be called from Python. In LeroiAir plates are discretised into cells, with the accuracy of the forward solver a function of the chosen cell-size. The forward solver is kept in a seperate Python package that is available here

Jacobian via finite differencing#

Parameter estimation methods frequently rely on the provision of a Jacobian for efficient optimisaiton. If an analytical Jacobian is not available it can be computed via finite differencing.

$ f’(x_0) = {f(x_0 +h)-f(x_0):raw-latex:`over `h} $

Care must be taken when choosing the step size \(h\) as a too small step size my result in a Jacobian that is affected by a limited accuracy of the forward solver and a too large step size \(h\) might result in a Jacobian that is not representative of the derivatives at location \(x_0\). In the following we use a relative step size \(q\) that is \(h=x0*(1.0+q)\). Further to this the gradient of the objective functions itself is affected by the noise on the data, thus for noisy data choosing a larger step size when computing the Jacobian can be advisable.

VTEM max Data#

Airborne electromagnetic systems can be categorised into either helicopter or fixed wing systems. The examples in this directory use a VTEM max system which is a helicopter based system developed and operated by Geotech.

https://geotech.ca/services/electromagnetic/vtem-versatile-time-domain-electromagnetic-system/

Successful inversion also relies on the objective function being smooth and predictable. For the data being inverted here it is advantageous to convert measurements to scale logarithmically to obtain a smoother and more predictable objective function when compared with using the unscaled data. Similarly plate orientation angels are converted into radians.

Further reading#

Hauser, J., Gunning, J., & Annetts, D. (2016). Probabilistic inversion of airborne electromagnetic data for basement conductors. Geophysics, 81(5), E389-E400.

Prikhodko, A., Morrison, E., Bagrianski, A., Kuzmin, P., Tishin, P., & Legault, J. (2010). Evolution of VTEM? technical solutions for effective exploration. ASEG Extended Abstracts, 2010(1), 1-4.

Raiche, A., Sugeng, F. and Wilson, G. (2007) Practical 3D EM inversion the P223F software suite, ASEG Extended Abstracts, 2007:1, 1-5

Wheelock, B., Constable, S., & Key, K. (2015). The advantages of logarithmically scaled data for electromagnetic inversion. Geophysical Journal International, 201(3), 1765–1780. https://doi.org/10.1093/GJI/GGV107

Problem definition#

survey_setup = {
    "tx": numpy.array([205.]),                  # transmitter easting/x-position
    "ty": numpy.array([100.]),                  # transmitter northing/y-position
    "tz": numpy.array([50.]),                   # transmitter height/z-position
    "tazi": numpy.deg2rad(numpy.array([90.])),  # transmitter azimuth
    "tincl": numpy.deg2rad(numpy.array([6.])),  # transmitter inclination
    "rx": numpy.array([205.]),                  # receiver easting/x-position
    "ry": numpy.array([100.]),                  # receiver northing/y-position
    "rz": numpy.array([50.]),                   # receiver height/z-position
    "trdx": numpy.array([0.]),                  # transmitter receiver separation inline
    "trdy": numpy.array([0.]),                  # transmitter receiver separation crossline
    "trdz": numpy.array([0.]),                  # transmitter receiver separation vertical
}

forward = ForwardWrapper(true_model, problem_setup, system_spec, survey_setup, ["pdip"])
true_param_value = numpy.array([60])

['pdip']

True model#

_, axes = plt.subplots(2, 2)
axes[1,1].axis("off")
plot_plate_faces(
    "plate_true", forward, true_param_value,
    axes[0,0], axes[0,1], axes[1,0], color="purple", label="True model"
)
plt.tight_layout()
point = Line2D([0], [0], label='Fiducial', marker='o', markersize=5,
         markeredgecolor='orange', markerfacecolor='orange', linestyle='')

handles, labels = axes[1,0].get_legend_handles_labels()
handles.extend([point])

axes[1,0].legend(handles=handles,bbox_to_anchor=(1.04, 0), loc="lower left")

<matplotlib.legend.Legend object at 0x7efdb4cfe350>

Generate synthetic data#

We use a simplfied noise model that assumes an absolute noise that is a standard deviation of \(0.05\) for the logarithms of the measured and observed data.

# The data
absolute_noise= 0.05

# create data and ad a realisation of the noise
data_pred_true = forward(true_param_value)
data_obs = data_pred_true + numpy.random.randn(len(data_pred_true))*absolute_noise

# define data covariance matrix
sigma=absolute_noise
Cdinv=numpy.identity(len(data_obs))*(1.0/(sigma*sigma))

Starting model#

Set an initial guess for the dip of the thin plate

init_param_value = numpy.array([45])

Define helper functions for CoFI#

def my_objective(model):
    dpred = forward(model)
    residual = dpred - data_obs
    return residual.T @ Cdinv @ residual

def my_gradient(model):
    dpred = forward(model)
    jacobian = forward.jacobian(model, relative_step=0.1)
    residual = dpred - data_obs
    return jacobian.T @ Cdinv @ residual

def my_hessian(model):
    jacobian = forward.jacobian(model)
    return jacobian.T @ Cdinv @ jacobian



class PerIterationCallbackFunction:
    def __init__(self):
        self.x = None
        self.i = 0

    def __call__(self, xk):
        print(f"Iteration #{self.i+1}")
        print(f"  objective value: {my_problem.objective(xk)}")
        self.x = xk
        self.i += 1

Sensitivity of the misfit function and usefulness of the gradient#

When setting up a new inverse problem prior to any inversion it makes sense to verify that the misfit function is sensitve to the parameter of itnerest and that the gradient of the objective function points in the right direction.

all_models = [numpy.array([pdip]) for pdip in range(40, 120, 5)]
all_misfits = []
all_gradients = []
for model in all_models:
    misfit = my_objective(model)
    gradient = my_gradient(model)
    all_misfits.append(misfit)
    all_gradients.append(gradient)
    print(f"pdip: {model}, data misfit: {misfit}, gradient: {gradient}")


fig, ax1 = plt.subplots()
color = 'tab:red'
ax1.plot(all_models, all_misfits,color=color)
ax1.tick_params(axis='y', labelcolor=color)
ax1.set_xlabel("pdip")
ax1.set_ylabel("Data misfit",color=color)

ax2 = ax1.twinx()
color = 'tab:blue'
ax2.plot(all_models, all_gradients,color=color)
ax2.set_ylabel('Gradient', color=color)
ax2.tick_params(axis='y', labelcolor=color)

fig.tight_layout()  # otherwise the right y-label is slightly clipped
plt.show()

pdip: [40], data misfit: 231.708835656254, gradient: [-8.38525334]
pdip: [45], data misfit: 161.52059630302878, gradient: [-5.17711959]
pdip: [50], data misfit: 115.89689554155724, gradient: [-3.41393833]
pdip: [55], data misfit: 91.08438559558219, gradient: [-1.79238028]
pdip: [60], data misfit: 78.13208270691166, gradient: [-0.40612889]
pdip: [65], data misfit: 79.47542276986738, gradient: [0.60701868]
pdip: [70], data misfit: 86.57794284916837, gradient: [0.77951888]
pdip: [75], data misfit: 99.82377503892143, gradient: [1.04776418]
pdip: [80], data misfit: 114.82121927725917, gradient: [1.88028555]
pdip: [85], data misfit: 135.80364266137616, gradient: [1.99672184]
pdip: [90], data misfit: 158.91312196239573, gradient: [2.12201565]
pdip: [95], data misfit: 190.36506697630017, gradient: [1.90338506]
pdip: [100], data misfit: 211.34574782298415, gradient: [3.0843773]
pdip: [105], data misfit: 248.95857134110955, gradient: [3.54243063]
pdip: [110], data misfit: 283.9454162116735, gradient: [3.73187192]
pdip: [115], data misfit: 338.73968428973234, gradient: [3.35857177]

Parameter estimation#

First we solve the inverse problem using optimisation that is we seek to find the minimum of the objective function given as

\[\chi^2 = (\mathbf{d} - \mathbf{f}(\mathbf{m}))^T\mathbf{C}_d^{-1}(\mathbf{d}-\mathbf{f}(\mathbf{m})),\]

with the full Newton step being

\[\begin{equation} \Delta \mathbf{m}= (\underbrace{\mathbf{J}^T \mathbf{C}_d^{-1} \mathbf{J}}_{\mathbf{Hessian}})^{-1} (\underbrace{ \mathbf{J}^T\mathbf{C}_d^{-1} (\mathbf{y}-\mathbf{f}(\mathbf{m}))}_\mathbf{Gradient}). \end{equation}\]
my_problem = cofi.BaseProblem()
my_problem.set_objective(my_objective)
my_problem.set_gradient(my_gradient)
my_problem.set_hessian(my_hessian)
my_problem.set_initial_model(init_param_value)

my_options = cofi.InversionOptions()
my_options.set_tool("scipy.optimize.minimize")
my_options.set_params(method="Newton-CG",callback=PerIterationCallbackFunction())

my_inversion = cofi.Inversion(my_problem, my_options)
my_result = my_inversion.run()
print(my_result.model)

Iteration #1
  objective value: 77.89136366904667
Iteration #2
  objective value: 77.67785905000744
[59.18480808]

_, (ax1, ax2) = plt.subplots(1, 2)
plot_transient(true_param_value, forward, "Data from true model", ax1, ax2, color="purple")
plot_transient(init_param_value, forward, "Data from starting model", ax1, ax2, color="green", linestyle=":")
plot_transient(my_result.model, forward, "Data from MAP model", ax1, ax2, color="red", linestyle="-.")
ax1.legend(loc="upper center")
ax2.legend(loc="upper center")
ax1.set_title("vertical")
ax2.set_title("inline")
plt.tight_layout()

_, axes = plt.subplots(2, 2)
axes[1,1].axis("off")
plot_plate_faces(
    "plate_true", forward, true_param_value,
    axes[0,0], axes[0,1], axes[1,0], color="purple", label="True model"
)
plot_plate_faces(
    "plate_init", forward, init_param_value,
    axes[0,0], axes[0,1], axes[1,0], color="green", label="Starting model"
)
plot_plate_faces(
    "plate_inverted", forward, my_result.model,
    axes[0,0], axes[0,1], axes[1,0], color="red", label="MAP solution", linestyle="dotted"
)
plt.tight_layout()
point = Line2D([0], [0], label='Fiducial', marker='o', markersize=5,
         markeredgecolor='orange', markerfacecolor='orange', linestyle='')

handles, labels = axes[1,0].get_legend_handles_labels()
handles.extend([point])

axes[1,0].legend(handles=handles,bbox_to_anchor=(1.04, 0), loc="lower left")

<matplotlib.legend.Legend object at 0x7efdb4cfd1d0>

Ensemble method#

Parameter estimation methods require an objective function while ensemble methods require a likelihood functions typically given in the form of a log likelidhood function. The objective function used for the parameter estimation consists of only a data misfit term and thus is closely related to the likelihood function, with the log likelihood being proportional to the value of the objective function multiplied by a factor of \(\frac{1}{2}\)

\[p({\mathbf d} | {\mathbf m}) \propto \exp \left\{- \frac{1}{2} ({\mathbf d}-{\mathbf f}({\mathbf m}))^T C_d^{-1} ({\mathbf d}-{\mathbf f}({\mathbf m})) \right\}\]

In the following we define a log likelihood function and log prior function. The prior distribution is a uniform distribution with an lower boundary of \(10 \degree\) and an upper boundary of \(80 \degree\)

def my_log_likelihood(model):
    return -0.5 * my_objective(model)

def my_log_prior(m):    # uniform distribution
    for i in range(len(m)):
        if m[i] < m_min[i] or m[i] > m_max[i]: return -numpy.inf
    return 0.0 # model lies within bounds -> return log(1)

m_min=numpy.array([10])
m_max=numpy.array([80])

To be able to use an ensemble method we need to augment our CoFI problem with a function providing the log of the prior probability and a second function that provides the log of the likelihood function.

my_problem.set_log_prior(my_log_prior)
my_problem.set_log_likelihood(my_log_likelihood)
my_problem.set_model_shape(len(init_param_value))

nwalkers = 3
ndim = len(init_param_value)
nsteps = 500
walkers_start = init_param_value + 1 * numpy.random.randn(nwalkers, ndim)

inv_options = cofi.InversionOptions()
inv_options.set_tool("emcee")
inv_options.set_params(nwalkers=nwalkers, nsteps=nsteps, initial_state=walkers_start, progress=True)

######## Run it
inv = cofi.Inversion(my_problem, inv_options)
inv_result = inv.run()

######## Check result
print(f"The inversion result from `emcee`:")
inv_result.summary()

  0%|          | 0/500 [00:00<?, ?it/s]
  0%|          | 1/500 [00:00<06:52,  1.21it/s]
  0%|          | 2/500 [00:01<06:52,  1.21it/s]
  1%|          | 3/500 [00:02<06:51,  1.21it/s]
  1%|          | 4/500 [00:03<06:49,  1.21it/s]
  1%|          | 5/500 [00:04<06:48,  1.21it/s]
  1%|          | 6/500 [00:04<06:51,  1.20it/s]
  1%|▏         | 7/500 [00:05<06:51,  1.20it/s]
  2%|▏         | 8/500 [00:06<06:51,  1.20it/s]
  2%|▏         | 9/500 [00:07<06:51,  1.19it/s]
  2%|▏         | 10/500 [00:08<06:49,  1.20it/s]
  2%|▏         | 11/500 [00:09<06:48,  1.20it/s]
  2%|▏         | 12/500 [00:09<06:46,  1.20it/s]
  3%|▎         | 13/500 [00:10<06:43,  1.21it/s]
  3%|▎         | 14/500 [00:11<06:41,  1.21it/s]
  3%|▎         | 15/500 [00:12<06:40,  1.21it/s]
  3%|▎         | 16/500 [00:13<06:38,  1.21it/s]
  3%|▎         | 17/500 [00:14<06:37,  1.21it/s]
  4%|▎         | 18/500 [00:14<06:36,  1.21it/s]
  4%|▍         | 19/500 [00:15<06:35,  1.22it/s]
  4%|▍         | 20/500 [00:16<06:34,  1.22it/s]
  4%|▍         | 21/500 [00:17<06:33,  1.22it/s]
  4%|▍         | 22/500 [00:18<06:32,  1.22it/s]
  5%|▍         | 23/500 [00:19<06:31,  1.22it/s]
  5%|▍         | 24/500 [00:19<06:31,  1.22it/s]
  5%|▌         | 25/500 [00:20<06:31,  1.21it/s]
  5%|▌         | 26/500 [00:21<06:31,  1.21it/s]
  5%|▌         | 27/500 [00:22<06:32,  1.21it/s]
  6%|▌         | 28/500 [00:23<06:31,  1.21it/s]
  6%|▌         | 29/500 [00:23<06:29,  1.21it/s]
  6%|▌         | 30/500 [00:24<06:27,  1.21it/s]
  6%|▌         | 31/500 [00:25<06:27,  1.21it/s]
  6%|▋         | 32/500 [00:26<06:26,  1.21it/s]
  7%|▋         | 33/500 [00:27<06:27,  1.21it/s]
  7%|▋         | 34/500 [00:28<06:29,  1.20it/s]
  7%|▋         | 35/500 [00:28<06:28,  1.20it/s]
  7%|▋         | 36/500 [00:29<06:27,  1.20it/s]
  7%|▋         | 37/500 [00:30<06:27,  1.19it/s]
  8%|▊         | 38/500 [00:31<06:25,  1.20it/s]
  8%|▊         | 39/500 [00:32<06:25,  1.19it/s]
  8%|▊         | 40/500 [00:33<06:25,  1.19it/s]
  8%|▊         | 41/500 [00:34<06:23,  1.20it/s]
  8%|▊         | 42/500 [00:34<06:21,  1.20it/s]
  9%|▊         | 43/500 [00:35<06:21,  1.20it/s]
  9%|▉         | 44/500 [00:36<06:20,  1.20it/s]
  9%|▉         | 45/500 [00:37<06:18,  1.20it/s]
  9%|▉         | 46/500 [00:38<06:19,  1.20it/s]
  9%|▉         | 47/500 [00:39<06:16,  1.20it/s]
 10%|▉         | 48/500 [00:39<06:15,  1.20it/s]
 10%|▉         | 49/500 [00:40<06:14,  1.20it/s]
 10%|█         | 50/500 [00:41<06:12,  1.21it/s]
 10%|█         | 51/500 [00:42<06:10,  1.21it/s]
 10%|█         | 52/500 [00:43<06:09,  1.21it/s]
 11%|█         | 53/500 [00:43<06:11,  1.20it/s]
 11%|█         | 54/500 [00:44<06:12,  1.20it/s]
 11%|█         | 55/500 [00:45<06:11,  1.20it/s]
 11%|█         | 56/500 [00:46<06:09,  1.20it/s]
 11%|█▏        | 57/500 [00:47<06:08,  1.20it/s]
 12%|█▏        | 58/500 [00:48<06:09,  1.20it/s]
 12%|█▏        | 59/500 [00:48<06:07,  1.20it/s]
 12%|█▏        | 60/500 [00:49<06:05,  1.20it/s]
 12%|█▏        | 61/500 [00:50<06:04,  1.20it/s]
 12%|█▏        | 62/500 [00:51<06:03,  1.20it/s]
 13%|█▎        | 63/500 [00:52<06:02,  1.20it/s]
 13%|█▎        | 64/500 [00:53<06:02,  1.20it/s]
 13%|█▎        | 65/500 [00:53<06:02,  1.20it/s]
 13%|█▎        | 66/500 [00:54<06:00,  1.20it/s]
 13%|█▎        | 67/500 [00:55<05:58,  1.21it/s]
 14%|█▎        | 68/500 [00:56<05:59,  1.20it/s]
 14%|█▍        | 69/500 [00:57<05:57,  1.21it/s]
 14%|█▍        | 70/500 [00:58<05:55,  1.21it/s]
 14%|█▍        | 71/500 [00:58<05:53,  1.21it/s]
 14%|█▍        | 72/500 [00:59<05:51,  1.22it/s]
 15%|█▍        | 73/500 [01:00<05:49,  1.22it/s]
 15%|█▍        | 74/500 [01:01<05:49,  1.22it/s]
 15%|█▌        | 75/500 [01:02<05:47,  1.22it/s]
 15%|█▌        | 76/500 [01:03<05:46,  1.22it/s]
 15%|█▌        | 77/500 [01:03<05:46,  1.22it/s]
 16%|█▌        | 78/500 [01:04<05:45,  1.22it/s]
 16%|█▌        | 79/500 [01:05<05:43,  1.23it/s]
 16%|█▌        | 80/500 [01:06<05:42,  1.23it/s]
 16%|█▌        | 81/500 [01:07<05:41,  1.23it/s]
 16%|█▋        | 82/500 [01:07<05:41,  1.22it/s]
 17%|█▋        | 83/500 [01:08<05:43,  1.21it/s]
 17%|█▋        | 84/500 [01:09<05:42,  1.22it/s]
 17%|█▋        | 85/500 [01:10<05:41,  1.22it/s]
 17%|█▋        | 86/500 [01:11<05:40,  1.22it/s]
 17%|█▋        | 87/500 [01:12<05:38,  1.22it/s]
 18%|█▊        | 88/500 [01:12<05:38,  1.22it/s]
 18%|█▊        | 89/500 [01:13<05:38,  1.21it/s]
 18%|█▊        | 90/500 [01:14<05:36,  1.22it/s]
 18%|█▊        | 91/500 [01:15<05:35,  1.22it/s]
 18%|█▊        | 92/500 [01:16<05:35,  1.22it/s]
 19%|█▊        | 93/500 [01:16<05:36,  1.21it/s]
 19%|█▉        | 94/500 [01:17<05:34,  1.21it/s]
 19%|█▉        | 95/500 [01:18<05:34,  1.21it/s]
 19%|█▉        | 96/500 [01:19<05:33,  1.21it/s]
 19%|█▉        | 97/500 [01:20<05:31,  1.21it/s]
 20%|█▉        | 98/500 [01:21<05:31,  1.21it/s]
 20%|█▉        | 99/500 [01:21<05:31,  1.21it/s]
 20%|██        | 100/500 [01:22<05:30,  1.21it/s]
 20%|██        | 101/500 [01:23<05:28,  1.21it/s]
 20%|██        | 102/500 [01:24<05:27,  1.21it/s]
 21%|██        | 103/500 [01:25<05:26,  1.22it/s]
 21%|██        | 104/500 [01:26<05:25,  1.22it/s]
 21%|██        | 105/500 [01:26<05:25,  1.21it/s]
 21%|██        | 106/500 [01:27<05:24,  1.21it/s]
 21%|██▏       | 107/500 [01:28<05:22,  1.22it/s]
 22%|██▏       | 108/500 [01:29<05:23,  1.21it/s]
 22%|██▏       | 109/500 [01:30<05:22,  1.21it/s]
 22%|██▏       | 110/500 [01:30<05:21,  1.21it/s]
 22%|██▏       | 111/500 [01:31<05:24,  1.20it/s]
 22%|██▏       | 112/500 [01:32<05:23,  1.20it/s]
 23%|██▎       | 113/500 [01:33<05:21,  1.21it/s]
 23%|██▎       | 114/500 [01:34<05:21,  1.20it/s]
 23%|██▎       | 115/500 [01:35<05:19,  1.21it/s]
 23%|██▎       | 116/500 [01:35<05:18,  1.21it/s]
 23%|██▎       | 117/500 [01:36<05:18,  1.20it/s]
 24%|██▎       | 118/500 [01:37<05:16,  1.21it/s]
 24%|██▍       | 119/500 [01:38<05:16,  1.21it/s]
 24%|██▍       | 120/500 [01:39<05:15,  1.20it/s]
 24%|██▍       | 121/500 [01:40<05:13,  1.21it/s]
 24%|██▍       | 122/500 [01:40<05:12,  1.21it/s]
 25%|██▍       | 123/500 [01:41<05:12,  1.21it/s]
 25%|██▍       | 124/500 [01:42<05:11,  1.21it/s]
 25%|██▌       | 125/500 [01:43<05:10,  1.21it/s]
 25%|██▌       | 126/500 [01:44<05:10,  1.20it/s]
 25%|██▌       | 127/500 [01:45<05:10,  1.20it/s]
 26%|██▌       | 128/500 [01:45<05:08,  1.21it/s]
 26%|██▌       | 129/500 [01:46<05:07,  1.21it/s]
 26%|██▌       | 130/500 [01:47<05:05,  1.21it/s]
 26%|██▌       | 131/500 [01:48<05:05,  1.21it/s]
 26%|██▋       | 132/500 [01:49<05:04,  1.21it/s]
 27%|██▋       | 133/500 [01:50<05:03,  1.21it/s]
 27%|██▋       | 134/500 [01:50<05:02,  1.21it/s]
 27%|██▋       | 135/500 [01:51<05:01,  1.21it/s]
 27%|██▋       | 136/500 [01:52<05:00,  1.21it/s]
 27%|██▋       | 137/500 [01:53<05:00,  1.21it/s]
 28%|██▊       | 138/500 [01:54<04:59,  1.21it/s]
 28%|██▊       | 139/500 [01:55<04:59,  1.21it/s]
 28%|██▊       | 140/500 [01:55<04:57,  1.21it/s]
 28%|██▊       | 141/500 [01:56<04:56,  1.21it/s]
 28%|██▊       | 142/500 [01:57<04:56,  1.21it/s]
 29%|██▊       | 143/500 [01:58<04:55,  1.21it/s]
 29%|██▉       | 144/500 [01:59<04:54,  1.21it/s]
 29%|██▉       | 145/500 [01:59<04:54,  1.21it/s]
 29%|██▉       | 146/500 [02:00<04:53,  1.21it/s]
 29%|██▉       | 147/500 [02:01<04:51,  1.21it/s]
 30%|██▉       | 148/500 [02:02<04:51,  1.21it/s]
 30%|██▉       | 149/500 [02:03<04:51,  1.20it/s]
 30%|███       | 150/500 [02:04<04:51,  1.20it/s]
 30%|███       | 151/500 [02:05<04:54,  1.18it/s]
 30%|███       | 152/500 [02:05<04:54,  1.18it/s]
 31%|███       | 153/500 [02:06<04:53,  1.18it/s]
 31%|███       | 154/500 [02:07<04:54,  1.17it/s]
 31%|███       | 155/500 [02:08<04:54,  1.17it/s]
 31%|███       | 156/500 [02:09<04:53,  1.17it/s]
 31%|███▏      | 157/500 [02:10<04:54,  1.17it/s]
 32%|███▏      | 158/500 [02:11<04:53,  1.17it/s]
 32%|███▏      | 159/500 [02:11<04:51,  1.17it/s]
 32%|███▏      | 160/500 [02:12<04:52,  1.16it/s]
 32%|███▏      | 161/500 [02:13<04:50,  1.17it/s]
 32%|███▏      | 162/500 [02:14<04:48,  1.17it/s]
 33%|███▎      | 163/500 [02:15<04:47,  1.17it/s]
 33%|███▎      | 164/500 [02:16<04:45,  1.18it/s]
 33%|███▎      | 165/500 [02:16<04:44,  1.18it/s]
 33%|███▎      | 166/500 [02:17<04:43,  1.18it/s]
 33%|███▎      | 167/500 [02:18<04:43,  1.17it/s]
 34%|███▎      | 168/500 [02:19<04:44,  1.17it/s]
 34%|███▍      | 169/500 [02:20<04:43,  1.17it/s]
 34%|███▍      | 170/500 [02:21<04:42,  1.17it/s]
 34%|███▍      | 171/500 [02:22<04:40,  1.17it/s]
 34%|███▍      | 172/500 [02:22<04:39,  1.18it/s]
 35%|███▍      | 173/500 [02:23<04:37,  1.18it/s]
 35%|███▍      | 174/500 [02:24<04:36,  1.18it/s]
 35%|███▌      | 175/500 [02:25<04:35,  1.18it/s]
 35%|███▌      | 176/500 [02:26<04:34,  1.18it/s]
 35%|███▌      | 177/500 [02:27<04:32,  1.18it/s]
 36%|███▌      | 178/500 [02:28<04:31,  1.19it/s]
 36%|███▌      | 179/500 [02:28<04:30,  1.18it/s]
 36%|███▌      | 180/500 [02:29<04:29,  1.19it/s]
 36%|███▌      | 181/500 [02:30<04:29,  1.18it/s]
 36%|███▋      | 182/500 [02:31<04:29,  1.18it/s]
 37%|███▋      | 183/500 [02:32<04:27,  1.18it/s]
 37%|███▋      | 184/500 [02:33<04:26,  1.19it/s]
 37%|███▋      | 185/500 [02:33<04:24,  1.19it/s]
 37%|███▋      | 186/500 [02:34<04:22,  1.20it/s]
 37%|███▋      | 187/500 [02:35<04:21,  1.19it/s]
 38%|███▊      | 188/500 [02:36<04:20,  1.20it/s]
 38%|███▊      | 189/500 [02:37<04:19,  1.20it/s]
 38%|███▊      | 190/500 [02:38<04:17,  1.20it/s]
 38%|███▊      | 191/500 [02:38<04:17,  1.20it/s]
 38%|███▊      | 192/500 [02:39<04:15,  1.20it/s]
 39%|███▊      | 193/500 [02:40<04:14,  1.20it/s]
 39%|███▉      | 194/500 [02:41<04:13,  1.21it/s]
 39%|███▉      | 195/500 [02:42<04:12,  1.21it/s]
 39%|███▉      | 196/500 [02:43<04:11,  1.21it/s]
 39%|███▉      | 197/500 [02:43<04:10,  1.21it/s]
 40%|███▉      | 198/500 [02:44<04:10,  1.21it/s]
 40%|███▉      | 199/500 [02:45<04:10,  1.20it/s]
 40%|████      | 200/500 [02:46<04:08,  1.21it/s]
 40%|████      | 201/500 [02:47<04:07,  1.21it/s]
 40%|████      | 202/500 [02:48<04:06,  1.21it/s]
 41%|████      | 203/500 [02:48<04:07,  1.20it/s]
 41%|████      | 204/500 [02:49<04:06,  1.20it/s]
 41%|████      | 205/500 [02:50<04:06,  1.19it/s]
 41%|████      | 206/500 [02:51<04:08,  1.18it/s]
 41%|████▏     | 207/500 [02:52<04:06,  1.19it/s]
 42%|████▏     | 208/500 [02:53<04:05,  1.19it/s]
 42%|████▏     | 209/500 [02:53<04:04,  1.19it/s]
 42%|████▏     | 210/500 [02:54<04:03,  1.19it/s]
 42%|████▏     | 211/500 [02:55<04:01,  1.20it/s]
 42%|████▏     | 212/500 [02:56<04:00,  1.20it/s]
 43%|████▎     | 213/500 [02:57<03:58,  1.20it/s]
 43%|████▎     | 214/500 [02:58<03:57,  1.21it/s]
 43%|████▎     | 215/500 [02:58<03:56,  1.21it/s]
 43%|████▎     | 216/500 [02:59<03:55,  1.20it/s]
 43%|████▎     | 217/500 [03:00<03:55,  1.20it/s]
 44%|████▎     | 218/500 [03:01<03:56,  1.19it/s]
 44%|████▍     | 219/500 [03:02<03:54,  1.20it/s]
 44%|████▍     | 220/500 [03:03<03:52,  1.21it/s]
 44%|████▍     | 221/500 [03:03<03:52,  1.20it/s]
 44%|████▍     | 222/500 [03:04<03:53,  1.19it/s]
 45%|████▍     | 223/500 [03:05<03:52,  1.19it/s]
 45%|████▍     | 224/500 [03:06<03:53,  1.18it/s]
 45%|████▌     | 225/500 [03:07<03:56,  1.16it/s]
 45%|████▌     | 226/500 [03:08<03:54,  1.17it/s]
 45%|████▌     | 227/500 [03:09<03:54,  1.16it/s]
 46%|████▌     | 228/500 [03:09<03:51,  1.18it/s]
 46%|████▌     | 229/500 [03:10<03:48,  1.18it/s]
 46%|████▌     | 230/500 [03:11<03:48,  1.18it/s]
 46%|████▌     | 231/500 [03:12<03:46,  1.19it/s]
 46%|████▋     | 232/500 [03:13<03:44,  1.19it/s]
 47%|████▋     | 233/500 [03:14<03:43,  1.19it/s]
 47%|████▋     | 234/500 [03:14<03:42,  1.19it/s]
 47%|████▋     | 235/500 [03:15<03:40,  1.20it/s]
 47%|████▋     | 236/500 [03:16<03:41,  1.19it/s]
 47%|████▋     | 237/500 [03:17<03:40,  1.19it/s]
 48%|████▊     | 238/500 [03:18<03:39,  1.20it/s]
 48%|████▊     | 239/500 [03:19<03:37,  1.20it/s]
 48%|████▊     | 240/500 [03:19<03:37,  1.19it/s]
 48%|████▊     | 241/500 [03:20<03:37,  1.19it/s]
 48%|████▊     | 242/500 [03:21<03:38,  1.18it/s]
 49%|████▊     | 243/500 [03:22<03:37,  1.18it/s]
 49%|████▉     | 244/500 [03:23<03:36,  1.18it/s]
 49%|████▉     | 245/500 [03:24<03:36,  1.18it/s]
 49%|████▉     | 246/500 [03:25<03:36,  1.17it/s]
 49%|████▉     | 247/500 [03:25<03:35,  1.17it/s]
 50%|████▉     | 248/500 [03:26<03:35,  1.17it/s]
 50%|████▉     | 249/500 [03:27<03:33,  1.18it/s]
 50%|█████     | 250/500 [03:28<03:30,  1.19it/s]
 50%|█████     | 251/500 [03:29<03:29,  1.19it/s]
 50%|█████     | 252/500 [03:30<03:29,  1.18it/s]
 51%|█████     | 253/500 [03:30<03:27,  1.19it/s]
 51%|█████     | 254/500 [03:31<03:27,  1.19it/s]
 51%|█████     | 255/500 [03:32<03:25,  1.19it/s]
 51%|█████     | 256/500 [03:33<03:24,  1.19it/s]
 51%|█████▏    | 257/500 [03:34<03:23,  1.20it/s]
 52%|█████▏    | 258/500 [03:35<03:22,  1.19it/s]
 52%|█████▏    | 259/500 [03:35<03:21,  1.20it/s]
 52%|█████▏    | 260/500 [03:36<03:21,  1.19it/s]
 52%|█████▏    | 261/500 [03:37<03:20,  1.19it/s]
 52%|█████▏    | 262/500 [03:38<03:19,  1.19it/s]
 53%|█████▎    | 263/500 [03:39<03:18,  1.19it/s]
 53%|█████▎    | 264/500 [03:40<03:18,  1.19it/s]
 53%|█████▎    | 265/500 [03:40<03:17,  1.19it/s]
 53%|█████▎    | 266/500 [03:41<03:17,  1.19it/s]
 53%|█████▎    | 267/500 [03:42<03:15,  1.19it/s]
 54%|█████▎    | 268/500 [03:43<03:14,  1.19it/s]
 54%|█████▍    | 269/500 [03:44<03:13,  1.20it/s]
 54%|█████▍    | 270/500 [03:45<03:12,  1.19it/s]
 54%|█████▍    | 271/500 [03:46<03:11,  1.20it/s]
 54%|█████▍    | 272/500 [03:46<03:10,  1.20it/s]
 55%|█████▍    | 273/500 [03:47<03:09,  1.19it/s]
 55%|█████▍    | 274/500 [03:48<03:08,  1.20it/s]
 55%|█████▌    | 275/500 [03:49<03:09,  1.19it/s]
 55%|█████▌    | 276/500 [03:50<03:10,  1.18it/s]
 55%|█████▌    | 277/500 [03:51<03:07,  1.19it/s]
 56%|█████▌    | 278/500 [03:51<03:06,  1.19it/s]
 56%|█████▌    | 279/500 [03:52<03:06,  1.18it/s]
 56%|█████▌    | 280/500 [03:53<03:04,  1.19it/s]
 56%|█████▌    | 281/500 [03:54<03:03,  1.20it/s]
 56%|█████▋    | 282/500 [03:55<03:02,  1.20it/s]
 57%|█████▋    | 283/500 [03:56<03:00,  1.20it/s]
 57%|█████▋    | 284/500 [03:56<02:59,  1.20it/s]
 57%|█████▋    | 285/500 [03:57<03:00,  1.19it/s]
 57%|█████▋    | 286/500 [03:58<02:58,  1.20it/s]
 57%|█████▋    | 287/500 [03:59<02:57,  1.20it/s]
 58%|█████▊    | 288/500 [04:00<02:56,  1.20it/s]
 58%|█████▊    | 289/500 [04:01<02:55,  1.20it/s]
 58%|█████▊    | 290/500 [04:01<02:53,  1.21it/s]
 58%|█████▊    | 291/500 [04:02<02:53,  1.21it/s]
 58%|█████▊    | 292/500 [04:03<02:52,  1.21it/s]
 59%|█████▊    | 293/500 [04:04<02:50,  1.21it/s]
 59%|█████▉    | 294/500 [04:05<02:50,  1.21it/s]
 59%|█████▉    | 295/500 [04:06<02:49,  1.21it/s]
 59%|█████▉    | 296/500 [04:06<02:48,  1.21it/s]
 59%|█████▉    | 297/500 [04:07<02:47,  1.21it/s]
 60%|█████▉    | 298/500 [04:08<02:46,  1.21it/s]
 60%|█████▉    | 299/500 [04:09<02:45,  1.21it/s]
 60%|██████    | 300/500 [04:10<02:45,  1.21it/s]
 60%|██████    | 301/500 [04:10<02:44,  1.21it/s]
 60%|██████    | 302/500 [04:11<02:43,  1.21it/s]
 61%|██████    | 303/500 [04:12<02:43,  1.21it/s]
 61%|██████    | 304/500 [04:13<02:42,  1.21it/s]
 61%|██████    | 305/500 [04:14<02:41,  1.21it/s]
 61%|██████    | 306/500 [04:15<02:40,  1.21it/s]
 61%|██████▏   | 307/500 [04:15<02:39,  1.21it/s]
 62%|██████▏   | 308/500 [04:16<02:38,  1.21it/s]
 62%|██████▏   | 309/500 [04:17<02:37,  1.21it/s]
 62%|██████▏   | 310/500 [04:18<02:37,  1.21it/s]
 62%|██████▏   | 311/500 [04:19<02:36,  1.21it/s]
 62%|██████▏   | 312/500 [04:20<02:35,  1.21it/s]
 63%|██████▎   | 313/500 [04:20<02:34,  1.21it/s]
 63%|██████▎   | 314/500 [04:21<02:33,  1.21it/s]
 63%|██████▎   | 315/500 [04:22<02:32,  1.21it/s]
 63%|██████▎   | 316/500 [04:23<02:32,  1.21it/s]
 63%|██████▎   | 317/500 [04:24<02:31,  1.21it/s]
 64%|██████▎   | 318/500 [04:25<02:30,  1.21it/s]
 64%|██████▍   | 319/500 [04:25<02:29,  1.21it/s]
 64%|██████▍   | 320/500 [04:26<02:28,  1.21it/s]
 64%|██████▍   | 321/500 [04:27<02:27,  1.21it/s]
 64%|██████▍   | 322/500 [04:28<02:27,  1.21it/s]
 65%|██████▍   | 323/500 [04:29<02:26,  1.21it/s]
 65%|██████▍   | 324/500 [04:29<02:25,  1.21it/s]
 65%|██████▌   | 325/500 [04:30<02:24,  1.21it/s]
 65%|██████▌   | 326/500 [04:31<02:23,  1.21it/s]
 65%|██████▌   | 327/500 [04:32<02:22,  1.22it/s]
 66%|██████▌   | 328/500 [04:33<02:22,  1.21it/s]
 66%|██████▌   | 329/500 [04:34<02:20,  1.21it/s]
 66%|██████▌   | 330/500 [04:34<02:19,  1.22it/s]
 66%|██████▌   | 331/500 [04:35<02:18,  1.22it/s]
 66%|██████▋   | 332/500 [04:36<02:18,  1.21it/s]
 67%|██████▋   | 333/500 [04:37<02:19,  1.20it/s]
 67%|██████▋   | 334/500 [04:38<02:18,  1.20it/s]
 67%|██████▋   | 335/500 [04:39<02:16,  1.21it/s]
 67%|██████▋   | 336/500 [04:39<02:15,  1.21it/s]
 67%|██████▋   | 337/500 [04:40<02:14,  1.21it/s]
 68%|██████▊   | 338/500 [04:41<02:13,  1.21it/s]
 68%|██████▊   | 339/500 [04:42<02:12,  1.21it/s]
 68%|██████▊   | 340/500 [04:43<02:12,  1.21it/s]
 68%|██████▊   | 341/500 [04:44<02:11,  1.21it/s]
 68%|██████▊   | 342/500 [04:44<02:10,  1.21it/s]
 69%|██████▊   | 343/500 [04:45<02:08,  1.22it/s]
 69%|██████▉   | 344/500 [04:46<02:08,  1.22it/s]
 69%|██████▉   | 345/500 [04:47<02:07,  1.22it/s]
 69%|██████▉   | 346/500 [04:48<02:06,  1.22it/s]
 69%|██████▉   | 347/500 [04:48<02:06,  1.21it/s]
 70%|██████▉   | 348/500 [04:49<02:05,  1.21it/s]
 70%|██████▉   | 349/500 [04:50<02:05,  1.21it/s]
 70%|███████   | 350/500 [04:51<02:04,  1.20it/s]
 70%|███████   | 351/500 [04:52<02:03,  1.21it/s]
 70%|███████   | 352/500 [04:53<02:02,  1.21it/s]
 71%|███████   | 353/500 [04:53<02:01,  1.21it/s]
 71%|███████   | 354/500 [04:54<02:00,  1.21it/s]
 71%|███████   | 355/500 [04:55<01:59,  1.21it/s]
 71%|███████   | 356/500 [04:56<01:58,  1.21it/s]
 71%|███████▏  | 357/500 [04:57<01:58,  1.21it/s]
 72%|███████▏  | 358/500 [04:58<01:57,  1.21it/s]
 72%|███████▏  | 359/500 [04:58<01:56,  1.21it/s]
 72%|███████▏  | 360/500 [04:59<01:55,  1.21it/s]
 72%|███████▏  | 361/500 [05:00<01:54,  1.21it/s]
 72%|███████▏  | 362/500 [05:01<01:53,  1.21it/s]
 73%|███████▎  | 363/500 [05:02<01:53,  1.21it/s]
 73%|███████▎  | 364/500 [05:03<01:52,  1.21it/s]
 73%|███████▎  | 365/500 [05:03<01:51,  1.21it/s]
 73%|███████▎  | 366/500 [05:04<01:51,  1.20it/s]
 73%|███████▎  | 367/500 [05:05<01:50,  1.20it/s]
 74%|███████▎  | 368/500 [05:06<01:50,  1.20it/s]
 74%|███████▍  | 369/500 [05:07<01:49,  1.19it/s]
 74%|███████▍  | 370/500 [05:08<01:49,  1.19it/s]
 74%|███████▍  | 371/500 [05:08<01:49,  1.17it/s]
 74%|███████▍  | 372/500 [05:09<01:48,  1.18it/s]
 75%|███████▍  | 373/500 [05:10<01:47,  1.18it/s]
 75%|███████▍  | 374/500 [05:11<01:46,  1.18it/s]
 75%|███████▌  | 375/500 [05:12<01:46,  1.17it/s]
 75%|███████▌  | 376/500 [05:13<01:45,  1.17it/s]
 75%|███████▌  | 377/500 [05:14<01:45,  1.17it/s]
 76%|███████▌  | 378/500 [05:14<01:43,  1.17it/s]
 76%|███████▌  | 379/500 [05:15<01:42,  1.18it/s]
 76%|███████▌  | 380/500 [05:16<01:41,  1.18it/s]
 76%|███████▌  | 381/500 [05:17<01:40,  1.18it/s]
 76%|███████▋  | 382/500 [05:18<01:39,  1.18it/s]
 77%|███████▋  | 383/500 [05:19<01:39,  1.18it/s]
 77%|███████▋  | 384/500 [05:19<01:38,  1.18it/s]
 77%|███████▋  | 385/500 [05:20<01:37,  1.18it/s]
 77%|███████▋  | 386/500 [05:21<01:36,  1.18it/s]
 77%|███████▋  | 387/500 [05:22<01:35,  1.18it/s]
 78%|███████▊  | 388/500 [05:23<01:34,  1.18it/s]
 78%|███████▊  | 389/500 [05:24<01:34,  1.18it/s]
 78%|███████▊  | 390/500 [05:25<01:34,  1.16it/s]
 78%|███████▊  | 391/500 [05:25<01:33,  1.16it/s]
 78%|███████▊  | 392/500 [05:26<01:33,  1.16it/s]
 79%|███████▊  | 393/500 [05:27<01:32,  1.16it/s]
 79%|███████▉  | 394/500 [05:28<01:30,  1.17it/s]
 79%|███████▉  | 395/500 [05:29<01:29,  1.17it/s]
 79%|███████▉  | 396/500 [05:30<01:29,  1.16it/s]
 79%|███████▉  | 397/500 [05:31<01:28,  1.16it/s]
 80%|███████▉  | 398/500 [05:31<01:27,  1.17it/s]
 80%|███████▉  | 399/500 [05:32<01:26,  1.17it/s]
 80%|████████  | 400/500 [05:33<01:24,  1.18it/s]
 80%|████████  | 401/500 [05:34<01:23,  1.18it/s]
 80%|████████  | 402/500 [05:35<01:23,  1.18it/s]
 81%|████████  | 403/500 [05:36<01:22,  1.17it/s]
 81%|████████  | 404/500 [05:37<01:22,  1.17it/s]
 81%|████████  | 405/500 [05:37<01:21,  1.16it/s]
 81%|████████  | 406/500 [05:38<01:20,  1.17it/s]
 81%|████████▏ | 407/500 [05:39<01:19,  1.18it/s]
 82%|████████▏ | 408/500 [05:40<01:18,  1.18it/s]
 82%|████████▏ | 409/500 [05:41<01:16,  1.19it/s]
 82%|████████▏ | 410/500 [05:42<01:15,  1.19it/s]
 82%|████████▏ | 411/500 [05:42<01:14,  1.20it/s]
 82%|████████▏ | 412/500 [05:43<01:13,  1.20it/s]
 83%|████████▎ | 413/500 [05:44<01:12,  1.20it/s]
 83%|████████▎ | 414/500 [05:45<01:11,  1.20it/s]
 83%|████████▎ | 415/500 [05:46<01:10,  1.20it/s]
 83%|████████▎ | 416/500 [05:47<01:09,  1.21it/s]
 83%|████████▎ | 417/500 [05:47<01:08,  1.20it/s]
 84%|████████▎ | 418/500 [05:48<01:08,  1.20it/s]
 84%|████████▍ | 419/500 [05:49<01:07,  1.21it/s]
 84%|████████▍ | 420/500 [05:50<01:06,  1.20it/s]
 84%|████████▍ | 421/500 [05:51<01:05,  1.21it/s]
 84%|████████▍ | 422/500 [05:52<01:04,  1.21it/s]
 85%|████████▍ | 423/500 [05:52<01:04,  1.19it/s]
 85%|████████▍ | 424/500 [05:53<01:03,  1.20it/s]
 85%|████████▌ | 425/500 [05:54<01:02,  1.20it/s]
 85%|████████▌ | 426/500 [05:55<01:01,  1.20it/s]
 85%|████████▌ | 427/500 [05:56<01:00,  1.20it/s]
 86%|████████▌ | 428/500 [05:57<01:00,  1.19it/s]
 86%|████████▌ | 429/500 [05:57<00:59,  1.20it/s]
 86%|████████▌ | 430/500 [05:58<00:58,  1.20it/s]
 86%|████████▌ | 431/500 [05:59<00:57,  1.20it/s]
 86%|████████▋ | 432/500 [06:00<00:56,  1.20it/s]
 87%|████████▋ | 433/500 [06:01<00:56,  1.20it/s]
 87%|████████▋ | 434/500 [06:02<00:55,  1.20it/s]
 87%|████████▋ | 435/500 [06:02<00:54,  1.20it/s]
 87%|████████▋ | 436/500 [06:03<00:53,  1.20it/s]
 87%|████████▋ | 437/500 [06:04<00:52,  1.20it/s]
 88%|████████▊ | 438/500 [06:05<00:51,  1.20it/s]
 88%|████████▊ | 439/500 [06:06<00:50,  1.20it/s]
 88%|████████▊ | 440/500 [06:07<00:50,  1.20it/s]
 88%|████████▊ | 441/500 [06:07<00:49,  1.19it/s]
 88%|████████▊ | 442/500 [06:08<00:48,  1.19it/s]
 89%|████████▊ | 443/500 [06:09<00:47,  1.20it/s]
 89%|████████▉ | 444/500 [06:10<00:46,  1.19it/s]
 89%|████████▉ | 445/500 [06:11<00:45,  1.20it/s]
 89%|████████▉ | 446/500 [06:12<00:45,  1.20it/s]
 89%|████████▉ | 447/500 [06:12<00:44,  1.20it/s]
 90%|████████▉ | 448/500 [06:13<00:43,  1.20it/s]
 90%|████████▉ | 449/500 [06:14<00:42,  1.21it/s]
 90%|█████████ | 450/500 [06:15<00:41,  1.21it/s]
 90%|█████████ | 451/500 [06:16<00:40,  1.21it/s]
 90%|█████████ | 452/500 [06:17<00:39,  1.21it/s]
 91%|█████████ | 453/500 [06:17<00:38,  1.21it/s]
 91%|█████████ | 454/500 [06:18<00:38,  1.21it/s]
 91%|█████████ | 455/500 [06:19<00:37,  1.21it/s]
 91%|█████████ | 456/500 [06:20<00:36,  1.21it/s]
 91%|█████████▏| 457/500 [06:21<00:35,  1.21it/s]
 92%|█████████▏| 458/500 [06:22<00:34,  1.21it/s]
 92%|█████████▏| 459/500 [06:22<00:33,  1.22it/s]
 92%|█████████▏| 460/500 [06:23<00:33,  1.20it/s]
 92%|█████████▏| 461/500 [06:24<00:32,  1.21it/s]
 92%|█████████▏| 462/500 [06:25<00:31,  1.21it/s]
 93%|█████████▎| 463/500 [06:26<00:30,  1.21it/s]
 93%|█████████▎| 464/500 [06:26<00:29,  1.21it/s]
 93%|█████████▎| 465/500 [06:27<00:28,  1.21it/s]
 93%|█████████▎| 466/500 [06:28<00:28,  1.21it/s]
 93%|█████████▎| 467/500 [06:29<00:27,  1.21it/s]
 94%|█████████▎| 468/500 [06:30<00:26,  1.20it/s]
 94%|█████████▍| 469/500 [06:31<00:25,  1.20it/s]
 94%|█████████▍| 470/500 [06:31<00:24,  1.20it/s]
 94%|█████████▍| 471/500 [06:32<00:24,  1.20it/s]
 94%|█████████▍| 472/500 [06:33<00:23,  1.20it/s]
 95%|█████████▍| 473/500 [06:34<00:22,  1.20it/s]
 95%|█████████▍| 474/500 [06:35<00:21,  1.20it/s]
 95%|█████████▌| 475/500 [06:36<00:20,  1.20it/s]
 95%|█████████▌| 476/500 [06:36<00:19,  1.20it/s]
 95%|█████████▌| 477/500 [06:37<00:19,  1.20it/s]
 96%|█████████▌| 478/500 [06:38<00:18,  1.20it/s]
 96%|█████████▌| 479/500 [06:39<00:17,  1.20it/s]
 96%|█████████▌| 480/500 [06:40<00:16,  1.20it/s]
 96%|█████████▌| 481/500 [06:41<00:15,  1.20it/s]
 96%|█████████▋| 482/500 [06:41<00:15,  1.20it/s]
 97%|█████████▋| 483/500 [06:42<00:14,  1.20it/s]
 97%|█████████▋| 484/500 [06:43<00:13,  1.21it/s]
 97%|█████████▋| 485/500 [06:44<00:12,  1.21it/s]
 97%|█████████▋| 486/500 [06:45<00:11,  1.21it/s]
 97%|█████████▋| 487/500 [06:46<00:10,  1.21it/s]
 98%|█████████▊| 488/500 [06:46<00:09,  1.21it/s]
 98%|█████████▊| 489/500 [06:47<00:09,  1.21it/s]
 98%|█████████▊| 490/500 [06:48<00:08,  1.21it/s]
 98%|█████████▊| 491/500 [06:49<00:07,  1.21it/s]
 98%|█████████▊| 492/500 [06:50<00:06,  1.22it/s]
 99%|█████████▊| 493/500 [06:51<00:05,  1.22it/s]
 99%|█████████▉| 494/500 [06:51<00:04,  1.21it/s]
 99%|█████████▉| 495/500 [06:52<00:04,  1.21it/s]
 99%|█████████▉| 496/500 [06:53<00:03,  1.21it/s]
 99%|█████████▉| 497/500 [06:54<00:02,  1.21it/s]
100%|█████████▉| 498/500 [06:55<00:01,  1.21it/s]
100%|█████████▉| 499/500 [06:56<00:00,  1.21it/s]
100%|██████████| 500/500 [06:56<00:00,  1.20it/s]
100%|██████████| 500/500 [06:56<00:00,  1.20it/s]
The inversion result from `emcee`:
============================
Summary for inversion result
============================
SUCCESS
----------------------------
sampler: <emcee.ensemble.EnsembleSampler object>
blob_names: ['log_likelihood', 'log_prior']

sampler = inv_result.sampler
arviz.style.use("arviz-variat")
var_names = [
    "plate dip (°)",
]
az_idata = inv_result.to_arviz(var_names=var_names)
true_values = [60]
pc = arviz.plot_trace_dist(
    az_idata.sel(draw=slice(100, None)),
    visuals={"xlabel_trace": False, "trace": {"color": "C0"}, "dist": {"color": "C0"}},
    figure_kwargs={"figsize": (12, 4), "constrained_layout": True},
)
var_list = list(az_idata.posterior.data_vars)
for i, vname in enumerate(var_list):
    ax_kde = pc.iget_target(i, 0)
    ax_trace = pc.iget_target(i, 1)
    ax_kde.set_title(vname)
    ax_trace.set_title(vname)
    ax_kde.axvline(true_values[i], color="green", linestyle="--", lw=1, alpha=0.5)
    ax_trace.axhline(true_values[i], color="green", linestyle="--", lw=1, alpha=0.5)
    ax_trace.margins(x=0)

Watermark#

watermark_list = ["cofi", "numpy", "scipy", "matplotlib"]
for pkg in watermark_list:
    pkg_var = __import__(pkg)
    print(pkg, getattr(pkg_var, "__version__"))

cofi 0.2.11+71.gb28b5b0
numpy 2.2.6
scipy 1.17.1
matplotlib 3.10.9

sphinx_gallery_thumbnail_number = -1

Total running time of the script: (7 minutes 37.191 seconds)

Gallery generated by Sphinx-Gallery

Contents