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 ofperturb_lognormal()
orperturb_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.
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
state (
ODESolverState
) – Current state of the ODE solver.dt (
Real
) – Step-size.
- Returns
New state.
- Return type
_odesolver_state.ODESolverState
- 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.
- solution_generator(ivp, stop_at=None, callbacks=None)¶
Generate ODE solver steps.
- solve(ivp, stop_at=None, callbacks=None)¶
Solve an IVP.
- Parameters
ivp (
InitialValueProblem
) – Initial value problem.stop_at (
Optional
[Iterable
[Real
]]) – Time-points through which the solver must step. Optional. Default is None.callbacks (
Union
[ODESolverCallback
,Iterable
[ODESolverCallback
],None
]) – Callbacks to happen after every accepted step.