OrdinaryDiffEqBDF

IMEX BDF (Implicit-Explicit Backward Differentiation Formula) methods for stiff differential equations that can be split into stiff and non-stiff components. These methods apply implicit BDF schemes to the stiff part while treating the non-stiff part explicitly, providing efficient handling of problems with mixed stiffness characteristics.

Key Properties

IMEX BDF methods provide:

Implicit-explicit splitting for mixed stiffness problems
BDF stability for the stiff component with A-stable and L-stable behavior
Explicit treatment of non-stiff terms avoiding unnecessary computational cost
High-order accuracy up to 4th order for both components
Efficient for large systems where full implicit treatment is expensive
Natural for operator splitting problems

When to Use IMEX BDF Methods

These methods are recommended for:

Reaction-diffusion systems where reaction terms are stiff and diffusion is moderate
Convection-diffusion problems with stiff source terms and explicit convection
Parabolic PDEs where diffusion operators are naturally split from other terms
Problems with natural stiffness separation where some terms require implicit treatment
Large-scale systems where full implicit methods are computationally prohibitive
Applications requiring operator splitting methodology

Mathematical Background

IMEX BDF methods split the ODE system du/dt = f(u,t) into: du/dt = f₁(u,t) + f₂(u,t)

where:

f₁(u,t) contains stiff terms (treated implicitly with BDF)
f₂(u,t) contains non-stiff terms (treated explicitly)

This splitting must be chosen carefully to ensure both stability and efficiency.

Problem Splitting Requirements

These methods require a SplitODEProblem formulation where:

First functionf₁ should contain stiff, implicit terms
Second functionf₂ should contain non-stiff, explicit terms
Splitting strategy significantly affects method performance
Stiffness characteristics should align with implicit/explicit treatment

Solver Selection Guide

IMEX Multistep Methods

SBDF2: Recommended - Second-order IMEX BDF method, good balance of accuracy and stability
SBDF3: Third-order method for higher accuracy requirements
SBDF4: Fourth-order method for maximum accuracy in IMEX BDF family
SBDF: Adaptive order method (experimental)

IMEX SDIRK Methods

IMEXEuler: First-order method for simple problems or debugging
IMEXEulerARK: Alternative first-order formulation

Performance Guidelines

When IMEX BDF methods excel

Natural stiffness separation where splitting is obvious
Large systems where full implicit treatment is expensive
Parabolic PDEs with natural operator splitting
Reaction-diffusion problems with well-separated timescales
Problems where implicit component has efficient linear algebra

Splitting strategy considerations

Identify stiff vs non-stiff terms based on eigenvalue analysis
Linear stiff terms work well in implicit component
Nonlinear non-stiff terms are suitable for explicit treatment
Test different splittings to optimize performance

Alternative Approaches

Consider these alternatives:

Full implicit methods (BDF, SDIRK) if splitting is unclear or ineffective
Standard IMEX Runge-Kutta methods for different accuracy/efficiency trade-offs
Exponential integrators for linear stiff problems with nonlinear non-stiff terms
Rosenbrock methods for moderately stiff problems without natural splitting

Usage Considerations

Careful splitting design is crucial for method effectiveness
Stability analysis should verify that explicit treatment doesn't introduce instabilities
Timestep restrictions may apply to the explicit component
Linear algebra efficiency in the implicit component affects overall performance

Installation

To be able to access the solvers in OrdinaryDiffEqBDF, you must first install them use the Julia package manager:

using Pkg
Pkg.add("OrdinaryDiffEqBDF")

This will only install the solvers listed at the bottom of this page. If you want to explore other solvers for your problem, you will need to install some of the other libraries listed in the navigation bar on the left.

Example usage

using OrdinaryDiffEqBDF

function lorenz!(du, u, p, t)
    du[1] = 10.0 * (u[2] - u[1])
    du[2] = u[1] * (28.0 - u[3]) - u[2]
    du[3] = u[1] * u[2] - (8 / 3) * u[3]
end
u0 = [1.0; 0.0; 0.0]
tspan = (0.0, 100.0)
prob = ODEProblem(lorenz!, u0, tspan)
sol = solve(prob, SBDF2())

Full list of solvers

IMEX Multistep

OrdinaryDiffEqBDF.SBDF — Type

SBDF(; chunk_size = Val{0}(),
       autodiff = AutoForwardDiff(),
       standardtag = Val{true}(),
       concrete_jac = nothing,
       linsolve = nothing,
       precs = DEFAULT_PRECS,
       κ = nothing,
       tol = nothing,
       nlsolve = NLNewton(),
       extrapolant = :linear,
       ark = false,
       order)

Multistep Method. Implicit-explicit (IMEX) method designed for SplitODEFunction equations, which reduce the size of the implicit handling to a subset of the equations. It's similar to the additive Runge-Kutta methods in splitting mode, like KenCarp4, but instead using a multistep BDF approach

Keyword Arguments

autodiff: Uses ADTypes.jl to specify whether to use automatic differentiation via ForwardDiff.jl or finite differencing via FiniteDiff.jl. Defaults to AutoForwardDiff() for automatic differentiation, which by default uses chunksize = 0, and thus uses the internal ForwardDiff.jl algorithm for the choice. To use FiniteDiff.jl, the AutoFiniteDiff() ADType can be used, which has a keyword argument fdtype with default value Val{:forward}(), and alternatives Val{:central}() and Val{:complex}().
standardtag: Specifies whether to use package-specific tags instead of the ForwardDiff default function-specific tags. For more information, see this blog post. Defaults to Val{true}().
concrete_jac: Specifies whether a Jacobian should be constructed. Defaults to nothing, which means it will be chosen true/false depending on circumstances of the solver, such as whether a Krylov subspace method is used for linsolve.
linsolve: Any LinearSolve.jl compatible linear solver. For example, to use KLU.jl, specify SBDF(linsolve = KLUFactorization()). When nothing is passed, uses DefaultLinearSolver.
precs: Any LinearSolve.jl-compatible preconditioner can be used as a left or right preconditioner. Preconditioners are specified by the Pl,Pr = precs(W,du,u,p,t,newW,Plprev,Prprev,solverdata) function where the arguments are defined as:
- W: the current Jacobian of the nonlinear system. Specified as either $I - \gamma J$ or $I/\gamma - J$ depending on the algorithm. This will commonly be a WOperator type defined by OrdinaryDiffEq.jl. It is a lazy representation of the operator. Users can construct the W-matrix on demand by calling convert(AbstractMatrix,W) to receive an AbstractMatrix matching the jac_prototype.
- du: the current ODE derivative
- u: the current ODE state
- p: the ODE parameters
- t: the current ODE time
- newW: a Bool which specifies whether the W matrix has been updated since the last call to precs. It is recommended that this is checked to only update the preconditioner when newW == true.
- Plprev: the previous Pl.
- Prprev: the previous Pr.
- solverdata: Optional extra data the solvers can give to the precs function. Solver-dependent and subject to change.
The return is a tuple (Pl,Pr) of the LinearSolve.jl-compatible preconditioners. To specify one-sided preconditioning, simply return nothing for the preconditioner which is not used. Additionally, precs must supply the dispatch:
```
Pl, Pr = precs(W, du, u, p, t, ::Nothing, ::Nothing, ::Nothing, solverdata)
```
which is used in the solver setup phase to construct the integrator type with the preconditioners (Pl,Pr). The default is precs=DEFAULT_PRECS where the default preconditioner function is defined as:
```
DEFAULT_PRECS(W, du, u, p, t, newW, Plprev, Prprev, solverdata) = nothing, nothing
```
/n- κ: TBD
tol: TBD
nlsolve: TBD
extrapolant: TBD
ark: TBD

References

@article{ascher1995implicit, title={Implicit-explicit methods for time-dependent partial differential equations}, author={Ascher, Uri M and Ruuth, Steven J and Wetton, Brian TR}, journal={SIAM Journal on Numerical Analysis}, volume={32}, number={3}, pages={797–823}, year={1995}, publisher={SIAM}}

OrdinaryDiffEqBDF.SBDF2 — Function

SBDF2(;kwargs...)

The two-step version of the IMEX multistep methods of

Uri M. Ascher, Steven J. Ruuth, Brian T. R. Wetton. Implicit-Explicit Methods for Time-Dependent Partial Differential Equations. Society for Industrial and Applied Mathematics. Journal on Numerical Analysis, 32(3), pp 797-823, 1995. doi: https://doi.org/10.1137/0732037

IMEX SDIRK

Note that Implicit Euler is the 1st order BDF method, and is thus implemented here using the same machinery.

OrdinaryDiffEqBDF.IMEXEuler — Function

IMEXEuler(;kwargs...)

The one-step version of the IMEX multistep methods of

Uri M. Ascher, Steven J. Ruuth, Brian T. R. Wetton. Implicit-Explicit Methods for Time-Dependent Partial Differential Equations. Society for Industrial and Applied Mathematics. Journal on Numerical Analysis, 32(3), pp 797-823, 1995. doi: https://doi.org/10.1137/0732037

When applied to a SplitODEProblem of the form

u'(t) = f1(u) + f2(u)

The default IMEXEuler() method uses an update of the form

unew = uold + dt * (f1(unew) + f2(uold))