# Solvers

ProbNumDiffEq.jl provides two solvers: the `EK1`

and the `EK0`

. Both based on extended Kalman filtering and smoothing, but the latter relies on evaluating the Jacobian of the vector field.

**Which solver should I use?**

- Use the
`EK1`

to get the best uncertainty quantification and to solve stiff problems. - Use the
`EK0`

to get the fastest runtimes and to solve high-dimensional problems.

All solvers are compatible with DAEs in mass-matrix ODE form. They also specialize on second-order ODEs: If the problem is of type `SecondOrderODEProblem`

, it solves the second-order problem directly; this is more efficient than solving the transformed first-order problem and provides more meaningful posteriors [1].

## API

`ProbNumDiffEq.EK1`

— Type```
EK1(; order=3,
smooth=true,
prior=IWP(order),
diffusionmodel=DynamicDiffusion(),
initialization=TaylorModeInit(num_derivatives(prior)),
kwargs...)
```

**Gaussian ODE filter with first-order vector field linearization.**

This is a *semi-implicit*, L-stable ODE solver so it can handle stiffness quite well [4], and it generally produces more expressive posterior covariances than the `EK0`

. However, as typical implicit ODE solvers it scales cubically with the ODE dimension [5], so if you're solving a high-dimensional non-stiff problem you might want to give the `EK0`

a try.

**Arguments**

`order::Integer`

: Order of the integrated Wiener process (IWP) prior.`smooth::Bool`

: Turn smoothing on/off; smoothing is required for dense output.`prior::AbstractGaussMarkovProcess`

: Prior to be used by the ODE filter. By default, uses a 3-times integrated Wiener process prior`IWP(3)`

. See also: Priors.`diffusionmodel::ProbNumDiffEq.AbstractDiffusion`

: See Diffusion models and calibration.`initialization::ProbNumDiffEq.InitializationScheme`

: See Initialization.

Some additional `kwargs`

relating to implicit solvers are supported; check out DifferentialEquations.jl's Extra Options page. Right now, we support `autodiff`

, `chunk_size`

, and `diff_type`

. In particular, `autodiff=false`

can come in handy to use finite differences instead of ForwardDiff.jl to compute Jacobians.

**Examples**

`julia> solve(prob, EK1())`

`ProbNumDiffEq.EK0`

— Type```
EK0(; order=3,
smooth=true,
prior=IWP(order),
diffusionmodel=DynamicDiffusion(),
initialization=TaylorModeInit(num_derivatives(prior)))
```

**Gaussian ODE filter with zeroth-order vector field linearization.**

This is an *explicit* ODE solver. It is fast and scales well to high-dimensional problems [5], but it is not L-stable [4]. So for stiff problems, use the `EK1`

.

Whenever possible this solver will use a Kronecker-factored implementation to achieve its linear scaling and to get the best runtimes. This can currently be done only with an `IWP`

prior (default), with a scalar diffusion model (either `DynamicDiffusion`

or `FixedDiffusion`

). *For other configurations the solver falls back to a dense implementation which scales cubically with the problem size.*

**Arguments**

`order::Integer`

: Order of the integrated Wiener process (IWP) prior.`smooth::Bool`

: Turn smoothing on/off; smoothing is required for dense output.`prior::AbstractGaussMarkovProcess`

: Prior to be used by the ODE filter. By default, uses a 3-times integrated Wiener process prior`IWP(3)`

. See also: Priors.`diffusionmodel::ProbNumDiffEq.AbstractDiffusion`

: See Diffusion models and calibration.`initialization::ProbNumDiffEq.InitializationScheme`

: See Initialization.

**Examples**

`julia> solve(prob, EK0())`

### Probabilistic Exponential Integrators

`ProbNumDiffEq.ExpEK`

— Function`ExpEK(; L, order=3, kwargs...)`

**Probabilistic exponential integrator**

Probabilistic exponential integrators are a class of integrators for semi-linear stiff ODEs that provide improved stability by essentially solving the linear part of the ODE exactly. In probabilistic numerics, this amounts to including the linear part into the prior model of the solver.

`ExpEK`

is therefore just a short-hand for `EK0`

with `IOUP`

prior:

`ExpEK(; order=3, L, kwargs...) = EK0(; prior=IOUP(order, L), kwargs...)`

See also `RosenbrockExpEK`

, `EK0`

, `EK1`

.

**Arguments**

See `EK0`

for available keyword arguments.

**Examples**

```
julia> prob = ODEProblem((du, u, p, t) -> (@. du = - u + sin(u)), [1.0], (0.0, 10.0))
julia> solve(prob, ExpEK(L=-1))
```

**Reference**

- [2] Bosch et al, "Probabilistic Exponential Integrators", arXiv (2021)

`ProbNumDiffEq.RosenbrockExpEK`

— Function`RosenbrockExpEK(; order=3, kwargs...)`

**Probabilistic Rosenbrock-type exponential integrator**

A probabilistic exponential integrator similar to `ExpEK`

, but with automatic linearization along the mean numerical solution. This brings the advantage that the linearity does not need to be specified manually, and the more accurate local linearization can sometimes also improve stability; but since the "prior" is adjusted at each step the probabilistic interpretation becomes more complicated.

`RosenbrockExpEK`

is just a short-hand for `EK1`

with locally-updated `IOUP`

prior:

`RosenbrockExpEK(; order=3, kwargs...) = EK1(; prior=IOUP(order, update_rate_parameter=true), kwargs...)`

**Arguments**

See `EK1`

for available keyword arguments.

**Examples**

```
julia> prob = ODEProblem((du, u, p, t) -> (@. du = - u + sin(u)), [1.0], (0.0, 10.0))
julia> solve(prob, RosenbrockExpEK())
```

**Reference**

- [2] Bosch et al, "Probabilistic Exponential Integrators", NeurIPS (2023)

## References

- [4]
- F. Tronarp, H. Kersting, S. Särkkä and P. Hennig.
*Probabilistic solutions to ordinary differential equations as nonlinear Bayesian filtering: a new perspective*. Statistics and Computing**29**, 1297–1315 (2019). - [5]
- N. Krämer, N. Bosch, J. Schmidt and P. Hennig.
*Probabilistic ODE Solutions in Millions of Dimensions*. In:*Proceedings of the 39th International Conference on Machine Learning*, Vol. 162 of*Proceedings of Machine Learning Research*, edited by K. Chaudhuri, S. Jegelka, L. Song, C. Szepesvari, G. Niu and S. Sabato (PMLR, 17–23 Jul 2022); pp. 11634–11649. - [2]
- N. Bosch, P. Hennig and F. Tronarp.
*Probabilistic Exponential Integrators*. In:*Thirty-seventh Conference on Neural Information Processing Systems*(2023).