OrdinaryDiffEqStabilizedIRK

Stabilized Implicit Runge-Kutta (IMEX) methods combine the benefits of stabilized explicit methods with implicit treatment of problematic eigenvalues. These IMEX schemes are designed for problems where the Jacobian has both large real eigenvalues (suitable for explicit stabilized methods) and large complex eigenvalues (requiring implicit treatment).

Key Properties

Stabilized IRK methods provide:

  • IMEX formulation treating different stiffness components appropriately
  • Large stability regions for real eigenvalues via explicit stabilized schemes
  • Implicit treatment of complex eigenvalues for unconditional stability
  • Efficient handling of mixed stiffness characteristics
  • Splitting-based approach requiring SplitODEProblem formulation

When to Use Stabilized IRK Methods

These methods are recommended for:

  • Mixed stiffness problems with both real and complex eigenvalues
  • Parabolic PDEs with convection where diffusion and advection have different scales
  • Reaction-diffusion systems with stiff reactions and moderate diffusion
  • Problems where pure explicit stabilized methods fail due to complex eigenvalues
  • Large-scale systems where full implicit methods are too expensive

Mathematical Background

Standard stabilized explicit methods (like RKC, ROCK) achieve large stability regions along the negative real axis but struggle with complex eigenvalues. Stabilized IRK methods address this by:

  1. Explicit stabilized treatment for large real eigenvalues
  2. Implicit treatment for complex eigenvalues
  3. IMEX coupling to maintain overall stability and accuracy

Problem Splitting Requirements

These methods require a SplitODEProblem where:

  • First component contains terms with large real eigenvalues (explicit treatment)
  • Second component contains terms with complex eigenvalues (implicit treatment)
  • Splitting design is crucial for method performance

Spectral Radius Estimation

Users can supply an upper bound on the spectral radius:

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

This bound applies to the explicit component of the split problem.

Solver Selection Guide

Available methods

  • IRKC: Implicit Runge-Kutta-Chebyshev method for mixed stiffness problems

Usage considerations

  • Requires careful splitting of the problem components
  • Spectral radius estimation needed for explicit component
  • Test splitting strategies for optimal performance
  • Compare with pure implicit or explicit stabilized alternatives

Performance Guidelines

When IMEX stabilized methods excel

  • Mixed eigenvalue distribution (both real and complex)
  • Moderate to large systems where splitting is natural
  • Problems where neither pure explicit nor implicit methods are ideal

Splitting strategy considerations

  • Identify dominant eigenvalue types in different terms
  • Real-dominated terms → explicit component
  • Complex-dominated terms → implicit component
  • Test different splittings for best performance

Alternative Approaches

Consider these alternatives:

  • Pure implicit methods (BDF, SDIRK, Rosenbrock) for highly stiff problems
  • Explicit stabilized methods (ROCK, RKC) if complex eigenvalues are small
  • Standard IMEX methods for natural explicit/implicit splitting

Installation

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

using Pkg
Pkg.add("OrdinaryDiffEqStabilizedIRK")

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 OrdinaryDiffEqStabilizedIRK

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

Full list of solvers

OrdinaryDiffEqStabilizedIRK.IRKCType
IRKC(; eigen_est = nothing)

Stabilized Implicit Runge Kutta method. Implicit Runge-Kutta-Chebyshev method.

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

REF TBD

source