.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/basicreservoir.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_basicreservoir.py: 2D Reservoir ESMDA example ========================== Ensemble Smoother Multiple Data Assimilation (ES-MDA) in Reservoir Simulation. .. GENERATED FROM PYTHON SOURCE LINES 8-18 .. code-block:: Python import numpy as np import matplotlib.pyplot as plt import resmda # For reproducibility, we instantiate a random number generator with a fixed # seed. For production, remove the seed! rng = np.random.default_rng(1848) .. GENERATED FROM PYTHON SOURCE LINES 20-22 Model parameters ---------------- .. GENERATED FROM PYTHON SOURCE LINES 22-49 .. code-block:: Python # Grid extension nx = 30 ny = 25 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)] nt = time.size # Assumed sandard deviation of our data dstd = 0.5 # Observation location indices (should be well locations) ox, oy = 1, 1 # Wells (if None, first and last cells are used with pressure 180 and 120) wells = None .. GENERATED FROM PYTHON SOURCE LINES 50-56 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 56-82 .. code-block:: Python # Get the model and ne prior models RP = resmda.RandomPermeability(nx, ny, perm_mean, perm_min, perm_max) perm_true = RP(1, random=rng) perm_prior = RP(ne, random=rng) # QC covariance, reference model, and first two random models pinp1 = {"origin": "lower", "vmin": perm_min, "vmax": perm_max} fig, axs = plt.subplots(2, 2, figsize=(6, 6), constrained_layout=True) axs[0, 0].set_title("Model") im = axs[0, 0].imshow(perm_true.T, **pinp1) axs[0, 1].set_title("Lower Covariance Matrix") im2 = axs[0, 1].imshow(RP.cov, cmap="plasma") axs[1, 0].set_title("Random Model 1") axs[1, 0].imshow(perm_prior[0, ...].T, **pinp1) axs[1, 1].set_title("Random Model 2") axs[1, 1].imshow(perm_prior[1, ...].T, **pinp1) fig.colorbar(im, ax=axs[1, :], orientation="horizontal", label="Log of Permeability (mD)") for ax in axs[1, :].ravel(): ax.set_xlabel("x-direction") for ax in axs[:, 0].ravel(): ax.set_ylabel("y-direction") .. image-sg:: /gallery/images/sphx_glr_basicreservoir_001.png :alt: Model, Lower Covariance Matrix, Random Model 1, Random Model 2 :srcset: /gallery/images/sphx_glr_basicreservoir_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 83-85 Run the prior models and the reference case ------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 85-111 .. code-block:: Python # Instantiate reservoir simulator RS = resmda.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)) # 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) # QC data and priors fig, ax = plt.subplots(1, 1, constrained_layout=True) ax.set_title("Observed and prior data") ax.plot(time*24*60*60, data_prior.T, color=".6", alpha=0.5) ax.plot(time*24*60*60, data_true, "ko", label="True data") ax.plot(time*24*60*60, data_obs, "C3o", label="Obs. data") ax.legend() ax.set_xlabel("Time (s)") ax.set_ylabel("Pressure (bar)") .. image-sg:: /gallery/images/sphx_glr_basicreservoir_002.png :alt: Observed and prior data :srcset: /gallery/images/sphx_glr_basicreservoir_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 112-114 ESMDA ----- .. GENERATED FROM PYTHON SOURCE LINES 114-133 .. code-block:: Python def restrict_permeability(x): """Restrict possible permeabilities.""" np.clip(x, perm_min, perm_max, out=x) perm_post, data_post = resmda.esmda( model_prior=perm_prior, forward=sim, data_obs=data_obs, sigma=dstd, alphas=4, data_prior=data_prior, callback_post=restrict_permeability, random=rng, ) .. rst-class:: sphx-glr-script-out .. code-block:: none ES-MDA step 1; α=4.0 ES-MDA step 2; α=4.0 ES-MDA step 3; α=4.0 ES-MDA step 4; α=4.0 .. GENERATED FROM PYTHON SOURCE LINES 134-140 Posterior Analysis ------------------ After running ESMDA, it's crucial to analyze the posterior ensemble of models. Here, we visualize the first three realizations from both the prior and posterior ensembles to see how the models have been updated. .. GENERATED FROM PYTHON SOURCE LINES 140-152 .. code-block:: Python # Plot posterior fig, ax = plt.subplots(1, 2, figsize=(8, 5), constrained_layout=True) pinp2 = {"origin": "lower", "vmin": 2.5, "vmax": 3.5} ax[0].set_title("Prior Mean") im = ax[0].imshow(perm_prior.mean(axis=0).T, **pinp2) ax[1].set_title("Post Mean") ax[1].imshow(perm_post.mean(axis=0).T, **pinp2) fig.colorbar(im, ax=ax, label="Log of Permeability (mD)", orientation="horizontal") .. image-sg:: /gallery/images/sphx_glr_basicreservoir_003.png :alt: Prior Mean, Post Mean :srcset: /gallery/images/sphx_glr_basicreservoir_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 153-158 Observing the monitored pressure at cell (1,1) for all realizations and the reference case, we can see that the ensemble of models after the assimilation steps (in blue) is closer to the reference case (in red) than the prior ensemble (in gray). This indicates that the ESMDA method is effectively updating the models to better represent the observed data. .. GENERATED FROM PYTHON SOURCE LINES 158-171 .. code-block:: Python # Compare posterior to prior and observed data fig, ax = plt.subplots(1, 1, constrained_layout=True) ax.set_title("Prior and posterior data") ax.plot(time*24*60*60, data_prior.T, color=".6", alpha=0.5) ax.plot(time*24*60*60, data_post.T, color="C0", alpha=0.5) ax.plot(time*24*60*60, data_true, "ko") ax.plot(time*24*60*60, data_obs, "C3o") ax.set_xlabel("Time (s)") ax.set_ylabel("Pressure (bar)") .. image-sg:: /gallery/images/sphx_glr_basicreservoir_004.png :alt: Prior and posterior data :srcset: /gallery/images/sphx_glr_basicreservoir_004.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 172-174 .. code-block:: Python resmda.Report() .. raw:: html
Wed Jul 24 14:04:06 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.11.9 (main, Jul 15 2024, 21:50:21) [GCC 11.4.0]
numpy 2.0.1 scipy 1.14.0 numba Module not found
resmda 0.2.0 matplotlib 3.9.1 IPython 8.26.0


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