Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fluvial example #11

Merged
merged 14 commits into from
Jul 11, 2024
3 changes: 3 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -2,12 +2,15 @@
__pycache__/
*.nc
*.pdf
*.zip

# Sphinx
docs/_build/
docs/api/resmda*
docs/savefig/
docs/gallery/*
download/
docs/sg_execution_times.rst

# Pytest and coverage related
htmlcov
Expand Down
7 changes: 3 additions & 4 deletions docs/conf.py
Original file line number Diff line number Diff line change
@@ -1,7 +1,6 @@
import time
import warnings
from resmda import __version__
from sphinx_gallery.sorting import FileNameSortKey

# ==== 1. Extensions ====

Expand Down Expand Up @@ -39,17 +38,17 @@
sphinx_gallery_conf = {
'examples_dirs': ['../examples', ],
'gallery_dirs': ['gallery', ],
'capture_repr': ('_repr_html_', '__repr__'),
'capture_repr': ('_repr_html_', ),
# Patter to search for example files
"filename_pattern": r"\.py",
# Sort gallery example by file name instead of number of lines (default)
"within_subsection_order": FileNameSortKey,
"within_subsection_order": "FileNameSortKey",
# Remove the settings (e.g., sphinx_gallery_thumbnail_number)
'remove_config_comments': True,
# Show memory
'show_memory': True,
# Custom first notebook cell
'first_notebook_cell': '%matplotlib notebook',
'first_notebook_cell': '%matplotlib widget',
}

# https://github.com/sphinx-gallery/sphinx-gallery/pull/521/files
Expand Down
7 changes: 3 additions & 4 deletions examples/basicESMDA.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,8 +19,8 @@
.. math::
x_{j,i+1} &= x_{j,i} + (C^e_{xy})_i \left((C^e_{yy})_i +
\alpha_iC^e_{dd}\right)^{-1} \left(d + \sqrt{\alpha_i}
\varepsilon_j - g(x_{j,i})\right) \\
y_{j,i+1} &= g(x_{j,i+1})
\varepsilon_j - g(x_{j,i})\right) \ , \\
y_{j,i+1} &= g(x_{j,i+1}) \ .

The model used for this example is

Expand Down Expand Up @@ -53,7 +53,7 @@ def forward(x, beta):
px = np.linspace(-5, 5, 301)
for i, b in enumerate([0.0, 0.2]):
axs[i].set_title(
f"{['Linear model', 'Nonlinear model'][i]}: $\\beta$ = {b}")
f"{['Linear model', 'Nonlinear model'][i]}: β = {b}")
axs[i].plot(px, forward(px, b))
axs[i].set_xlabel('x')
axs[i].set_ylabel('y')
Expand Down Expand Up @@ -109,7 +109,6 @@ def plot_result(mpost, dpost, dobs, title, ylim):
ax2.set_xlabel('y')
ax2.set_xlim([-5, 8])
ax2.legend()
fig.show()


###############################################################################
Expand Down
85 changes: 40 additions & 45 deletions examples/basicreservoir.py
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@
# ----------------

# Grid extension
nx = 25
nx = 30
ny = 25
nc = nx*ny

Expand All @@ -32,7 +32,6 @@

# ESMDA parameters
ne = 100 # Number of ensembles
# dt = np.logspace(-5, -3, 10)
dt = np.zeros(10)+0.0001 # Time steps (could be irregular, e.g., increasing!)
time = np.r_[0, np.cumsum(dt)]
nt = time.size
Expand Down Expand Up @@ -61,25 +60,23 @@
perm_prior = RP(ne, random=rng)


# TODO: change scale in imshow to represent meters

# QC covariance, reference model, and first two random models
fig, axs = plt.subplots(2, 2, constrained_layout=True)
axs[0, 0].set_title('Model')
im = axs[0, 0].imshow(perm_true.T, vmin=perm_min, vmax=perm_max)
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, vmin=perm_min, vmax=perm_max)
axs[1, 1].set_title('Random Model 2')
axs[1, 1].imshow(perm_prior[1, ...].T, vmin=perm_min, vmax=perm_max)
fig.colorbar(im, ax=axs[1, :], orientation='horizontal',
label='Log of Permeability (mD)')
for ax in axs[1, :]:
ax.set_xlabel('x-direction')
for ax in axs[:, 0]:
ax.set_ylabel('y-direction')
fig.show()
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")


###############################################################################
Expand All @@ -101,15 +98,14 @@ def sim(x):
data_obs = rng.normal(data_true, dstd)

# QC data and priors
fig, ax = plt.subplots(1, 1)
ax.set_title('Observed and prior data')
ax.plot(time, data_prior.T, color='.6', alpha=0.5)
ax.plot(time, data_true, 'ko', label='True data')
ax.plot(time, data_obs, 'C3o', label='Obs. data')
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 (???)')
ax.set_ylabel('Pressure (???)')
fig.show()
ax.set_xlabel("Time (s)")
ax.set_ylabel("Pressure (bar)")


###############################################################################
Expand Down Expand Up @@ -143,14 +139,14 @@ def restrict_permeability(x):
# and posterior ensembles to see how the models have been updated.

# Plot posterior
fig, ax = plt.subplots(1, 3, figsize=(12, 5))
ax[0].set_title('Prior Mean')
ax[0].imshow(perm_prior.mean(axis=0).T)
ax[1].set_title('Post Mean')
ax[1].imshow(perm_post.mean(axis=0).T)
ax[2].set_title('"Truth"')
ax[2].imshow(perm_true.T)
fig.show()
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")


###############################################################################
Expand All @@ -162,15 +158,14 @@ def restrict_permeability(x):


# Compare posterior to prior and observed data
fig, ax = plt.subplots(1, 1)
ax.set_title('Prior and posterior data')
ax.plot(time, data_prior.T, color='.6', alpha=0.5)
ax.plot(time, data_post.T, color='C0', alpha=0.5)
ax.plot(time, data_true, 'ko')
ax.plot(time, data_obs, 'C3o')
ax.set_xlabel('Time (???)')
ax.set_ylabel('Pressure (???)')
fig.show()
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)")


###############################################################################
Expand Down
Loading