PerturbedStepSolver

class probnum.diffeq.perturbed.step.PerturbedStepSolver(rng, solver, noise_scale, perturb_function)[source]

Bases: probnum.diffeq.ODESolver

Probabilistic ODE solver based on random perturbation of the step-sizes.

Perturbs the steps accordingly and projects the solution back to the originally proposed time points. Proposed by Abdulle and Garegnani (2020) 1.

Parameters
  • rng (Generator) – Random number generator.

  • solver (WrappedScipyRungeKutta) – Currently this has to be a Runge-Kutta method based on SciPy.

  • noise-scale – Scales the amount of noise that is introduced.

  • perturb_function (Callable) – Defines how the stepsize is distributed. This can be either one of perturb_lognormal() or perturb_uniform() or any other perturbation function with the same signature.

References

1

Abdulle, A. and Garegnani, G. Random time step probabilistic methods for uncertainty quantification in chaotic and geometric numerical integration. Statistics and Computing. 2020.

Methods Summary

attempt_step(state, dt)

Perturb the original stopping point.

construct_with_lognormal_perturbation(rng, …)

construct_with_uniform_perturbation(rng, …)

initialize(ivp)

Initialise and reset the solver.

method_callback(state)

Call dense output after each step and store the interpolants.

perform_full_step(state, initial_dt)

Perform a full ODE solver step.

postprocess(odesol)

Process the ODESolution object before returning.

rvlist_to_odesol(times, rvs)

Create an ODESolution object.

solution_generator(ivp[, stop_at, callbacks])

Generate ODE solver steps.

solve(ivp[, stop_at, callbacks])

Solve an IVP.

Methods Documentation

attempt_step(state, dt)[source]

Perturb the original stopping point.

Perform one perturbed step and project the solution back to the original stopping point.

Parameters
Returns

New state.

Return type

_odesolver_state.ODESolverState

classmethod construct_with_lognormal_perturbation(rng, solver, noise_scale)[source]
classmethod construct_with_uniform_perturbation(rng, solver, noise_scale)[source]
initialize(ivp)[source]

Initialise and reset the solver.

method_callback(state)[source]

Call dense output after each step and store the interpolants.

perform_full_step(state, initial_dt)

Perform a full ODE solver step.

This includes the acceptance/rejection decision as governed by error estimation and steprule.

postprocess(odesol)[source]

Process the ODESolution object before returning.

rvlist_to_odesol(times, rvs)[source]

Create an ODESolution object.

solution_generator(ivp, stop_at=None, callbacks=None)

Generate ODE solver steps.

solve(ivp, stop_at=None, callbacks=None)

Solve an IVP.

Parameters