GitHub

OrdinaryDiffEqStabilizedRK

Stabilized Runge-Kutta methods are explicit schemes designed to handle moderately stiff problems by extending the stability region through careful tableau construction. These methods use an upper bound on the spectral radius of the Jacobian to achieve much larger stable timesteps than conventional explicit methods. These methods are good for large real eigenvalue problems, but not for problems where the complex eigenvalues are large.

Key Properties

Stabilized RK methods provide:

  • Extended stability regions for moderately stiff problems
  • Explicit formulation avoiding nonlinear solvers
  • Large stable timestep sizes compared to standard explicit methods
  • Automatic spectral radius estimation or user-supplied bounds
  • Efficient for parabolic PDEs with moderate stiffness
  • Good performance on problems with well-separated timescales

When to Use Stabilized RK Methods

These methods are recommended for:

  • Moderately stiff problems where implicit methods are overkill
  • Parabolic PDEs with diffusion-dominated behavior
  • Problems with large spatial grids where implicit methods become expensive
  • Systems with well-separated timescales but not extreme stiffness
  • Cases where explicit is preferred but standard methods are unstable
  • Large-scale problems where linear algebra cost of implicit methods is prohibitive

Mathematical Background

Stabilized methods achieve extended stability by constructing tableaus with enlarged stability regions, often using Chebyshev polynomials or orthogonal polynomial techniques. The stable timestep is determined by the spectral radius bound rather than the CFL condition. Important: These methods extend stability primarily along the negative real axis, making them effective for large real eigenvalues but ineffective when complex eigenvalues dominate the stiffness.

Spectral Radius Estimation

Users can supply an upper bound on the spectral radius using:

eigen_est = (integrator) -> integrator.eigen_est = upper_bound

If not provided, the methods include automatic estimation procedures.

Solver Selection Guide

  • ROCK2: Second-order ROW-type stabilized method with extended stability
  • ROCK4: Fourth-order stabilized method for higher accuracy requirements

Performance Guidelines

When stabilized methods excel

  • Large real eigenvalue problems where stiffness comes from real eigenvalues
  • Moderate stiffness ratio (10³ to 10⁶) dominated by real eigenvalues
  • Large spatial discretizations where implicit solver cost is high
  • Very large systems where stabilized RK methods are more efficient than BDF methods due to no linear algebra requirements
  • Parabolic PDEs with diffusion-dominated (real eigenvalue) stiffness
  • Problems where spectral radius can be estimated reliably

When to use alternatives

  • Complex eigenvalue dominated problems: Use implicit methods (BDF, SDIRK, Rosenbrock)
  • Non-stiff problems: Use standard explicit methods (Tsit5, Verner)

Usage Considerations

  • Spectral radius estimation is crucial for performance
  • Method efficiency depends on stiffness ratio
  • Test against implicit methods for highly stiff problems
  • Consider adaptive spectral radius estimation for varying stiffness

Installation

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

using Pkg
Pkg.add("OrdinaryDiffEqStabilizedRK")

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 OrdinaryDiffEqStabilizedRK

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

Full list of solvers

OrdinaryDiffEqStabilizedRK.ROCK2 — Type
ROCK2(; min_stages = 0,
        max_stages = 200,
        eigen_est = nothing)

Stabilized Explicit Method. High stability for real eigenvalues. Second order method. Exhibits high stability for real eigenvalues and is smoothened to allow for moderate sized complex eigenvalues.

Keyword Arguments

  • min_stages: The minimum degree of the Chebyshev polynomial.
  • max_stages: The maximumdegree of the Chebyshev polynomial.
  • eigen_est: function of the form (integrator) -> integrator.eigen_est = upper_bound, where upper_bound is an estimated upper bound on the spectral radius of the Jacobian matrix. If eigen_est is not provided, upper_bound will be estimated using the power iteration.

References

Assyr Abdulle, Alexei A. Medovikov. Second Order Chebyshev Methods based on Orthogonal Polynomials. Numerische Mathematik, 90 (1), pp 1-18, 2001. doi: https://dx.doi.org/10.1007/s002110100292

source
OrdinaryDiffEqStabilizedRK.ROCK4 — Type
ROCK4(; min_stages = 0,
        max_stages = 152,
        eigen_est = nothing)

Stabilized Explicit Method. High stability for real eigenvalues. Fourth order method. Exhibits high stability for real eigenvalues and is smoothened to allow for moderate sized complex eigenvalues.

Keyword Arguments

  • min_stages: The minimum degree of the Chebyshev polynomial.
  • max_stages: The maximumdegree of the Chebyshev polynomial.
  • eigen_est: function of the form (integrator) -> integrator.eigen_est = upper_bound, where upper_bound is an estimated upper bound on the spectral radius of the Jacobian matrix. If eigen_est is not provided, upper_bound will be estimated using the power iteration.

References

Assyr Abdulle. Fourth Order Chebyshev Methods With Recurrence Relation. 2002 Society for Industrial and Applied Mathematics Journal on Scientific Computing, 23(6), pp 2041-2054, 2001. doi: https://doi.org/10.1137/S1064827500379549

source
OrdinaryDiffEqStabilizedRK.RKC — Type
RKC(; eigen_est = nothing)

Stabilized Explicit Method. Second order method. Exhibits high stability for real eigenvalues.

Keyword Arguments

  • eigen_est: function of the form (integrator) -> integrator.eigen_est = upper_bound, where upper_bound is an estimated upper bound on the spectral radius of the Jacobian matrix. If eigen_est is not provided, upper_bound will be estimated using the power iteration.

References

B. P. Sommeijer, L. F. Shampine, J. G. Verwer. RKC: An Explicit Solver for Parabolic PDEs, Journal of Computational and Applied Mathematics, 88(2), pp 315-326, 1998. doi: https://doi.org/10.1016/S0377-0427(97)00219-7

source
OrdinaryDiffEqStabilizedRK.SERK2 — Type
SERK2(; controller = :PI
        eigen_est = nothing)

Stabilized Explicit Method. Second order method.

Keyword Arguments

  • controller: TBD
  • eigen_est: function of the form (integrator) -> integrator.eigen_est = upper_bound, where upper_bound is an estimated upper bound on the spectral radius of the Jacobian matrix. If eigen_est is not provided, upper_bound will be estimated using the power iteration.

References

@article{kleefeld2013serk2v2, title={SERK2v2: A new second-order stabilized explicit Runge-Kutta method for stiff problems}, author={Kleefeld, B and Martin-Vaquero, J}, journal={Numerical Methods for Partial Differential Equations}, volume={29}, number={1}, pages={170–185}, year={2013}, publisher={Wiley Online Library}}

source
OrdinaryDiffEqStabilizedRK.ESERK4 — Type
ESERK4(; eigen_est = nothing)

Stabilized Explicit Method. Fourth order method. Exhibits high stability for real eigenvalues and is smoothened to allow for moderate sized complex eigenvalues.

Keyword Arguments

  • eigen_est: function of the form (integrator) -> integrator.eigen_est = upper_bound, where upper_bound is an estimated upper bound on the spectral radius of the Jacobian matrix. If eigen_est is not provided, upper_bound will be estimated using the power iteration.

References

J. Martín-Vaquero, B. Kleefeld. Extrapolated stabilized explicit Runge-Kutta methods, Journal of Computational Physics, 326, pp 141-155, 2016. doi: https://doi.org/10.1016/j.jcp.2016.08.042.

source
OrdinaryDiffEqStabilizedRK.ESERK5 — Type
ESERK5(; eigen_est = nothing)

Stabilized Explicit Method. Fifth order method. Exhibits high stability for real eigenvalues and is smoothened to allow for moderate sized complex eigenvalues.

Keyword Arguments

  • eigen_est: function of the form (integrator) -> integrator.eigen_est = upper_bound, where upper_bound is an estimated upper bound on the spectral radius of the Jacobian matrix. If eigen_est is not provided, upper_bound will be estimated using the power iteration.

References

J. Martín-Vaquero, A. Kleefeld. ESERK5: A fifth-order extrapolated stabilized explicit Runge-Kutta method, Journal of Computational and Applied Mathematics, 356, pp 22-36, 2019. doi: https://doi.org/10.1016/j.cam.2019.01.040.

source