iterative_ensemble_smoother.SIES#

class iterative_ensemble_smoother.SIES(parameters: npt.NDArray[np.double], covariance: npt.NDArray[np.double], observations: npt.NDArray[np.double], *, inversion: str = 'subspace_exact', truncation: float = 1.0, seed: None | int | np.random._generator.Generator = None)[source]#

Bases: object

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 Evensen et al. [2019].

Methods

`evaluate_objective`(W, Y)	Evaluate the objective function in Equation (18), taking the mean over each ensemble member.
`evaluate_objective_elementwise`(W, Y)	Equation (18) elementwise and termwise, returning a tuple of vectors (prior, likelihood).
`propose_W`(responses[, step_length])	Return a proposal for W_i, without updating the internal W.
`propose_W_masked`(responses, ensemble_mask[, ...])	Return a proposal for W_i, without updating the internal W.
`sies_iteration`(responses[, step_length, ...])	Perform a single SIES iteration (Gauss-Newton step).

__init__(parameters: npt.NDArray[np.double], covariance: npt.NDArray[np.double], observations: npt.NDArray[np.double], *, inversion: str = 'subspace_exact', truncation: float = 1.0, seed: None | int | np.random._generator.Generator = None) → None[source]#

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.

Methods definition

evaluate_objective(W: npt.NDArray[np.double], Y: npt.NDArray[np.double]) → npt.NDArray[np.double][source]#

Evaluate the objective function in Equation (18), taking the mean over each ensemble member.

Parameters:

W (npt.NDArray[np.double]) – Current weight matrix W, which represents the current X_i as a linear combination of the prior. See equation (17) and line 9 in Algorithm 1.
Y (npt.NDArray[np.double]) – Responses when evaluating the model at X_i. In other words, Y = g(X_i).

Returns:

The objective function value. Lower objective is better.

Return type:

float

evaluate_objective_elementwise(W: npt.NDArray[np.double], Y: npt.NDArray[np.double]) → tuple[npt.NDArray[np.double], npt.NDArray[np.double]][source]#

Equation (18) elementwise and termwise, returning a tuple of vectors (prior, likelihood).

Parameters:

W (npt.NDArray[np.double]) – Current weight matrix W, which represents the current X_i as a linear combination of the prior. See equation (17) and line 9 in Algorithm 1.
Y (npt.NDArray[np.double]) – Responses when evaluating the model at X_i. In other words, Y = g(X_i).

Returns:

prior (npt.NDArray[np.double]) – A 1D array representing distance from prior for each ensemble member.
likelihood (npt.NDArray[np.double]) – A 1D array representing distance from observations for each ensemble member.

propose_W(responses: npt.NDArray[np.double], step_length: float = 0.5) → npt.NDArray[np.double][source]#

Return a proposal for W_i, without updating the internal W.

This is an implementation of lines 4-8 in Algorithm 1.

Parameters:

responses (npt.NDArray[np.double]) – The model evaluated at X_i. In other words, responses = g(X_i). This is Y in the paper.
step_length (float, optional) – Step length for Gauss-Newton. The default is 0.5.

Returns:

W_i – A proposal for a new W in the algorithm.

Return type:

npt.NDArray[np.double]

propose_W_masked(responses: npt.NDArray[np.double], ensemble_mask: npt.NDArray[np.bool_], step_length: float = 0.5) → npt.NDArray[np.double][source]#

Return a proposal for W_i, without updating the internal W.

This is an implementation of lines 4-8 in Algorithm 1.

Parameters:

responses (npt.NDArray[np.double]) – The model evaluated at X_i. In other words, responses = g(X_i). This is Y in the paper.
step_length (float, optional) – Step length for Gauss-Newton. The default is 0.5.

Returns:

W_i – A proposal for a new W in the algorithm.

Return type:

npt.NDArray[np.double]

sies_iteration(responses: npt.NDArray[np.double], step_length: float = 0.5, ensemble_mask: npt.NDArray[np.bool_] | None = None) → npt.NDArray[np.double][source]#

Perform a single SIES iteration (Gauss-Newton step).

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]) – The model evaluated at X_i. In other words, responses = g(X_i). This is Y in the paper.
step_length (float, optional) – Step length for Gauss-Newton. The default is 0.5.

Returns:

Updated parameter ensemble.

Return type:

npt.NDArray[np.double]