.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "examples_gallery/plot_future_simulations.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_examples_gallery_plot_future_simulations.py>`
        to download the full example code.

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_examples_gallery_plot_future_simulations.py:


==================
Future simulations
==================

Shows how to Monte Carlo simulate the future.

.. GENERATED FROM PYTHON SOURCE LINES 9-114



.. image-sg:: /examples_gallery/images/sphx_glr_plot_future_simulations_001.png
   :alt: GYDA - Simulations
   :srcset: /examples_gallery/images/sphx_glr_plot_future_simulations_001.png
   :class: sphx-glr-single-img


.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    {'p': np.float64(1.1750182019164204), 'sigma': np.float64(0.2658888425595641), 'phi': np.float64(0.019681858945647165), 'theta': array([ 4.02310131,  3.95211774, -1.69933228])}






|

.. code-block:: Python


    import matplotlib.pyplot as plt
    import numpy as np
    import pandas as pd

    from dca.adca.well import Well
    from dca.datasets import load_monthly_sodir_production

    FIELD = "GYDA"

    # Read data from a public source (SODIR - sokkeldirektoratet)
    df = (
        load_monthly_sodir_production()
        .loc[lambda df: df.prfInformationCarrier == FIELD, :]
        .rename(
            columns={
                "prfInformationCarrier": "field",
                "prfYear": "year",
                "prfMonth": "month",
                "prfPrdOilNetMillSm3": "production",
            }
        )
        .assign(day=1)
        .assign(
            time=lambda df: pd.to_datetime(df[["year", "month", "day"]]).dt.to_period(
                freq="M"
            )
        )
        .loc[lambda df: df.production > 0]
    )


    well = Well(
        production=df.production.to_numpy(),
        time=df.time,
        id=FIELD,
        preprocessing="producing_time",
        curve_model="arps",
    )


    # Fit the well
    SPLIT = 0.3
    HALF_LIFE = 12  # The uncertainty is also very dependent on half-life
    well = well.fit(half_life=HALF_LIFE, prior_strength=1, split=SPLIT)
    print(well.parameters_)


    # Compute pre_eta, needed to initialize the simulations
    (x_train, y_train, w_train), (x_test, y_test, w_test) = well.get_train_test(split=SPLIT)
    preds_train = well.predict(x_train, time_on=w_train, q=0.5)
    prev_eta = (np.log(y_train) - np.log(preds_train))[-1]

    # Simulate a few production rates
    simulations = well.simulate(
        x_test,
        time_on=w_test,
        prev_eta=prev_eta,
        seed=abs(hash(well.id)) % 2**8,
        simulations=9,
    )


    fig, (ax1, ax2) = plt.subplots(figsize=(6, 3), nrows=2, ncols=1, sharex=True)

    # ===================== SIMULATE FUTURE PRODUCTION RATES =====================
    ax1.set_title(f"{FIELD} - Simulations")
    ax1.plot(x_train, y_train, color="black")
    ax1.plot(x_test, y_test, color="green")


    # Plot a single simulation (change to plot several)
    for j in range(3):
        ax1.plot(x_test, simulations[j, :], zorder=5, color="blue", alpha=0.2)

    ax1.set_ylabel("log(Prodrate)")
    ax1.grid(True, ls="--", alpha=0.33)

    # ===================== SIMULATE FUTURE PRODUCTION =====================

    # Simuate future cumulative production
    simulations = well.simulate_cumprod(
        x=x_test,
        prev_eta=prev_eta,
        seed=42,
        simulations=999,
    ) + np.sum(y_train * w_train)
    q10, q50, q90 = np.percentile(simulations, q=[10, 50, 90], axis=0)


    # Simulation results
    ax2.plot(x_test, q10, lw=1, ls="--", color="black")
    ax2.plot(x_test, q50, ls="-", color="black")
    ax2.plot(x_test, q90, lw=1, ls="--", color="black")

    # Train and test cumulatives
    ax2.plot(x_train, np.cumsum(y_train * w_train))
    ax2.plot(x_test, np.sum(y_train * w_train) + np.cumsum(y_test * w_test))


    ax2.set_ylabel("Production")
    ax2.set_xlabel("Months since production start")
    ax2.grid(True, ls="--", alpha=0.33)
    fig.tight_layout()
    plt.show()


.. rst-class:: sphx-glr-timing

   **Total running time of the script:** (0 minutes 0.474 seconds)


.. _sphx_glr_download_examples_gallery_plot_future_simulations.py:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example

    .. container:: sphx-glr-download sphx-glr-download-jupyter

      :download:`Download Jupyter notebook: plot_future_simulations.ipynb <plot_future_simulations.ipynb>`

    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download Python source code: plot_future_simulations.py <plot_future_simulations.py>`

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: plot_future_simulations.zip <plot_future_simulations.zip>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_