All posts
·7 min read

Implementing IMEX and Multirate Integrators in OrdinaryDiffEq.jl

  • Julia
  • SciML
  • Numerical Methods
  • ODE

A big chunk of my OrdinaryDiffEq.jl work has been adding IMEX and multirate integrators. Both come from the same observation: a single system often contains terms with wildly different character, and forcing one integration strategy on all of them is wasteful.

IMEX: implicit where you must, explicit where you can

In many PDE discretizations you can split the right-hand side:

u' = f_explicit(u) + f_implicit(u)

The stiff part (say, diffusion) needs an implicit solve for stability; the non-stiff part (advection, reactions) is cheaper done explicitly. IMEX Runge–Kutta methods do exactly that in one scheme.

I implemented several families:

  • ARS schemes (ARS222, ARS232, ARS443) - Ascher–Ruuth–Spiteri.
  • BHR(5,5,3) - Boscarino–Russo.
  • The IMEX-SSP (strong-stability-preserving) family.

The unifying idea was an ESDIRKIMEXTableau abstraction so all of these share one generic perform_step! instead of each method hand-rolling its own stepping code.

# One generic stepper, driven by a tableau, instead of N copies
function perform_step!(integrator, cache::ESDIRKIMEXCache)
    # explicit stages + implicit Newton solves, read from the tableau
end

Multirate: different clocks for different terms

Multirate methods go further: when part of your system evolves much faster than the rest, you take small steps for the fast part and large steps for the slow part. I added:

  • MRI-GARK (Multirate Infinitesimal Generalized-structure ARK) - ERK22a/22b, Sandu 2019.
  • MIS2 - Multirate Infinitesimal Step (Wensch–Knoth–Galant).
  • MRAB - a multirate Adams–Bashforth method.
  • MREEF - a multirate extrapolated explicit Euler.

The fun (and the hard part) is the slow–fast coupling: the slow method produces a "forcing" function that the fast integrator sees as it sub-steps. Get the accumulation of that forcing wrong and you silently lose an order of accuracy - which only a careful convergence-order test catches.

Why this matters

These methods aren't academic curiosities - they're what make large stiff simulations (climate, combustion, circuits) tractable. Contributing them taught me to:

  • Read papers as specifications. A Butcher tableau in a 2009 paper becomes a tested method in a library.
  • Trust convergence tests over intuition. "It looks right" is not a proof; an order-of-accuracy plot is.
  • Design for the next method. The best PR isn't one new solver - it's an abstraction that makes the next ten solvers easy.