.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/localization.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_localization.py: Localization ========================== Example using several well doublets and compares ESMDA with and without localization. The example follows contextually :ref:`sphx_glr_gallery_basicreservoir.py`. .. GENERATED FROM PYTHON SOURCE LINES 11-21 .. code-block:: Python import numpy as np import matplotlib.pyplot as plt import dageo # For reproducibility, we instantiate a random number generator with a fixed # seed. For production, remove the seed! rng = np.random.default_rng(2020) .. GENERATED FROM PYTHON SOURCE LINES 23-25 Model parameters ---------------- .. GENERATED FROM PYTHON SOURCE LINES 25-59 .. code-block:: Python # Grid extension nx = 35 ny = 30 nc = nx*ny # Permeabilities perm_mean = 3.0 perm_min = 0.5 perm_max = 5.0 # ESMDA parameters ne = 100 # Number of ensembles dt = np.zeros(10)+0.0001 # Time steps (could be irregular, e.g., increasing!) time = np.r_[0, np.cumsum(dt)] # Assumed sandard deviation of our data dstd = 0.5 # Measurement points ox = (5, 15, 24) oy = (5, 10, 24) # Number of data points nd = time.size * len(ox) # Wells wells = np.array([ [ox[0], oy[0], 180], [5, 12, 120], [ox[1], oy[1], 180], [20, 5, 120], [ox[2], oy[2], 180], [24, 17, 120] ]) .. GENERATED FROM PYTHON SOURCE LINES 60-66 Create permeability maps for ESMDA ---------------------------------- We will create a set of permeability maps that will serve as our initial guess (prior). These maps are generated using a Gaussian random field and are constrained by certain statistical properties. .. GENERATED FROM PYTHON SOURCE LINES 66-73 .. code-block:: Python # Get the model and ne prior models RP = dageo.RandomPermeability(nx, ny, perm_mean, perm_min, perm_max) perm_true = RP(1, random=rng) perm_prior = RP(ne, random=rng) .. GENERATED FROM PYTHON SOURCE LINES 74-76 Run the prior models and the reference case ------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 76-93 .. code-block:: Python # Instantiate reservoir simulator RS = dageo.Simulator(nx, ny, wells=wells) def sim(x): """Custom fct to use exp(x), and specific dt & location.""" return RS(np.exp(x), dt=dt, data=(ox, oy)).reshape((x.shape[0], -1)) # Simulate data for the prior and true fields data_prior = sim(perm_prior) data_true = sim(perm_true) data_obs = rng.normal(data_true, dstd) data_obs[0, :3] = data_true[0, :3] .. GENERATED FROM PYTHON SOURCE LINES 94-96 Localization Matrix ------------------- .. GENERATED FROM PYTHON SOURCE LINES 96-114 .. code-block:: Python # Vector of nd length with the well x and y position for each nd data point nd_positions = np.tile(np.array([ox, oy]), time.size).T # Create matrix loc_mat = dageo.localization_matrix(RP.cov, nd_positions, (nx, ny)) # QC localization matrix fig, ax = plt.subplots(1, 1, constrained_layout=True) im = ax.imshow(loc_mat.sum(axis=2).T/time.size, origin='lower', cmap="plasma") ax.contour(loc_mat.sum(axis=2).T/time.size, levels=[0.2, ], colors='w') ax.set_xlabel('x-direction') ax.set_ylabel('y-direction') for well in wells: ax.plot(well[0], well[1], ['C3v', 'C1^'][int(well[2] == 120)]) fig.colorbar(im, ax=ax, label="Localization Weight (-)") .. image-sg:: /gallery/images/sphx_glr_localization_001.png :alt: localization :srcset: /gallery/images/sphx_glr_localization_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 115-117 ESMDA ----- .. GENERATED FROM PYTHON SOURCE LINES 117-136 .. code-block:: Python def restrict_permeability(x): """Restrict possible permeabilities.""" np.clip(x, perm_min, perm_max, out=x) inp = { 'model_prior': perm_prior, 'forward': sim, 'data_obs': data_obs, 'sigma': dstd, 'alphas': 4, 'data_prior': data_prior, 'callback_post': restrict_permeability, 'random': rng, } .. GENERATED FROM PYTHON SOURCE LINES 137-139 Without localization '''''''''''''''''''' .. GENERATED FROM PYTHON SOURCE LINES 139-143 .. code-block:: Python nl_perm_post, nl_data_post = dageo.esmda(**inp) .. rst-class:: sphx-glr-script-out .. code-block:: none ESMDA step 1; α=4.0 ESMDA step 2; α=4.0 ESMDA step 3; α=4.0 ESMDA step 4; α=4.0 .. GENERATED FROM PYTHON SOURCE LINES 144-146 With localization ''''''''''''''''' .. GENERATED FROM PYTHON SOURCE LINES 146-150 .. code-block:: Python wl_perm_post, wl_data_post = dageo.esmda(**inp, localization_matrix=loc_mat) .. rst-class:: sphx-glr-script-out .. code-block:: none ESMDA step 1; α=4.0 ESMDA step 2; α=4.0 ESMDA step 3; α=4.0 ESMDA step 4; α=4.0 .. GENERATED FROM PYTHON SOURCE LINES 151-153 Compare permeabilities ---------------------- .. GENERATED FROM PYTHON SOURCE LINES 153-182 .. code-block:: Python # Plot posterior fig, axs = plt.subplots( 1, 3, figsize=(8, 4), sharex=True, sharey=True, constrained_layout=True) par = {"vmin": perm_min, "vmax": perm_max, "origin": "lower"} axs[0].set_title("Prior Mean") im = axs[0].imshow(perm_prior.mean(axis=0).T, **par) axs[1].set_title("Post Mean; No localization") axs[1].imshow(nl_perm_post.mean(axis=0).T, **par) axs[2].set_title("Post Mean: Localization") axs[2].imshow(wl_perm_post.mean(axis=0).T, **par) axs[2].contour(loc_mat.sum(axis=2).T, levels=[2.0, ], colors="w") fig.colorbar(im, ax=axs, label="Log Permeabilities (mD)", orientation="horizontal") for ax in axs: for well in wells: ax.plot(well[0], well[1], ["C3v", "C1^"][int(well[2] == 120)]) for ax in axs: ax.set_xlabel('x-direction') axs[0].set_ylabel('y-direction') .. image-sg:: /gallery/images/sphx_glr_localization_002.png :alt: Prior Mean, Post Mean; No localization, Post Mean: Localization :srcset: /gallery/images/sphx_glr_localization_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 183-185 Compare data ------------ .. GENERATED FROM PYTHON SOURCE LINES 185-207 .. code-block:: Python # QC data and priors fig, axs = plt.subplots( 2, 3, figsize=(8, 5), sharex=True, sharey=True, constrained_layout=True) fig.suptitle('Prior and posterior data') for i, ax in enumerate(axs[0, :]): ax.set_title(f"Well ({ox[i]}; {oy[i]})") ax.plot(time*24*60*60, data_prior[..., i::3].T, color='.6', alpha=0.5) ax.plot(time*24*60*60, nl_data_post[..., i::3].T, color='C0', alpha=0.5) ax.plot(time*24*60*60, data_obs[0, i::3], 'C3o') for i, ax in enumerate(axs[1, :]): ax.plot(time*24*60*60, data_prior[..., i::3].T, color='.6', alpha=0.5) ax.plot(time*24*60*60, wl_data_post[..., i::3].T, color='C0', alpha=0.5) ax.plot(time*24*60*60, data_obs[0, i::3], 'C3o') ax.set_xlabel("Time (s)") for i, ax in enumerate(axs[:, 0]): ax.set_ylabel("Pressure (bar)") for i, txt in enumerate(["No l", "L"]): axs[i, 2].yaxis.set_label_position("right") axs[i, 2].set_ylabel(f"{txt}ocalization") .. image-sg:: /gallery/images/sphx_glr_localization_003.png :alt: Prior and posterior data, Well (5; 5), Well (15; 10), Well (24; 24) :srcset: /gallery/images/sphx_glr_localization_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 208-210 .. code-block:: Python dageo.Report() .. raw:: html
Fri Dec 06 15:59:20 2024 UTC
OS Linux (Ubuntu 22.04) CPU(s) 4 Machine x86_64
Architecture 64bit RAM 15.6 GiB Environment Python
File system ext4
Python 3.12.7 (main, Oct 1 2024, 15:17:55) [GCC 11.4.0]
dageo 1.1.0 numpy 2.1.3 scipy 1.14.1
matplotlib 3.9.3 IPython 8.30.0


.. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 20.389 seconds) **Estimated memory usage:** 129 MB .. _sphx_glr_download_gallery_localization.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: localization.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: localization.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: localization.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_