Skip to content

update_schemes⚓︎

Ensemble analysis/conditioning/inversion schemes.

AnalysisBindingMixin ⚓︎

Resolve an analysis flavour to an analysis object and delegate to it.

bind_analysis(analysis) ⚓︎

Bind the analysis for analysis, unless a mixin already supplies one.

Nothing shipped in this repository takes that path today (see the module docstring); it remains for a scheme that mixes an analysis directly into its bases instead of listing it in COMPATIBLE_ANALYSES, in which case it keeps the inherited implementation and binds nothing.

resolve_analysis(analysis=None, keys_da=None) ⚓︎

Decide the flavour: explicit argument, else the config, else "approx".

update(*args, **kwargs) ⚓︎

Delegate the analysis step to the bound analysis.

Only reached when nothing else in the MRO defines update; a mixed-in flavour takes precedence and never gets here.

AssimilationResult ⚓︎

Bases: OptimizeResult

Result of an assimilation run.

A dict subclass with attribute access, mirroring :class:scipy.optimize.OptimizeResult so that PIPT and POPT results can be handled the same way. Typical fields:

nit Number of accepted iterations. success Whether the run stopped on a convergence criterion rather than by exhausting maxiter. message Human-readable reason the run stopped. why_stop Mapping of criterion name to whether it fired. data_misfit / prior_data_misfit Final and initial mean data misfit.

AssimilationScheme ⚓︎

Bases: AnalysisBindingMixin, RestartMixin, ABC

What every iterative ensemble data-assimilation scheme inherits.

Subclasses implement :meth:update_step, which performs one iteration and reports what it produced. Everything else is here: the loop, convergence bookkeeping, restart files, the run table, the result object, and the diagnostics and artifact saving that surround a run.

Those last two used to be a separate AssimilationWorkflowMixin that a combined class mixed in ahead of the loop. The split bought nothing -- every shipped scheme wanted both halves -- and cost the reader two classes and one load-bearing MRO order, in which listing the mixin second silently stopped a run from saving anything.

RESTART_ATTRIBUTES: tuple = () ⚓︎

Attributes a scheme needs restored to resume mid-run: what its iterations change and what it drew at construction (perturbed observations, a damping parameter). The loop's own bookkeeping and the ensemble's state are covered by the base state; a subclass only names what it adds. Missing names are skipped, so a scheme that has not yet set one of them checkpoints fine.

__init__(ensemble, **options) ⚓︎

Parameters:

Name Type Description Default
ensemble object

Collaborator satisfying the ensemble protocol described in the module docstring. Owns the state, the observed data and the forward simulator.

 required
**options

Scheme configuration.

  • maxiter: Maximum number of accepted iterations (default: 100).
  • misfit_tol: Relative data-misfit tolerance for convergence (default: 0.01). The assimilation counterpart of an optimizer's ftol.
  • step_tol: Absolute tolerance on the norm of the state update (default: 1e-8). Counterpart of an optimizer's xtol.
  • restart: Restore from a restart file on startup (default: False).
  • restartsave: Write a restart file after the prior forecast and each accepted iteration (default: False).
  • restart_file: Path for the restart file (default: '{scheme_name}_restart.pkl'). Config-driven schemes take these three from the [dataassim] block via :func:restart_options.
 {}

after_accepted_iteration() ⚓︎

Persist iteration artifacts and run QA/QC after an accepted update.

after_analysis() ⚓︎

Between analysis and forecast.

The odd one out: it marks a point inside :meth:update_step, and this class does not dictate the shape of a step, so a scheme calls it itself. The rest of the hooks here are called by :meth:run_assimilation. Nothing runs here at present; it used to refresh QA/QC's variance after data screening, which is no longer supported.

after_forecast(state) ⚓︎

Between forecast and scoring: replace outlier members.

Ordering matters -- outliers are replaced before the misfit is scored, so the replacement feeds into the number the scheme sees. The resampled state is returned rather than written back, so the caller keeps ownership of what it is forecasting.

after_loop(converged) ⚓︎

Save the posterior and the reason the run stopped.

after_prior_forecast() ⚓︎

Handle the prior forecast: prior QA, saved artifacts.

Outlier replacement is not done here. The prior goes through :meth:after_forecast like every other forecast, so it has already happened by the time this runs -- and before :meth:record_prior_score computes the misfit, which is the order that matters.

assimilate(*args, **options) ⚓︎

Construct the scheme and run it to completion.

The assimilation counterpart of scipy.optimize.minimize: one call that builds the scheme, runs every iteration, and returns the outcome. Use it when the scheme object itself is not needed afterwards; when it is, construct the class and call :meth:run_assimilation instead.

Every argument is forwarded verbatim to the constructor, so this accepts whatever the scheme accepts rather than imposing a second signature.

Parameters:

Name Type Description Default
*args

Positional arguments for the constructor. For the shipped PIPT schemes that is (keys_da, keys_en, sim) -- the parsed data-assimilation config, the parsed ensemble config, and the forward simulator -- from which the scheme builds its own ensemble. A scheme written directly against the collaborator protocol is handed its ensemble here instead.

 ()
**options

Keyword arguments for the constructor, such as analysis to override the flavour named in the config.

 {}

Returns:

Type Description
AssimilationResult

Outcome of the run. x is the posterior state ensemble, nit the number of accepted iterations, data_misfit and prior_data_misfit the final and initial mean misfits, and message the reason the run stopped.

Examples:

>>> keys_da, keys_sim, keys_en = read_config.read("case.toml")
>>> result = ESMDA.assimilate(keys_da, keys_en, flow(keys_sim))
>>> result.prior_data_misfit, result.data_misfit
(539.2, 70.1)

Overriding the flavour named in the config:

>>> result = ESMDA.assimilate(keys_da, keys_en, sim, analysis="subspace")
Notes

success reports whether the run stopped on a convergence criterion rather than by exhausting maxiter. Schemes with a fixed iteration schedule -- ES-MDA in particular -- therefore finish normally with success=False, which is expected rather than a failure.

See Also

run_assimilation : Run an already-constructed scheme.

build_ensemble(keys_da, keys_en, sim, ensemble=None) ⚓︎

The collaborator to run on: ensemble if given, else a fresh ENSEMBLE_CLASS.

Handing one in lets two schemes share a prior and its forecasts, and lets a test substitute a stand-in without the config, data files and simulator a real ensemble needs.

check_convergence() ⚓︎

Check scheme-specific convergence criteria.

Returns:

Type Description
bool

True if a subclass-specific stopping criterion is satisfied. The default implementation never stops the loop.

check_misfit_convergence() ⚓︎

Check convergence on the relative change in mean data misfit.

check_state_convergence() ⚓︎

Check convergence on the norm of the state update.

The counterpart of :meth:popt.optimization_methods.optimizer_base. OptimizerBase.check_state_convergence, which compares xk against xk_old. enX_old is snapshotted by :meth:run_assimilation before each attempt, but only when step_tol > 0 -- see there for why.

Opt-in in practice: every shipped scheme passes step_tol=0.0, because ‖Δx‖₂ over a state that mixes variables on different scales (log-permeability alongside saturations, say) has no tolerance that is meaningful across cases. The base default of 1e-8 is small enough to mean "the state did not move at all" rather than being a guess at a scale.

log_columns(prior_run=False) ⚓︎

Trailing columns for the run table -- typically the scheme's control parameter, e.g. {"λ": self.lam}. Empty by default.

log_update(success=None, prior_run=False) ⚓︎

Log one row of the run table.

Called by :meth:run_assimilation -- once for the prior and once per accepted iteration -- so a scheme gets its rows without asking, and the attempts it makes inside :meth:update_step stay its own business. The row is the same for every scheme apart from its control parameter, which :meth:log_columns supplies.

propose_state(result, step_scale=1.0) ⚓︎

The trial state an analysis result implies.

Parameters:

Name Type Description Default
result AnalysisResult or array - like

What self.update(...) returned. A plain array is a state-space step.

 required
step_scale float

Step length applied to the step (GN-EnRML's gamma); 1 for schemes without one.

 1.0

Returns:

Type Description
ndarray

The state to forecast. Weight-space results also advance self.W from self.current_W; the scheme commits W to current_W when it accepts the step.

record_prior_score() ⚓︎

Score the prior forecast and record it, before any iteration.

Sets prior_data_misfit_mean, data_misfit_mean and the per-realisation ensemble_misfit, so the prior is described by the same attributes as every later iteration -- and early enough that the iteration-0 artifacts written by :meth:after_prior_forecast can capture them.

This used to be a score_prior() hook that each scheme implemented, which meant every scheme spelled out both the misfit expression and the five assignments around it. The expression is now :meth:score and the bookkeeping is here; a scheme customises the former.

Does nothing when :meth:score reports no misfit, which is how a scheme with nothing to score opts out.

run_assimilation() ⚓︎

Run this scheme's assimilation to completion.

Named for the job rather than the mechanism, and matching the run_forecast already on this class. The counterpart in popt is OptimizerBase.run_optimization.

Restores a checkpoint if configured, forecasts and scores the prior, then calls :meth:update_step until a convergence criterion fires or maxiter iterations have been taken. One call is one iteration: a scheme that retries -- re-damping, backtracking a step length -- does so inside :meth:update_step, so a report coming back rejected means it has run out of attempts, and the run stops rather than asking again for a step it just said it could not find.

Convergence is checked on rejected reports too, before that stop takes effect: a scheme's :meth:check_convergence can legitimately fire on a step it is about to reject (a stalled misfit that did not actually improve), and that verdict decides how the run is reported.

Returns:

Type Description
AssimilationResult

Populated result object, also stored on self.results.

run_forecast(state) ⚓︎

Forecast state, then run the post-forecast step.

Returns the state to carry forward -- the same one unless :meth:after_forecast replaced members in it.

score(pred_data=None) ⚓︎

Per-realisation data misfit of a forecast.

Called every time a new state has been forecast and needs a number: once for the prior, by :meth:record_prior_score, and then by each scheme for every attempt it takes inside :meth:update_step. One definition per scheme, rather than the same expression repeated in a prior-scoring hook and again in the step.

Parameters:

Name Type Description Default
pred_data optional

The forecast to score -- a PETDataFrame or an (nd, ne) matrix. Defaults to self.pred_data, which is what the ensemble's most recent forecast produced, so the usual call is self.score() straight after run_forecast. Pass one explicitly to score a forecast the ensemble no longer holds.

 None

Returns:

Type Description
ndarray or None

(ne,) misfit per realisation, or None when the scheme has no observation ensemble bound -- a scheme that scores some other way overrides this, and one that reports no misfit at all (the base's own tests) leaves the loop's misfit bookkeeping alone.
Notes

The default is the objective function every shipped scheme uses,

.. math::

\Phi_j = (g(m_j) - d_j)^{\mathsf T} C_d^{-1} (g(m_j) - d_j),

against the perturbed observations enObs and the data covariance cov_data. ES-MDA overrides it to score against an un-inflated copy of the perturbations (enObs_conv); the multilevel scheme to score all fidelity levels at once.

update_step() ⚓︎

Perform one scheme-specific analysis step.

Implementations compute the analysis update, apply it to the ensemble state, run the resulting forecast, and score the result. How they do that is entirely theirs -- the base calls this and nothing inside it.

Returns:

Type Description
StepReport

accepted decides whether the loop advances or gives the scheme another attempt at the same iteration number, which is how the Levenberg-Marquardt schemes back off by increasing their damping parameter. misfit is the per-realisation data misfit as of now; the loop derives data_misfit and data_misfit_std from it.

ES ⚓︎

Bases: EnKF

Ensemble Smoother (ES).

Assimilates all observations simultaneously in a single update, rather than sequentially in time as the filter does. It is :class:EnKF specialised to one data group, and shares its analysis step; only the iteration budget and the misfit bookkeeping differ.

A single conditioning step is cheap but can over-correct when the model is strongly non-linear. :class:ESMDA addresses this by spreading the same update over several inflated steps.

Parameters:

Name Type Description Default
keys_da dict

Parsed dataassim configuration. Besides the keys every scheme reads -- data, datavar, obsname, truedataindex -- the ones this scheme acts on are listed under Notes.

 required
keys_en dict

Parsed ensemble configuration: ensemble size ne, the state variable names, and the prior_<name> blocks describing each.

 required
sim object

Forward simulator instance, e.g. simulator.opm.flow.

 required
analysis (approx, full, subspace)

Analysis flavour, i.e. how the ensemble-approximated sensitivity is inverted. Defaults to the analysis key in keys_da, falling back to 'approx'. The flavours differ in cost and in how they handle a rank-deficient ensemble; they solve the same update equation.

 'approx'

Attributes:

Name Type Description
ensemble AssimilationEnsemble

Collaborator holding the state realisations, observed data and simulator. Its state is exposed as properties on the scheme, so scheme.enX and scheme.keys_da read straight through.
analysis AnalysisBase

The bound analysis object. Note the constructor takes analysis as a name and this attribute holds the resulting object, the way Model(optimizer="adam").optimizer is an optimizer instance.
analysis_name str

The flavour name that was resolved, e.g. 'approx'.
iteration int

Accepted iterations completed so far.
data_misfit, prior_data_misfit float

Current and initial mean data misfit.
Notes

assimindex is flattened to a single group at construction, so the ordering that matters for :class:EnKF has no effect here.

Because there is only one step, the full flavour coincides with approx -- the prior-increment term they differ over is only reached when iterating -- so :attr:EnKF.COMPATIBLE_ANALYSES, inherited unchanged here, points "full" at the cheaper approx analysis.

Examples:

>>> result = ES.assimilate(keys_da, keys_en, flow(keys_sim))
>>> result.nit
1
References

Evensen, Data Assimilation: The Ensemble Kalman Filter evensen2009a.

See Also

EnKF : Sequential form of the same update. ESMDA : Spreads the conditioning over several inflated steps.

__init__(keys_da, keys_en, sim, analysis=None, ensemble=None) ⚓︎

Build the ensemble from the config (or take the one given) and bind the analysis.

See the class docstring for the parameters.

check_convergence() ⚓︎

ES takes a single all-data-at-once step; nothing stops early.

score_and_commit() ⚓︎

Calculate the "convergence" of the method. Important to

ESMDA ⚓︎

Bases: AssimilationScheme

Ensemble Smoother with Multiple Data Assimilation (ES-MDA).

An iterative ensemble smoother that assimilates all data repeatedly over a fixed number of steps, inflating the data-error covariance at each one so that the repeated conditioning does not over-fit. With inflation factors :math:\alpha_i satisfying :math:\sum_i 1/\alpha_i = 1, each step applies

.. math::

m \leftarrow m + C_{md} (C_{dd} + \alpha_i C_d)^{-1} (d_{obs} - g(m))

with the observations re-perturbed as :math:d_{obs} = d_{true} + \sqrt{\alpha_i} C_d^{1/2} Z.

The schedule is fixed rather than convergence-driven, so a run normally ends by exhausting its steps and reports success=False. That is the expected outcome, not a failure.

Parameters:

Name Type Description Default
keys_da dict

Parsed dataassim configuration. Besides the keys every scheme reads -- data, datavar, obsname, truedataindex -- the ones this scheme acts on are listed under Notes.

 required
keys_en dict

Parsed ensemble configuration: ensemble size ne, the state variable names, and the prior_<name> blocks describing each.

 required
sim object

Forward simulator instance, e.g. simulator.opm.flow.

 required
analysis (approx, full, subspace)

Analysis flavour, i.e. how the ensemble-approximated sensitivity is inverted. Defaults to the analysis key in keys_da, falling back to 'approx'. The flavours differ in cost and in how they handle a rank-deficient ensemble; they solve the same update equation.

 'approx'

Attributes:

Name Type Description
ensemble AssimilationEnsemble

Collaborator holding the state realisations, observed data and simulator. Its state is exposed as properties on the scheme, so scheme.enX and scheme.keys_da read straight through.
analysis AnalysisBase

The bound analysis object. Note the constructor takes analysis as a name and this attribute holds the resulting object, the way Model(optimizer="adam").optimizer is an optimizer instance.
analysis_name str

The flavour name that was resolved, e.g. 'approx'.
iteration int

Accepted iterations completed so far.
data_misfit, prior_data_misfit float

Current and initial mean data misfit.
Notes

Configured through the mda block of keys_da:

tot_assim_steps Number of assimilation steps, e.g. 3. inflation_param Inflation factors, one per step, e.g. [3, 3, 3]. Their reciprocals must sum to 1, which is asserted at construction. Defaults to tot_assim_steps repeated, which satisfies the constraint.

Examples:

>>> result = ESMDA.assimilate(keys_da, keys_en, flow(keys_sim))
>>> result.nit
3
References

Emerick and Reynolds, Ensemble smoother with multiple data assimilation emerick2013a.

See Also

ES : Single-step smoother; ES-MDA with one assimilation step. LMEnRML : Iterates to convergence instead of on a fixed schedule.

__init__(keys_da, keys_en, sim, analysis=None, ensemble=None) ⚓︎

Build the ensemble from the config (or take the one given) and bind the analysis.

See the class docstring for the parameters; ensemble is a ready-made collaborator to run on instead of building one.

calc_analysis() ⚓︎

Analysis step of ES-MDA. The analysis algorithm is similar to EnKF analysis, only difference is that the data covariance matrix is inflated with an inflation parameter alpha. The update is done as an iterative smoother where all data is assimilated at once.

Notes

ES-MDA is an iterative ensemble smoother with a predefined number of iterations, where the updates is done with the EnKF update equations but where the data covariance matrix have been inflated:

\[ \begin{align} d_{obs} &= d_{true} + \sqrt{\alpha}C_d^{1/2}Z \\ m &= m_{prior} + C_{md}(C_g + \alpha C_d)^{-1}(g(m) - d_{obs}) \end{align} \]

where \(d_{true}\) is the true observed data, \(\alpha\) is the inflation factor, \(C_d\) is the data covariance matrix, \(Z\) is a standard normal random variable, \(C_{md}\) and \(C_{g}\) are sample covariance matrices, \(m\) is the model parameter, and \(g(\)\) is the predicted data. Note that \(\alpha\) can have a different value in each assimilation step and must fulfill:

\[ \sum_{i=1}^{N_a} \frac{1}{\alpha} = 1 \]

where \(N_a\) being the total number of assimilation steps.

check_convergence() ⚓︎

ES-MDA runs its full schedule of inflated steps; nothing stops early.

log_columns(prior_run=False) ⚓︎

ES-MDA reports the inflation factor for the step just taken.

score(pred_data=None) ⚓︎

Data misfit against the un-inflated perturbed observations.

enObs is redrawn each step with the covariance inflated by alpha[iteration], so scoring against it would compare every iteration to a different yardstick. enObs_conv is the copy taken before any inflation, which is what makes the misfit trajectory comparable across the schedule.

score_and_commit() ⚓︎

Score the forecast that followed the analysis, then commit the step.

Was the second half of check_convergence: ES-MDA never actually tested for convergence there, it recomputed the misfit, logged the iteration and promoted enX_temp. Under the new contract the convergence question lives in :meth:check_convergence and this keeps the bookkeeping.

Returns:

Type Description
dict

The why_stop record, also stored on self.why_stop.

update_step() ⚓︎

Run one ES-MDA assimilation step.

Computes the inflated analysis, forecasts the trial state, then scores the resulting misfit and promotes the state. Scoring after the forecast is what lets outlier replacement, which runs in between, feed into the number the scheme sees.

Returns:

Type Description
bool

Always True. ES-MDA takes a fixed number of inflated steps and never rejects one. The success flag it logs compares the misfit against the previous iteration and is a reporting signal only -- returning it here would make the base class discard accepted steps.

EnKF ⚓︎

Bases: AssimilationScheme

Ensemble Kalman Filter (EnKF).

Assimilates data sequentially, updating the state once per group of observations in the order given by assimindex. Each update applies the Kalman equations with the covariances approximated from the ensemble:

.. math::

m \leftarrow m + C_{md} (C_{dd} + C_d)^{-1} (d_{obs} - g(m))

There is no damping and no rejection: every step is accepted, and the run ends once the data groups are exhausted.

Parameters:

Name Type Description Default
keys_da dict

Parsed dataassim configuration. Besides the keys every scheme reads -- data, datavar, obsname, truedataindex -- the ones this scheme acts on are listed under Notes.

 required
keys_en dict

Parsed ensemble configuration: ensemble size ne, the state variable names, and the prior_<name> blocks describing each.

 required
sim object

Forward simulator instance, e.g. simulator.opm.flow.

 required
analysis (approx, full, subspace)

Analysis flavour, i.e. how the ensemble-approximated sensitivity is inverted. Defaults to the analysis key in keys_da, falling back to 'approx'. The flavours differ in cost and in how they handle a rank-deficient ensemble; they solve the same update equation.

 'approx'

Attributes:

Name Type Description
ensemble AssimilationEnsemble

Collaborator holding the state realisations, observed data and simulator. Its state is exposed as properties on the scheme, so scheme.enX and scheme.keys_da read straight through.
analysis AnalysisBase

The bound analysis object. Note the constructor takes analysis as a name and this attribute holds the resulting object, the way Model(optimizer="adam").optimizer is an optimizer instance.
analysis_name str

The flavour name that was resolved, e.g. 'approx'.
iteration int

Accepted iterations completed so far.
data_misfit, prior_data_misfit float

Current and initial mean data misfit.
Notes

assimindex determines the grouping and ordering of the sequential updates. If all data are to be assimilated in a single step, use :class:ES, which is this scheme specialised to one group.

energy sets the fraction of singular values retained in the truncated SVD (default 0.98); values above 1 are read as percentages.

Every data group is assimilated exactly once, so the prior-increment term that distinguishes full from approx is never reached: "full" is pointed at the same class as "approx" in :attr:COMPATIBLE_ANALYSES. :class:ES inherits this.

Examples:

>>> result = EnKF.assimilate(keys_da, keys_en, flow(keys_sim))
References

Evensen, Data Assimilation: The Ensemble Kalman Filter evensen2009a.

See Also

ES : All-data-at-once form of the same update.

__init__(keys_da, keys_en, sim, analysis=None, ensemble=None) ⚓︎

Build the ensemble from the config and bind the analysis.

See the class docstring for the parameters.

calc_analysis() ⚓︎

Calculate the analysis step of the EnKF procedure. The updating is done using the Kalman filter equations, using svd for numerical stability. Localization is available.

check_convergence() ⚓︎

The EnKF runs its full sweep of data groups; nothing stops early.

score_and_commit() ⚓︎

Calculate the "convergence" of the method. Important to

update_step() ⚓︎

Run one EnKF step: analysis, forecast, then score and commit.

Returns:

Type Description
bool

Always True. The EnKF applies one update per data group and has no rejection path.

GNEnRML ⚓︎

Bases: IterativeEnRML

Gauss-Newton Ensemble Randomized Maximum Likelihood (GN-EnRML).

Solves the same randomized maximum likelihood problem as :class:LMEnRML, but takes undamped Gauss-Newton steps scaled by a step length :math:\gamma \in (0, 1] rather than inflating the Hessian:

.. math::

m \leftarrow m + \gamma \, C_{md} (C_d + C_{dd})^{-1}
(d_{obs} - g(m))

Steps are accepted or rejected on the mean data misfit as in LM-EnRML. On acceptance :math:\gamma is relaxed towards gamma_max; on rejection it is divided by gamma_factor and the step re-solved, in the same within-:meth:update_step loop LM-EnRML uses for :math:\lambda.

Parameters:

Name Type Description Default
keys_da dict

Parsed dataassim configuration. Besides the keys every scheme reads -- data, datavar, obsname, truedataindex -- the ones this scheme acts on are listed under Notes.

 required
keys_en dict

Parsed ensemble configuration: ensemble size ne, the state variable names, and the prior_<name> blocks describing each.

 required
sim object

Forward simulator instance, e.g. simulator.opm.flow.

 required
analysis (approx, full, subspace)

Analysis flavour, i.e. how the ensemble-approximated sensitivity is inverted. Defaults to the analysis key in keys_da, falling back to 'approx'. The flavours differ in cost and in how they handle a rank-deficient ensemble; they solve the same update equation.

 'approx'

Attributes:

Name Type Description
ensemble AssimilationEnsemble

Collaborator holding the state realisations, observed data and simulator. Its state is exposed as properties on the scheme, so scheme.enX and scheme.keys_da read straight through.
analysis AnalysisBase

The bound analysis object. Note the constructor takes analysis as a name and this attribute holds the resulting object, the way Model(optimizer="adam").optimizer is an optimizer instance.
analysis_name str

The flavour name that was resolved, e.g. 'approx'.
iteration int

Accepted iterations completed so far.
data_misfit, prior_data_misfit float

Current and initial mean data misfit.
Notes

Configured through the iteration block of keys_da:

max_iter Maximum accepted iterations. gamma Initial step length (default 0.2). gamma_max Value the step length relaxes towards on success (default 0.5). gamma_factor Divisor applied to the step length on rejection (default 2.5). max_inner_iter Step-length attempts one iteration may make before the run gives up (default 10). There is no gamma_min, so this is what bounds it. data_misfit_tol Relative misfit change treated as converged (default 0.01).

The margis flavour is backed by margIS_update, ported from an older layout. It returns a matrix-form ensemble transform step (AnalysisResult(W_step=...), starting from W = I) rather than the weight step most other flavours use; propose_state reconstructs the state for either. Run against real data it produces a large, sensible misfit reduction, but is still one run on one case with no committed reference pinning it -- see its module docstring (:mod:pipt.update_schemes.analysis.margis) for what was fixed in the port and what remains a modelling choice rather than a bug.

Examples:

>>> result = GNEnRML.assimilate(keys_da, keys_en, flow(keys_sim))
References

Chen and Oliver chen2013; see also Raanes, Stordal and Evensen, Revising the stochastic iterative ensemble smoother raanes2019, and Evensen et al. evensen2019.

See Also

IterativeEnRML : The loop, scoring and bookkeeping both schemes share. LMEnRML : Levenberg-Marquardt form, damped via the Hessian.

log_columns(prior_run=False) ⚓︎

GN-EnRML reports the step length the logged iteration took.

IterativeEnRML ⚓︎

Bases: AssimilationScheme

What LM-EnRML and GN-EnRML share: everything but the control parameter.

Both solve the randomized maximum likelihood problem by repeated linearisation, accept or reject each step on the mean data misfit, retry a rejected step from the same state inside :meth:update_step, and stop on the relative misfit change, on max_inner_iter failed attempts in one iteration, or on max_iter. They differ only in the control parameter that reacts to an attempt: LM-EnRML's damping :math:\lambda inflates the Hessian and grows on rejection; GN-EnRML's step length :math:\gamma scales the step and shrinks on rejection. A subclass supplies that behaviour through the hooks below and nothing else.

Hooks

_read_damping_options(options) Read the control parameter(s) from the iteration block. _step_scale() Factor applied to the analysis step: 1 for LM-EnRML, :math:\gamma for GN-EnRML. _record_control() Remember the control the attempt ran with, for the run table. _control_exhausted() and _exhausted_message() Whether the control itself says stop (LM-EnRML: :math:\lambda \ge lambda_max), and the stop reason to report then. _why_stop_control() The control's entries in why_stop. _on_improved() Accepted with a smaller misfit spread: relax the control. _on_rejected() Rejected: tighten the control. _give_up_message(attempt) Stop reason when max_inner_iter attempts all failed. log_columns() The control's column in the run table.

__init__(keys_da, keys_en, sim, analysis=None, ensemble=None) ⚓︎

Build the ensemble from the config (or take the one given) and bind the analysis.

See the subclass docstrings for the parameters; ensemble is a ready-made collaborator to run on instead of building one.

calc_analysis() ⚓︎

Compute the trial state: the analysis step, scaled and clipped.

check_convergence() ⚓︎

Report the verdict reached by the preceding :meth:score_and_commit.

score_and_commit() ⚓︎

Score the forecast, decide on the attempt, and adjust the control.

Returns:

Name Type Description
why_stop dict

The convergence criteria with their values, including the control's own entries.

update_step() ⚓︎

Run one iteration, retrying until an attempt improves the misfit.

The retry loop is here rather than in the base loop: one call is one iteration, and the attempts it took to get there are this scheme's business. That mirrors popt, where EnOpt.update_step backtracks over its own step length and returns only once it has an improving step or has run out of attempts.

Each attempt re-solves the analysis with the current control parameter, forecasts the proposal and scores it. A worse misfit tightens the control (:meth:_on_rejected) and tries again from the same state -- nothing was committed -- so the retries cost forecasts, not correctness.

Returns:

Type Description
StepReport

accepted is whether an attempt improved the misfit. It is False only when the scheme has also decided to stop, which :meth:check_convergence then reports to the loop.

LMEnRML ⚓︎

Bases: IterativeEnRML

Levenberg-Marquardt Ensemble Randomized Maximum Likelihood (LM-EnRML).

An iterative ensemble smoother that solves the randomized maximum likelihood problem by repeated linearisation, with a Levenberg-Marquardt damping parameter :math:\lambda controlling the step size. The damped update inflates the Hessian approximation:

.. math::

m \leftarrow m + C_{md} \big((1 + \lambda) C_d + C_{dd}\big)^{-1}
(d_{obs} - g(m))

Unlike ES-MDA, steps are accepted or rejected. A step that increases the mean data misfit is discarded, :math:\lambda is multiplied by lambda_factor and the step re-solved from the same state; one that decreases it is kept and :math:\lambda reduced. That retry loop lives inside :meth:update_step, so one iteration is one call however many attempts it takes -- the shape popt's optimizers have. The run stops when the relative misfit change falls below data_misfit_tol, when :math:\lambda reaches lambda_max, when a single iteration exhausts max_inner_iter attempts, or on max_iter.

Parameters:

Name Type Description Default
keys_da dict

Parsed dataassim configuration. Besides the keys every scheme reads -- data, datavar, obsname, truedataindex -- the ones this scheme acts on are listed under Notes.

 required
keys_en dict

Parsed ensemble configuration: ensemble size ne, the state variable names, and the prior_<name> blocks describing each.

 required
sim object

Forward simulator instance, e.g. simulator.opm.flow.

 required
analysis (approx, full, subspace)

Analysis flavour, i.e. how the ensemble-approximated sensitivity is inverted. Defaults to the analysis key in keys_da, falling back to 'approx'. The flavours differ in cost and in how they handle a rank-deficient ensemble; they solve the same update equation.

 'approx'

Attributes:

Name Type Description
ensemble AssimilationEnsemble

Collaborator holding the state realisations, observed data and simulator. Its state is exposed as properties on the scheme, so scheme.enX and scheme.keys_da read straight through.
analysis AnalysisBase

The bound analysis object. Note the constructor takes analysis as a name and this attribute holds the resulting object, the way Model(optimizer="adam").optimizer is an optimizer instance.
analysis_name str

The flavour name that was resolved, e.g. 'approx'.
iteration int

Accepted iterations completed so far.
data_misfit, prior_data_misfit float

Current and initial mean data misfit.
Notes

Configured through the iteration block of keys_da:

max_iter Maximum accepted iterations. lambda Initial damping parameter (default 100). 'auto' derives it from the prior data misfit. lambda_factor Factor by which damping grows on rejection and shrinks on acceptance (default 5). Held as lam_factor -- not gamma, which is GN-EnRML's step length, a different quantity entirely. lambda_max, lambda_min Bounds on the damping parameter. max_inner_iter Damping attempts one iteration may make before the run gives up (default 10). lambda_max normally stops it first. data_misfit_tol Relative misfit change treated as converged (default 0.01).

Examples:

>>> result = LMEnRML.assimilate(keys_da, keys_en, flow(keys_sim))
>>> result.message
'Maximum number of iterations reached'

success distinguishes the two ways a run can end: True when a convergence criterion fired, False when max_iter was reached first. Both are ordinary outcomes -- check prior_data_misfit against data_misfit to judge whether the run achieved anything.

References

Chen and Oliver, Levenberg-Marquardt forms of the iterative ensemble smoother for efficient history matching and uncertainty quantification chen2013.

See Also

IterativeEnRML : The loop, scoring and bookkeeping both schemes share. GNEnRML : Gauss-Newton form, damped by a step length instead. ESMDA : Fixed schedule rather than convergence-driven iteration.

log_columns(prior_run=False) ⚓︎

LM-EnRML reports the damping the logged iteration ran with.

score(pred_data=None) ⚓︎

Data misfit, sizing lambda='auto' the first time there is one.

:math:\lambda_0 = \Phi_{prior} / 2 N_d is defined against the prior misfit, so it cannot be settled in __init__. The first score of a run is the prior's, which makes this the earliest point it can be resolved -- and everything downstream needs a number: the prior row reports λ, and the prior QA/QC pass computes with it.

MultilevelEnsemble ⚓︎

Bases: AssimilationEnsemble

Ensemble whose state is partitioned into fidelity levels.

enX is a list of matrices, one per level, rather than a single (nx, ne) matrix, and the simulator is configured to run each level. Everything else is the ordinary assimilation ensemble.

Attributes:

Name Type Description
enX list of ndarray

State ensemble per level; enX[l] has shape (nx, ml_ne[l]).
tot_level int

Number of fidelity levels.
ml_ne list of int

Ensemble size at each level.

reorganize_ml_prior(enX) ⚓︎

Reorganize prior ensemble to multilevel structure (list of matrices).

StepReport ⚓︎

What one attempt produced. Returned by :meth:update_step.

The base does not dictate how a scheme takes its step; this is what it needs back afterwards, to score convergence, log, and build the result. Required fields are positional, so forgetting one is a TypeError at construction rather than a None surfacing several iterations later.

accepted: bool ⚓︎

Keep this step? False says the scheme found no improving step and has exhausted the attempts it makes inside :meth:update_step, so the loop stops rather than asking for the same step again.

misfit: np.ndarray ⚓︎

Per-realisation data misfit as of now. The loop derives data_misfit and data_misfit_std from it, so the three can no longer drift apart the way separately-assigned attributes could.

"As of now" matters for a scheme that gives up: LM-EnRML restores the last accepted misfit when it backs off, and returns that, so the value the loop records and logs is the one the run actually reached.

state: Any ⚓︎

The state this attempt produced, committed by the loop when accepted. A scheme still writes it to ensemble.enX_temp first, because that is what the forecast predicts on -- but handing it back here is what lets the loop own the commit, rather than every scheme remembering the same two lines. Forgetting them used to give a run that iterated and logged normally while returning the prior untouched.

why_stop: dict | None = None ⚓︎

Criterion record, merged into result.why_stop.

approx_update ⚓︎

Bases: AnalysisBase

Approximate LM Update scheme as defined in "Chen, Y., & Oliver, D. S. (2013). Levenberg–Marquardt forms of the iterative ensemble smoother for efficient history matching and uncertainty quantification. Computational Geosciences, 17(4), 689–703. https://doi.org/10.1007/s10596-013-9351-5". Note that for a EnKF or ES update, or for update within GN scheme, lambda = 0.

update(enX, enY, enE, **kwargs) ⚓︎

Perform the approximate LM update.

Parameters: 
enX : np.ndarray
    State ensemble matrix (nx, ne)

enY : np.ndarray
    Predicted data ensemble matrix (nd, ne)

enE : np.ndarray
    Ensemble of perturbed observations (nd, ne)

esmda_hybrid ⚓︎

Bases: ESMDA

A multilevel implementation of the ES-MDA algorithm with the hybrid gain.

Composes a :class:MultilevelEnsemble and binds hybrid_update for the per-level gain, the same way :class:~pipt.update_schemes.esmda.ESMDA binds approx_update and friends. It is not just ESMDA with an extra flavour, though: its own COMPATIBLE_ANALYSES offers only "hybrid", deliberately narrower than ESMDA's -- approx_update et al. expect a single enX/proj matrix, and this scheme's state is partitioned into one such matrix per level, which those analyses were never written to handle.

Notes

Requires a multilevel block in keys_en giving levels, en_size per level and ml_weights.

calc_analysis() ⚓︎

The ES-MDA analysis over every fidelity level: per-level predictions, redrawn observations, the hybrid update, clipped proposals.

score(pred_data=None) ⚓︎

Data misfit over every fidelity level at once.

pred_data is one frame per level here, so the levels are concatenated along the ensemble axis and scored as a single ensemble against the un-inflated perturbations, as :meth:pipt.update_schemes.esmda.ESMDA.score does for one level.

score_and_commit() ⚓︎

Score the forecast that followed the analysis, then commit the step.

Was the second half of check_convergence. ES-MDA never tested for convergence there; it recomputed the misfit, logged the iteration and promoted enX_temp.

Returns:

Type Description
dict

The why_stop record, also stored on self.why_stop.

subspace_update ⚓︎

Bases: AnalysisBase

Ensemble subspace update (weight-space IES).

The update is formulated in the ensemble weight space W (shape ne × ne) rather than model space, making it efficient when ne ≪ nx. The caller checks self.w_step (not self.step) to apply the update.

References

Raanes, P. N., Stordal, A. S., & Evensen, G. (2019). Revising the stochastic iterative ensemble smoother. Nonlinear Processes in Geophysics, 26(3), 325-338. https://doi.org/10.5194/npg-26-325-2019

Evensen, G., Raanes, P. N., Stordal, A. S., & Hove, J. (2019). Efficient implementation of an iterative ensemble smoother for data assimilation and reservoir history matching. Frontiers in Applied Mathematics and Statistics, 5, 47. https://doi.org/10.3389/fams.2019.00047

update(enX, enY, enE, **kwargs) ⚓︎

Perform the subspace (weight-space) LM update.

Sets self.scheme.w_step (shape ne × ne) and returns None -- the caller applies the weight update, not a state-space step.

Parameters:

Name Type Description Default
enX (ndarray, shape(nx, ne))

State ensemble matrix (unused directly; included for interface parity).

 required
enY (ndarray, shape(nd, ne))

Predicted data ensemble matrix.

 required
enE (ndarray, shape(nd, ne))

Perturbed observations ensemble.

 required

Returns:

Type Description
AnalysisResult

The weight-space step w_step (ne, ne).

gen_real(mean, var, number, rng=None, limits=None, return_chol=False) ⚓︎

Realisations of a Gaussian with the given mean and (co)variance.

Draw for draw the same as geostat.decomp.Cholesky.gen_real -- the same shapes drawn in the same order with the same arithmetic -- so runs are bit-identical to what geostat produced; only the stream is a parameter now.

Parameters:

Name Type Description Default
mean (array - like, shape(n))

Mean vector.

 required
var array - like

Variance vector (n,), covariance matrix (n, n), or a scalar when mean has one element.

 required
number int

Number of realisations.

 required
rng RandomState - like

The stream to draw from; the global one by default.

 None
limits dict

{'lower': ..., 'upper': ...} to clip the realisations to.

 None
return_chol bool

Also return the factor used: sqrt(var) for a diagonal, the upper Cholesky factor otherwise.

 False

Returns:

Type Description
ndarray, shape (n, number), and the factor when ``return_chol``.

restart_options(keys_da) ⚓︎

The checkpoint settings of a config's [dataassim] block, as scheme options.

restart (resume from the checkpoint), restartsave (write one after the prior forecast and every accepted iteration) and restart_file (default <scheme>_restart.pkl). Legacy yes/no strings are accepted. Schemes pass **restart_options(keys_da) to the base so the keys reach :class:~ensemble.checkpoint.RestartMixin; they used to stop at the ensemble, which loaded a pickle of itself and left the scheme's own state -- iteration, damping, misfit history -- at its initial values.