From 29c4a359cf50d2b43de869e5bea31992fbf9711c Mon Sep 17 00:00:00 2001 From: antoinecollet5 Date: Tue, 3 Oct 2023 23:31:02 +0200 Subject: [PATCH] DOCS: improved documentation using citet and citep for references, makefiles for unix and windows to avoid having to use tox, use pydata theme rather than read the docs and improve docstrings for classes SIES and ES. Add badges for project visibility. --- .gitignore | 3 + README.md | 35 ++-- docs/source/_templates/autosummary/base.rst | 8 + docs/source/_templates/autosummary/class.rst | 18 ++ docs/source/_templates/autosummary/module.rst | 64 +++++++ docs/source/api.rst | 13 +- docs/source/bibliography.rst | 2 + docs/source/conf.py | 160 +++++++++++++++++- docs/source/examples.rst | 2 + docs/source/glossary.rst | 4 +- docs/source/index.rst | 26 +-- docs/source/refs.bib | 15 ++ pyproject.toml | 11 +- src/iterative_ensemble_smoother/__init__.py | 32 +++- src/iterative_ensemble_smoother/esmda.py | 7 +- .../esmda_inversion.py | 23 ++- src/iterative_ensemble_smoother/sies.py | 101 +++++------ src/iterative_ensemble_smoother/utils.py | 11 +- 18 files changed, 394 insertions(+), 141 deletions(-) create mode 100644 docs/source/_templates/autosummary/base.rst create mode 100644 docs/source/_templates/autosummary/class.rst create mode 100644 docs/source/_templates/autosummary/module.rst diff --git a/.gitignore b/.gitignore index 8b9dee82..bc1db546 100644 --- a/.gitignore +++ b/.gitignore @@ -176,4 +176,7 @@ poetry.toml # LSP config files pyrightconfig.json +# For the documentation +_autosummary + # End of https://www.toptal.com/developers/gitignore/api/python diff --git a/README.md b/README.md index e066cec8..554b13c2 100644 --- a/README.md +++ b/README.md @@ -1,8 +1,15 @@ Iterative Ensemble Smoother =========================== +[![License: GPL v3](https://img.shields.io/badge/License-GPLv3-blue.svg)](https://github.com/equinor/iterative_ensemble_smoother/blob/main/COPYING) +[![Stars](https://img.shields.io/github/stars/equinor/iterative_ensemble_smoother.svg?style=social&label=Star&maxAge=2592000)](https://github.com/equinor/iterative_ensemble_smoother/stargazers) +[![Python](https://img.shields.io/pypi/pyversions/iterative_ensemble_smoother.svg)](https://pypi.org/pypi/iterative_ensemble_smoother) +[![PyPI](https://img.shields.io/pypi/v/iterative_ensemble_smoother.svg)](https://pypi.org/pypi/iterative_ensemble_smoother) +[![Downloads](https://static.pepy.tech/badge/iterative_ensemble_smoother)](https://pepy.tech/project/iterative_ensemble_smoother) +[![Build Status](https://github.com/equinor/iterative_ensemble_smoother/actions/workflows/upload_to_pypi.yml/badge.svg)](https://github.com/equinor/iterative_ensemble_smoother/actions/workflows/main.yml) [![Precommit: enabled](https://img.shields.io/badge/pre--commit-enabled-brightgreen?logo=pre-commit)](https://github.com/pre-commit/pre-commit) [![Ruff](https://img.shields.io/endpoint?url=https://raw.githubusercontent.com/astral-sh/ruff/main/assets/badge/v2.json)](https://github.com/astral-sh/ruff) +[![Mypy](https://www.mypy-lang.org/static/mypy_badge.svg)](https://mypy-lang.org/) [![Code style: black](https://img.shields.io/badge/code%20style-black-000000.svg)](https://github.com/psf/black) [![docs](https://readthedocs.org/projects/iterative_ensemble_smoother/badge/?version=latest&style=plastic)](https://iterative-ensemble-smoother.readthedocs.io/) @@ -21,31 +28,9 @@ pip install iterative_ensemble_smoother ## Usage -The following illustrates usage but does not actually create the inputs and can not be run. +**iterative_ensemble_smoother** mainly implements `SIES` and `ESMDA` classes. Check out +the examples section to see how to use them. -The example below illustrates the usage of the package. -Note that it does not create the necessary inputs and cannot be run directly. -For complete information and runnable examples, please refer to [the documentation](http://iterative_ensemble_smoother.rtfd.io). - -```python -import iterative_ensemble_smoother as ies - -# `prior` is a matrix of size (num_params, ensemble_size). -# In geostatistical applications, it typically consists of Gaussian random fields -# that model rock properties like porosities and permeabilities. - -# `response_ensemble` is a matrix of size (num_obs, ensemble_size) and is -# generated by running a dynamical model, such as a flow simulator. -# In geostatistical applications, this is typically a porous flow simulator like Eclipse or OPM flow. - -# `obs_errors` and `obs_values` are 1D array of size `num_obs` that hold observations from real-life -# sensors and their uncertainty specifications. -# In geostatistical applications, the observations are typically pressures, temperatures, production rates etc. - -smoother = ies.ES() -smoother.fit(response_ensemble, obs_errors, obs_values) -posterior = smoother.update(prior) -``` ## Building from source @@ -57,7 +42,7 @@ cd iterative_ensemble_smoother pip install . ``` -### Building the documentation +## Building the documentation ```bash apt install pandoc # Pandoc is required to build the documentation. diff --git a/docs/source/_templates/autosummary/base.rst b/docs/source/_templates/autosummary/base.rst new file mode 100644 index 00000000..d1e5d94f --- /dev/null +++ b/docs/source/_templates/autosummary/base.rst @@ -0,0 +1,8 @@ +.. + Template for the html base rendering + +{{ fullname | escape | underline}} + +.. currentmodule:: {{ module }} + +.. auto{{ objtype }}:: {{ objname }} diff --git a/docs/source/_templates/autosummary/class.rst b/docs/source/_templates/autosummary/class.rst new file mode 100644 index 00000000..96adfa2d --- /dev/null +++ b/docs/source/_templates/autosummary/class.rst @@ -0,0 +1,18 @@ +.. + Template for the html class rendering + +{{ fullname | escape | underline}} + +.. currentmodule:: {{ module }} + +.. autoclass:: {{ objname }} + :members: + :private-members: + :show-inheritance: + :inherited-members: + + {% block methods %} + .. automethod:: __init__ + {% endblock %} + + .. rubric:: {{ _('Methods definition') }} diff --git a/docs/source/_templates/autosummary/module.rst b/docs/source/_templates/autosummary/module.rst new file mode 100644 index 00000000..47fa73ec --- /dev/null +++ b/docs/source/_templates/autosummary/module.rst @@ -0,0 +1,64 @@ +.. + Template for the html module rendering + +{{ fullname | escape | underline}} + +.. automodule:: {{ fullname }} + + {% block attributes %} + {% if attributes %} + .. rubric:: {{ _('Module Attributes') }} + + .. autosummary:: + {% for item in attributes %} + {{ item }} + {%- endfor %} + {% endif %} + {% endblock %} + + {% block functions %} + {% if functions %} + .. rubric:: {{ _('Functions') }} + + .. autosummary:: + {% for item in functions %} + {{ item }} + {%- endfor %} + {% endif %} + {% endblock %} + + {% block classes %} + {% if classes %} + .. rubric:: {{ _('Classes') }} + + .. autosummary:: + {% for item in classes %} + :template: class.rst + {{ item }} + {%- endfor %} + {% endif %} + {% endblock %} + + {% block exceptions %} + {% if exceptions %} + .. rubric:: {{ _('Exceptions') }} + + .. autosummary:: + {% for item in exceptions %} + {{ item }} + {%- endfor %} + {% endif %} + {% endblock %} + +{% block modules %} +{% if modules %} +.. rubric:: Modules + +.. autosummary:: + :toctree: + :recursive: +{% for item in modules %} + {{ item }} +{%- endfor %} +{% endif %} +{% endblock %} diff --git a/docs/source/api.rst b/docs/source/api.rst index a00951cb..c2896973 100644 --- a/docs/source/api.rst +++ b/docs/source/api.rst @@ -1,10 +1,11 @@ -API reference +.. _api_ref: + +============= +API Reference ============= -.. currentmodule:: iterative_ensemble_smoother +.. automodule:: iterative_ensemble_smoother -.. autoclass:: iterative_ensemble_smoother.SIES - :members: +.. raw:: latex -.. autoclass:: iterative_ensemble_smoother.ESMDA - :members: + \clearpage diff --git a/docs/source/bibliography.rst b/docs/source/bibliography.rst index 0d21440d..eb35d1fa 100644 --- a/docs/source/bibliography.rst +++ b/docs/source/bibliography.rst @@ -2,3 +2,5 @@ Bibliography ============ .. bibliography:: + :all: + :labelprefix: Bib- diff --git a/docs/source/conf.py b/docs/source/conf.py index fd4b3d99..53b5801d 100644 --- a/docs/source/conf.py +++ b/docs/source/conf.py @@ -1,26 +1,168 @@ +""" +iterative_ensemble_smoother documentation build configuration file. +""" +import datetime +import os +import re +import sys from subprocess import check_output +import iterative_ensemble_smoother as ies + +package_path = os.path.abspath("..") +sys.path.insert(0, package_path) + project = "iterative_ensemble_smoother" -copyright = "2022, Equinor" author = "Equinor" -release = "0.1.1" - +copyright = f"2022-{datetime.datetime.today().year}, {author}" +# The default replacements for |version| and |release|, also used in various +# other places throughout the built documents. +version = re.sub(r"\.dev.*$", r".dev", ies.__version__) +release = version +# convert the python file to a notebook check_output(["jupytext", "Polynomial.py", "-o", "Polynomial.ipynb"]) +# Do the same for this file check_output(["jupytext", "Oscillator.py", "-o", "Oscillator.ipynb"]) - +# Add any Sphinx extension module names here, as strings. They can be +# extensions coming with Sphinx (named 'sphinx.ext.*') or your custom ones. extensions = [ - "sphinx.ext.autodoc", - "sphinx.ext.doctest", - "nbsphinx", + "sphinx.ext.napoleon", # autodoc understands numpy docstrings + "sphinx.ext.autodoc", # Core library for html generation from docstrings + "sphinx.ext.autosummary", # Create neat summary tables + "sphinx.ext.doctest", # Test snippets in the documentation + "sphinx.ext.viewcode", # Add links to highlighted source code + "sphinx.ext.intersphinx", # Link to other projects’ documentation + "sphinx.ext.autosectionlabel", # Allow reference sections using its title + # allows BibTeX citations to be inserted into documentation generated by Sphinx "sphinxcontrib.bibtex", + "jupyter_sphinx", "numpydoc", + "m2r2", # compatibility with markdown + "myst_nb", # to include jupyter notebooks ] -bibtex_bibfiles = ["refs.bib"] + +# ----------------------------------------------------------------------------- +# Autosummary +# ----------------------------------------------------------------------------- + +# autosummaries from source-files +autosummary_generate = True +# dont show __init__ docstring +autoclass_content = "class" +# sort class members +autodoc_member_order = "groupwise" +# autodoc_member_order = 'bysource' + +autosummary_generate = True # Turn on sphinx.ext.autosummary +# This is because numpydoc screew up with autosummary +# numpydoc_show_class_members=False + +# Napoleon settings +napoleon_google_docstring = True +napoleon_numpy_docstring = True +napoleon_include_init_with_doc = False +napoleon_include_private_with_doc = True +napoleon_include_special_with_doc = False +# napoleon_use_admonition_for_examples = False +# napoleon_use_admonition_for_notes = False +# napoleon_use_admonition_for_references = False +napoleon_use_ivar = True +# napoleon_use_param = True +# napoleon_use_rtype = True + +# Add any paths that contain templates here, relative to this directory. +templates_path = ["_templates"] + +# Add any external modules you want to refer to in the docs here. +intersphinx_mapping = { + "python": ("https://docs.python.org/3", None), + "numpy": ("https://numpy.org/doc/stable/", None), +} + +# The suffix(es) of source filenames. +# You can specify multiple suffix as a list of string: +# source_suffix = ['.rst', '.md'] +source_suffix = [".rst", ".md", ".ipynb"] + +# The encoding of source files. +# source_encoding = 'utf-8-sig' + +# The master toctree document. +master_doc = "index" + +# The language for content autogenerated by Sphinx. Refer to documentation +# for a list of supported languages. +# +# This is also used if you do content translation via gettext catalogs. +# Usually you set "language" from the command line for these cases. language = "en" -html_theme = "sphinx_rtd_theme" +# List of patterns, relative to source directory, that match files and +# directories to ignore when looking for source files. +exclude_patterns = ["build", "_templates/*.rst'", "Thumbs.db", ".DS_Store"] + +# The name of the Pygments (syntax highlighting) style to use. +pygments_style = "sphinx" + +# ----------------------------------------------------------------------------- +# HTML output +# ----------------------------------------------------------------------------- + +# The theme to use for HTML and HTML Help pages. See the documentation for +# a list of builtin themes. +# html_theme = "sphinx_rtd_theme" +html_theme = "pydata_sphinx_theme" + +# TODO: when ies gets a nice logo +# html_logo = '_static/logo.svg' +# html_favicon = '_static/favicon.ico' + +# Theme options are theme-specific and customize the look and feel of a theme +# further. For a list of options available for each theme, see the +# documentation. +# +html_theme_options = { + # "google_analytics_id": "UA-140243896-1", + "show_prev_next": False, + "github_url": "https://github.com/equinor/iterative_ensemble_smoother", + "icon_links": [ + { + "name": "Support", + "url": "https://github.com/equinor/iterative_ensemble_smoother/issues", + "icon": "fa fa-comment fa-fw", + }, + # { + # "name": "The Paper", + # "url": "https://doi.org/10.21105/joss.01450", + # "icon": "fa fa-file-text fa-fw", + # }, + ], +} + +# Title displayed on the html page +html_title = "Iterative Ensemble Smoother" + +# ----------------------------------------------------------------------------- +# Bibliography +# ----------------------------------------------------------------------------- + +bibtex_bibfiles = ["./refs.bib"] +bibtex_default_style = "unsrt" +bibtex_reference_style = "author_year" +suppress_warnings = ["bibtex.duplicate_citation", "autosectionlabel.*"] numpydoc_class_members_toctree = False + +# ----------------------------------------------------------------------------- +# Myst configuration +# ----------------------------------------------------------------------------- + +# see https://myst-parser.readthedocs.io/en/latest/syntax/optional.html#syntax-mathjax +# for math parsing in jupyter notebooks +myst_enable_extensions = [ + "amsmath", + "dollarmath", +] diff --git a/docs/source/examples.rst b/docs/source/examples.rst index 55edb7b0..d6506b3e 100644 --- a/docs/source/examples.rst +++ b/docs/source/examples.rst @@ -7,3 +7,5 @@ Example Usage Polynomial Oscillator + +* :ref:`genindex` diff --git a/docs/source/glossary.rst b/docs/source/glossary.rst index 4e85de9d..8551ee35 100644 --- a/docs/source/glossary.rst +++ b/docs/source/glossary.rst @@ -17,8 +17,8 @@ Glossary forward model A forward model maps model parameters to predicted measurements - (:term:`prediction`). See, for instance, :cite:`evensen2018analysis` and - :cite:`evensen2019efficient`. + (:term:`prediction`). See, for instance, :cite:t:`evensen2018analysis` and + :cite:t:`evensen2019efficient`. error function The error function, or erf_, is a function that for a diff --git a/docs/source/index.rst b/docs/source/index.rst index 37976b64..fc61b3f1 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -1,6 +1,6 @@ -=========================== -Iterative Ensemble Smoother -=========================== +**Date**: |today| **Version**: |version| + +.. mdinclude:: ../../README.md .. toctree:: :hidden: @@ -12,26 +12,6 @@ Iterative Ensemble Smoother bibliography -A library for the iterative ensemble smoothers. - -.. currentmodule:: iterative_ensemble_smoother - -Currently, two main algorithms are implemented: - -* :class:`SIES` - Subspace Iterative Ensemble Smoother - based on the method developed in :cite:`evensen2018analysis`. -* :class:`ESMDA` - Ensemble Smoother with Multiple Data Assimilation - based on the method developed in :cite:`EMERICK2013`. - - -Installation -============ - -Install from `PyPI `_: - -.. code-block:: console - - pip install iterative_ensemble_smoother - - Indices and tables ================== diff --git a/docs/source/refs.bib b/docs/source/refs.bib index aa587b04..74daf105 100644 --- a/docs/source/refs.bib +++ b/docs/source/refs.bib @@ -43,3 +43,18 @@ @article{EMERICK2013 url = {https://www.sciencedirect.com/science/article/pii/S0098300412000994}, keywords = {History matching, Ensemble smoother, Ensemble Kalman filter, Multiple data assimilation}, } + +@article{emerickHistoryMatchingTimelapse2012, + title = {History Matching Time-Lapse Seismic Data Using the Ensemble {{Kalman}} Filter with Multiple Data Assimilations}, + author = {Emerick, Alexandre A. and Reynolds, Albert C.}, + year = {2012}, + month = jun, + journal = {Computational Geosciences}, + volume = {16}, + number = {3}, + pages = {639--659}, + issn = {1573-1499}, + doi = {10.1007/s10596-012-9275-5}, + langid = {english}, + keywords = {Ensemble Kalman filter,Multiple data assimilations,Time-lapse seismic}, +} diff --git a/pyproject.toml b/pyproject.toml index 95e99431..8fea3bfe 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -39,19 +39,22 @@ text = "GPL-3.0" [project.optional-dependencies] doc = [ "sphinx", - "sphinx-rtd-theme", - "nbsphinx", + "pydata_sphinx_theme", + "jupyter_sphinx", "sphinxcontrib.bibtex", "matplotlib", + "pygments", "jupytext", "pandas", - "scipy", "p_tqdm", "ipykernel", "ipywidgets", "numpydoc", + "scipy", + "jinja2", + "m2r2", + "myst_nb" ] - dev = [ "pytest", "pytest-snapshot", diff --git a/src/iterative_ensemble_smoother/__init__.py b/src/iterative_ensemble_smoother/__init__.py index 6c9f1215..6464a67c 100644 --- a/src/iterative_ensemble_smoother/__init__.py +++ b/src/iterative_ensemble_smoother/__init__.py @@ -1,7 +1,31 @@ -""" Implementation of the iterative ensemble smoother history matching algorithms -from Evensen et al. "Efficient Implementation of an Iterative Ensemble Smoother -for Data Assimilation and Reservoir History Matching" -https://www.frontiersin.org/articles/10.3389/fams.2019.00047/full +""" +Purpose +======= + +**iterative_ensemble_smoother** is an open-source, pure python, +and object-oriented library that provides +a user friendly implementation of history matching algorithms +from :cite:t:`evensen2019efficient`. + +The following functionalities are directly provided on module-level. + +Classes +======= + +.. autosummary:: + :toctree: _autosummary + + SIES + ESMDA + +Functions +========= + +.. autosummary:: + :toctree: _autosummary + + steplength_exponential + """ try: from ._version import version as __version__ diff --git a/src/iterative_ensemble_smoother/esmda.py b/src/iterative_ensemble_smoother/esmda.py index 8255d5f3..2707bfea 100644 --- a/src/iterative_ensemble_smoother/esmda.py +++ b/src/iterative_ensemble_smoother/esmda.py @@ -21,6 +21,7 @@ https://helper.ipam.ucla.edu/publications/oilws3/oilws3_14147.pdf """ + import numbers from typing import Optional, Union @@ -36,9 +37,10 @@ class ESMDA: - """Initialize Ensemble Smoother with Multiple Data Assimilation (ES-MDA). + """ + Implement an Ensemble Smoother with Multiple Data Assimilation (ES-MDA). - The implementation follows the 2013 paper by Emerick et al. + The implementation follows :cite:t:`EMERICK2013`. Parameters ---------- @@ -89,6 +91,7 @@ def __init__( seed: Union[np.random._generator.Generator, int, None] = None, inversion: str = "exact", ) -> None: + """Initialize the instance.""" # Validate inputs if not (isinstance(covariance, np.ndarray) and covariance.ndim in (1, 2)): raise TypeError( diff --git a/src/iterative_ensemble_smoother/esmda_inversion.py b/src/iterative_ensemble_smoother/esmda_inversion.py index 5d1586d9..549a9c85 100644 --- a/src/iterative_ensemble_smoother/esmda_inversion.py +++ b/src/iterative_ensemble_smoother/esmda_inversion.py @@ -77,7 +77,7 @@ def empirical_cross_covariance( def normalize_alpha(alpha: npt.NDArray[np.double]) -> npt.NDArray[np.double]: """Assure that sum_i (1/alpha_i) = 1. - This is Eqn (22) in the 2013 Emerick paper. + This is Eqn (22) in :cite:t:`EMERICK2013`. Examples -------- @@ -135,10 +135,10 @@ def singular_values_to_keep( # # C_MD @ inv(C_DD + alpha * C_D) @ (D - Y) # -# where C_MD = empirical_cross_covariance(X, Y) = -# center(X) @ center(Y).T / (X.shape[1] - 1) -# C_DD = empirical_cross_covariance(Y, Y) = -# center(Y) @ center(Y).T / (Y.shape[1] - 1) +# where C_MD = empirical_cross_covariance(X, Y) = center(X) @ center(Y).T +# / (X.shape[1] - 1) +# C_DD = empirical_cross_covariance(Y, Y) = center(Y) @ center(Y).T +# / (Y.shape[1] - 1) # # The methods can be classified as # - exact : with truncation=1.0, these methods compute the exact solution @@ -264,7 +264,9 @@ def inversion_exact_rescaled( ) -> npt.NDArray[np.double]: """Compute a rescaled inversion. - See Appendix A.1 in Emerick et al (2012) for details regarding this approach.""" + See Appendix A.1 in :cite:t:`emerickHistoryMatchingTimelapse2012` + for details regarding this approach. + """ C_DD = empirical_cross_covariance(Y, Y) if C_D.ndim == 2: @@ -338,9 +340,6 @@ def inversion_exact_subspace_woodbury( (A + U @ U.T)^-1 = A^-1 - A^-1 @ U @ (1 + U.T @ A^-1 @ U )^-1 @ U.T @ A^-1 to compute inv(C_DD + alpha * C_D). - - - """ # Woodbury: @@ -397,7 +396,7 @@ def inversion_subspace( X: npt.NDArray[np.double], truncation: float = 1.0, ) -> npt.NDArray[np.double]: - """See Appendix A.2 in Emerick et al (2012) + """See Appendix A.2 in :cite:t:`emerickHistoryMatchingTimelapse2012`. This is an approximate solution. The approximation is that when U, w, V.T = svd(D_delta) @@ -484,10 +483,10 @@ def inversion_rescaled_subspace( X: npt.NDArray[np.double], truncation: float = 1.0, ) -> npt.NDArray[np.double]: - """See Appendix A.2 in Emerick et al (2012) + """ + See Appendix A.2 in :cite:t:`emerickHistoryMatchingTimelapse2012`. Subspace inversion with rescaling. - """ # TODO: I don't understand why this approach is not approximate, when # `inversion_subspace` is approximate diff --git a/src/iterative_ensemble_smoother/sies.py b/src/iterative_ensemble_smoother/sies.py index 21bc8f61..2ba5a1e1 100644 --- a/src/iterative_ensemble_smoother/sies.py +++ b/src/iterative_ensemble_smoother/sies.py @@ -17,50 +17,13 @@ class SIES: - """ - Initialize a Subspace Iterative Ensemble Smoother (SIES) instance. - - This is an implementation of the algorithm described in the paper: - Efficient Implementation of an Iterative Ensemble Smoother for - Data Assimilation and Reservoir History Matching - written by Evensen et al (2019), - URL: https://www.frontiersin.org/articles/10.3389/fams.2019.00047/full - - Parameters - ---------- - parameters : npt.NDArray[np.double] - A 2D array of shape (num_parameters, ensemble_size). Each row corresponds - to a parameter in the model, and each column corresponds to an ensemble - member (realization). This is X in Evensen (2019). - covariance : npt.NDArray[np.double] - Either a 1D array of diagonal covariances, or a 2D covariance matrix. - The shape is either (num_observations,) or (num_observations, num_observations). - This is C_dd in Evensen (2019), and represents observation or measurement - errors. We observe d from the real world, y from the model g(x), and - assume that d = y + e, where the error e is multivariate normal with - covariance given by `covariance`. - observations : npt.NDArray[np.double] - A 1D array of observations, with shape (num_observations,). - This is d in Evensen (2019). - inversion : str - The type of inversion used in the algorithm. Every inversion method - scales the variables. The default is `subspace.` - The options are: - - - `direct`: Solve Eqn (42) directly, which involves inverting a - matrix of shape (num_parameters, num_parameters). - - `subspace_exact` : Solve Eqn (42) using Eqn (50), i.e., the Woodbury - lemma to invert a matrix of size (ensemble_size, ensemble_size). - - `subspace_projected` : Solve Eqn (42) using Section 3.3, i.e., - by projecting the covariance onto S. This approach utilizes the - truncation factor `truncation`. - - truncation : float - How much of the total energy (singular values squared) to keep in the - SVD when `inversion` equals `subspace_projected`. Choosing 1.0 - retains all information, while 0.0 removes all information. - seed : Union[None, int, np.random._generator.Generator], optional - Integer used to seed the random number generator. The default is None. + r""" + Implement a Sequential Iterative Ensemble Smoother (SIES) for Data Assimilation. + + This is an implementation of the algorithm described in the paper: "Efficient + Implementation of an Iterative Ensemble Smoother for Data Assimilation + and Reservoir History Matching", written by :cite:t:`evensen2019efficient`. + """ inversion_funcs = { @@ -78,7 +41,50 @@ def __init__( inversion: str = "subspace_exact", truncation: float = 1.0, seed: Union[None, int, np.random._generator.Generator] = None, - ): + ) -> None: + """ + Initialize the instance. + + Parameters + ---------- + parameters : npt.NDArray[np.double] + A 2D array of shape (num_parameters, ensemble_size). Each row corresponds + to a parameter in the model, and each column corresponds to an ensemble + member (realization). This is X in Evensen (2019). + covariance : npt.NDArray[np.double] + Either a 1D array of diagonal covariances, or a 2D covariance matrix. + The shape is either (num_observations,) or (num_observations, + num_observations). + This is C_dd in Evensen (2019), and represents observation or measurement + errors. We observe d from the real world, y from the model g(x), and + assume that d = y + e, where the error e is multivariate normal with + covariance given by `covariance`. + observations : npt.NDArray[np.double] + A 1D array of observations, with shape (num_observations,). + This is d in Evensen (2019). + inversion : str + The type of inversion used in the algorithm. Every inversion method + scales the variables. The default is `subspace.` + The options are: + + * `direct`: + Solve Eqn (42) directly, which involves inverting a + matrix of shape (num_parameters, num_parameters). + * `subspace_exact` : + Solve Eqn (42) using Eqn (50), i.e., the Woodbury + lemma to invert a matrix of size (ensemble_size, ensemble_size). + * `subspace_projected` : + Solve Eqn (42) using Section 3.3, i.e., by projecting the covariance + onto S. This approach utilizes the truncation factor `truncation`. + + truncation : float + How much of the total energy (singular values squared) to keep in the + SVD when `inversion` equals `subspace_projected`. Choosing 1.0 + retains all information, while 0.0 removes all information. + The default is 1.0. + seed : Union[None, int, np.random._generator.Generator], optional + Integer used to seed the random number generator. The default is None. + """ _validate_inputs( parameters=parameters, covariance=covariance, observations=observations ) @@ -203,7 +209,6 @@ def sies_iteration( This method implements lines 4-9 in Algorithm 1. It returns an updated X and updates the internal state W. - Parameters ---------- responses : npt.NDArray[np.double] @@ -240,7 +245,7 @@ def sies_iteration( def propose_W( self, responses: npt.NDArray[np.double], step_length: float = 0.5 ) -> npt.NDArray[np.double]: - """Returns a proposal for W_i, without updating the internal W. + """Return a proposal for W_i, without updating the internal W. This is an implementation of lines 4-8 in Algorithm 1. @@ -322,7 +327,7 @@ def propose_W_masked( ensemble_mask: npt.NDArray[np.bool_], step_length: float = 0.5, ) -> npt.NDArray[np.double]: - """Returns a proposal for W_i, without updating the internal W. + """Return a proposal for W_i, without updating the internal W. This is an implementation of lines 4-8 in Algorithm 1. diff --git a/src/iterative_ensemble_smoother/utils.py b/src/iterative_ensemble_smoother/utils.py index b2e55040..ab2190a3 100644 --- a/src/iterative_ensemble_smoother/utils.py +++ b/src/iterative_ensemble_smoother/utils.py @@ -14,13 +14,12 @@ def steplength_exponential( max_steplength: float = 0.6, halflife: float = 1.5, ) -> float: - """ - This is an implementation of Eq. (49), which calculates a suitable step length for - the update step, from the book: + r""" + Compute a suitable step length for the update step. - Geir Evensen, Formulating the history matching problem with - consistent error statistics, - Computational Geosciences (2021) 25:945 –970 + This is an implementation of Eq. (49), which calculates a suitable step length for + the update step, from the book: \"Formulating the history matching problem with + consistent error statistics", written by :cite:t:`evensen2021formulating`. Examples --------