Stochastic doublet

Stochastic Doublet Performance Calculations

This is a simple example of how to perform stochastic doublet simulations. In a stochastic simulation, multiple p-values of transmissivity, thickness and permeability are used. Transmissivity is defined through the combination of a Normal distribution of thickness and a Log-Normal distribution of permeability. For a predetermined set of p-values simulations are run across these transmissivity samples resulting in a output Dataset with the dimension p_value. This is how ThermoGIS currently handles uncertainty in the reservoir properties.

xarray is at the heart of data handling in pythermogis

Xarray is a very powerfull tool. Pythermogis is build using xarray, and leverages its powers. For basic 1D usage, no xarray knowledge is required. For 2D the Aquifer class expects xarray DataArrays. For more information on how to use xarray, see the xarray documentation.

🧪 Basic 1D Example

This example demonstrates how to run a stochastic doublet simulation using the StochasticAquifer and ThermoGISDoublet classes. For this example only the P50 values are calculated. The outcomes are stochastic, meaning that different output values can occur with each run.

Permeability vs Log-permeability

The StochasticAquifer requires permeability in natural space, while the stochastic method defines permeability as a Log-normal distribution with a mean and a standard deviation. This means that permeability mean and sd is in log space. A log-permeability of 5 is equal to 148.6mD natural-permability (for example).

from pythermogis.aquifer import StochasticAquifer
from pythermogis.doublet import ThermoGISDoublet

aquifer = StochasticAquifer(
    thickness_mean=300,
    thickness_sd=50,
    ntg=0.5,
    porosity=0.2,
    depth=2000,
    ln_permeability_mean=4,
    ln_permeability_sd=0.5,
)
results = ThermoGISDoublet(aquifer=aquifer).simulate(p_values=[50]).to_dataset()

print(results)

🧪 P-values Example

This example is an extension of the deterministic doublet simulation using the StochasticAquifer class. The outcomes correspond to the p-value range of the transmissivity are presented as an expectation curve of transmissivity. This example corresponds to test case test_example5 in test_doc_examples.py in the tests directory of the repository.

from pythermogis.aquifer import StochasticAquifer
from pythermogis.doublet import ThermoGISDoublet

from matplotlib import pyplot as plt
from pathlib import Path
output_data_path = Path(__file__).parent.parent / "resources" / "test_output" / "example5"
output_data_path.mkdir(parents=True, exist_ok=True)

aquifer = StochasticAquifer(
    thickness_mean=300,
    thickness_sd=50,
    ntg=0.5,
    porosity=0.2,
    depth=2000,
    ln_permeability_mean=5,
    ln_permeability_sd=0.5,
)
doublet = ThermoGISDoublet(aquifer=aquifer)
results = doublet.simulate(
    p_values=[1, 10, 20, 30, 40, 50, 60, 70, 80, 90, 99]
).to_dataset()

fig, axe = plt.subplots(figsize=(10, 5))
kh = results.transmissivity_with_ntg * 1.0
kh.plot(y="p_value", ax=axe)
axe.set_title(f"Aquifer Net transmissivity\n kH")
temp =  aquifer.temperature
inj_temp = results.inj_temp.sel(p_value=50, method="nearest")
prd_temp = results.prd_temp.sel(p_value=50, method="nearest")
axe.axhline(50.0, label=f"TEMP res: {temp:.1f} inj: {inj_temp:.1f} prd: {prd_temp:.1f} ",
            ls="--", c="tab:orange")
# axes.axhline(pos, label=f"probability of success: {pos:.1f}%",ls="--", c="tab:orange")
kh50 = kh.sel(p_value=50, method="nearest")
axe.axvline(kh50, ls="--", c="tab:orange" )
axe.legend()
axe.set_xlabel("kH [Dm]")
axe.set_ylabel("p-value [%]")
plt.savefig(output_data_path / "kh.png")

The other examples from deterministic doublet can easily be adapted to this stochastic framework.