OrdinaryDiffEqRosenbrock

Rosenbrock methods are semi-implicit Runge-Kutta methods designed for small to medium-sized stiff systems of differential equations. These methods combine the stability properties of implicit methods with computational efficiency by using approximate Jacobians. They work particularly well for strict tolerance requirements.

Key Properties

Rosenbrock methods provide:

Excellent efficiency for small to medium systems (<1000 ODEs)
L-stable and A-stable variants for stiff problems
W-method structure making them robust to inaccurate Jacobians
Automatic differentiation compatibility for Jacobian computation
High-order accuracy with embedded error estimation
Excellent performance for strict tolerance requirements
Stiffly accurate variants for enhanced stability

When to Use Rosenbrock Methods

Rosenbrock methods are recommended for:

Small to medium stiff systems (10 to 1000 equations)
Problems where Jacobians are available or can be computed efficiently
Medium tolerance requirements (1e-8 to 1e-2)
Stiff ODEs arising from reaction kinetics and chemical systems
Moderately stiff PDEs after spatial discretization
Problems requiring reliable error control for stiff systems

Solver Selection Guide

Low tolerance (>1e-2)

Rosenbrock23: Second/third-order method, recommended for low tolerance requirements

Medium tolerance (1e-8 to 1e-2)

Rodas5P: Fifth-order method, most efficient for many problems
Rodas4P: Fourth-order method, more reliable than Rodas5P
Rodas5Pe: Enhanced fifth-order variant with stiffly accurate embedded estimate for better adaptivity on highly stiff equations
Rodas5Pr: Fifth-order variant with residual test for robust error estimation that guarantees accuracy on interpolation

Performance Guidelines

Rodas5P: Best overall efficiency at medium tolerances
Rodas4P: Most reliable when Rodas5P fails
Rosenbrock23: Fastest at high tolerances (>1e-2)

When to Choose Alternatives

Consider other methods when:

System size > 1000: Use BDF methods (QNDF, FBDF)
Matrix-free methods: If using Krylov solvers for the linear solver (matrix-free methods), SDIRK or BDF methods are preferred

Advantages of Rosenbrock Methods

Very low tolerances: Rosenbrock23 performs well even at very low tolerances
Very high stiffness: Rosenbrock methods are often more stable than SDIRK or BDF methods because other methods can diverge due to bad initial guesses for the Newton method (leading to nonlinear solver divergence)

Installation

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

using Pkg
Pkg.add("OrdinaryDiffEqRosenbrock")

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 OrdinaryDiffEqRosenbrock

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, Rodas5P())

Full list of solvers

OrdinaryDiffEqRosenbrock.Rosenbrock23 — Type

Rosenbrock23(; chunk_size = Val{0}(),
               standardtag = Val{true}(),
               autodiff = AutoForwardDiff(),
               concrete_jac = nothing,
               diff_type = Val{:forward}(),
               linsolve = nothing,
               precs = DEFAULT_PRECS,
               step_limiter! = OrdinaryDiffEq.trivial_limiter!)

Rosenbrock-Wanner-W(olfbrandt) Method. An Order 2/3 L-Stable Rosenbrock-W method which is good for very stiff equations with oscillations at low tolerances. 2nd order stiff-aware interpolation.

Keyword Arguments

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}().
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}().
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 Rosenbrock23(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
```
step_limiter!: function of the form limiter!(u, integrator, p, t)

References

Shampine L.F. and Reichelt M., (1997) The MATLAB ODE Suite, SIAM Journal of

Scientific Computing, 18 (1), pp. 1-22.