Implicit Methods for Stiff SDEs

When SDEs have stiff drift terms, explicit methods may require impractically small time steps for stability. Implicit methods treat the drift term implicitly while keeping the diffusion explicit, providing excellent stability properties.

StochasticDiffEq.SKenCarp — Type
SKenCarp(;chunk_size=0, autodiff=true, diff_type=Val{:central}, 
         standardtag=Val{true}(), concrete_jac=nothing, precs=DEFAULT_PRECS,
         linsolve=nothing, nlsolve=NLNewton(), smooth_est=true, 
         extrapolant=:min_correct, new_jac_conv_bound=1e-3, 
         controller=:Predictive, ode_error_est=true)

SKenCarp: Stochastic KenCarp Method (Stiff) - Highly Recommended for Stiff Problems

Adaptive L-stable drift-implicit method with strong order 1.5. Highly recommended for stiff problems with additive noise.

Method Properties

  • Strong Order: 1.5 (for additive noise)
  • Weak Order: 2.0
  • Time stepping: Adaptive
  • Noise types: Additive noise (diagonal, non-diagonal, and scalar)
  • SDE interpretation: Both Itô and Stratonovich
  • Stability: L-stable (excellent for stiff problems)
  • Implicit: Drift-implicit (handles stiffness in drift term)

When to Use

  • Highly recommended for stiff additive noise problems
  • When the drift term f(u,p,t) is stiff
  • For problems requiring high accuracy with stiff dynamics
  • When implicit treatment of the drift is necessary for stability
  • Best choice for stiff problems with additive noise structure

Algorithm Description

SKenCarp applies implicit treatment to the drift term while keeping the diffusion explicit. This provides excellent stability for stiff SDEs with additive noise.

Stiffness and Stability

  • L-stable: Excellent for stiff problems
  • Handles large negative eigenvalues in the drift term
  • Maintains accuracy while providing stability

Configuration Options

  • Linear solver options via linsolve parameter
  • Nonlinear solver options via nlsolve parameter
  • Jacobian computation control via autodiff and related parameters
  • Step size control via controller parameter

References

  • Based on KenCarp methods from OrdinaryDiffEq.jl
  • Adapted for stochastic problems with additive noise

Basic Implicit Methods

ImplicitEM - Implicit Euler-Maruyama

StochasticDiffEq.ImplicitEM — Type
ImplicitEM(;chunk_size=0, autodiff=true, diff_type=Val{:central},
           standardtag=Val{true}(), concrete_jac=nothing, precs=DEFAULT_PRECS,
           linsolve=nothing, nlsolve=NLNewton(), extrapolant=:constant,
           theta=1, symplectic=false, new_jac_conv_bound=1e-3, 
           controller=:Predictive)

ImplicitEM: Implicit Euler-Maruyama Method (Stiff)

Drift-implicit version of the Euler-Maruyama method with theta-method treatment of the drift term.

Method Properties

  • Strong Order: 0.5 (Itô sense)
  • Weak Order: 1.0
  • Time stepping: Adaptive (1.0/1.5 heuristic)
  • Noise types: All forms (non-diagonal, scalar, colored noise)
  • SDE interpretation: Itô
  • Implicit treatment: Drift term only (diffusion remains explicit)

Parameters

  • theta::Real = 1: Implicitness parameter (0=explicit, 1=fully implicit, 0.5=trapezoidal)
  • symplectic::Bool = false: When true and theta=0.5, uses symplectic implicit midpoint
  • Linear/nonlinear solver options via linsolve and nlsolve

When to Use

  • For mildly stiff SDEs where drift term causes stability issues
  • When explicit methods require very small time steps
  • As a robust fallback for difficult problems
  • When all noise types need to be supported

Theta Method Variants

  • theta = 0: Explicit Euler (not recommended, use EM instead)
  • theta = 0.5: Trapezoidal rule (second order accurate for deterministic part)
  • theta = 1: Backward Euler (default, maximum stability)

Symplectic Option

When symplectic=true and theta=0.5, the method preserves the symplectic structure in distribution for appropriate problems.

Algorithm Description

Treats the SDE du = f(u,t)dt + g(u,t)dW using:

u_{n+1} = u_n + theta*f(u_{n+1},t_{n+1})*dt + (1-theta)*f(u_n,t_n)*dt + g(u_n,t_n)*dW_n

References

  • Standard implicit methods adapted for SDEs

ImplicitEulerHeun - Implicit Euler-Heun (Stratonovich)

StochasticDiffEq.ImplicitEulerHeun — Type
ImplicitEulerHeun(;chunk_size=0, autodiff=true, diff_type=Val{:central},
                  standardtag=Val{true}(), concrete_jac=nothing, precs=DEFAULT_PRECS,
                  linsolve=nothing, nlsolve=NLNewton(), extrapolant=:constant,
                  theta=1, symplectic=false, new_jac_conv_bound=1e-3, 
                  controller=:Predictive)

ImplicitEulerHeun: Implicit Euler-Heun Method (Stiff)

Drift-implicit version of the Euler-Heun method for Stratonovich SDEs with stiff drift terms.

Method Properties

  • Strong Order: 0.5 (Stratonovich sense)
  • Weak Order: 1.0
  • Time stepping: Adaptive (1.0/1.5 heuristic)
  • Noise types: All forms (non-diagonal, scalar, colored noise)
  • SDE interpretation: Stratonovich
  • Implicit treatment: Drift term only (diffusion remains explicit)

Parameters

  • theta::Real = 1: Implicitness parameter (0=explicit, 1=fully implicit, 0.5=trapezoidal)
  • symplectic::Bool = false: When true and theta=1, uses symplectic implicit midpoint
  • Linear/nonlinear solver options via linsolve and nlsolve

When to Use

  • Stiff Stratonovich SDEs where drift term causes stability issues
  • When working in Stratonovich interpretation with stiff dynamics
  • Alternative to ImplicitEM for Stratonovich problems
  • When all noise types need to be supported in Stratonovich form

Theta Method Variants

  • theta = 0.5: Trapezoidal rule (default, good accuracy/stability balance)
  • theta = 1: Backward Euler (maximum stability)

Symplectic Option

When symplectic=true and theta=1, preserves symplectic structure in distribution for appropriate Stratonovich problems.

References

  • Implicit methods for stiff SDEs in Stratonovich interpretation

ImplicitRKMil - Implicit Runge-Kutta Milstein

StochasticDiffEq.ImplicitRKMil — Type
ImplicitRKMil(;chunk_size=0, autodiff=true, diff_type=Val{:central},
              standardtag=Val{true}(), concrete_jac=nothing, precs=DEFAULT_PRECS,
              linsolve=nothing, nlsolve=NLNewton(), extrapolant=:constant,
              theta=1, symplectic=false, new_jac_conv_bound=1e-3, 
              controller=:Predictive, interpretation=AlgorithmInterpretation.Ito)

ImplicitRKMil: Implicit Runge-Kutta Milstein Method (Stiff)

Drift-implicit Runge-Kutta Milstein method achieving order 1.0 for stiff problems with diagonal/scalar noise.

Method Properties

  • Strong Order: 1.0
  • Weak Order: Depends on tableau
  • Time stepping: Adaptive (1.5/2.0 heuristic)
  • Noise types: Diagonal and scalar noise only
  • SDE interpretation: Configurable (Itô or Stratonovich)
  • Implicit treatment: Drift term only (diffusion remains explicit)

Parameters

  • theta::Real = 1: Implicitness parameter (0.5=trapezoidal, 1=backward Euler)
  • symplectic::Bool = false: When true and theta=0.5, uses symplectic implicit midpoint
  • interpretation: Choose AlgorithmInterpretation.Ito (default) or AlgorithmInterpretation.Stratonovich
  • Linear/nonlinear solver options via linsolve and nlsolve

When to Use

  • Stiff problems requiring higher accuracy than ImplicitEM
  • When strong order 1.0 is needed with implicit stability
  • Diagonal or scalar noise problems with stiff drift
  • Alternative to SKenCarp for non-additive noise

Restrictions

  • Only works with diagonal or scalar noise
  • For non-diagonal noise, use ISSEM/ISSEulerHeun
  • For additive noise, prefer SKenCarp

Algorithm Features

  • Higher order accuracy than ImplicitEM
  • Milstein correction for improved strong convergence
  • Configurable interpretation (Itô/Stratonovich)

References

  • Implicit Milstein methods for stiff SDEs

Split-Step Implicit Methods

ISSEM - Implicit Split-Step Euler-Maruyama

StochasticDiffEq.ISSEM — Type
ISSEM(;chunk_size=0, autodiff=true, diff_type=Val{:central},
      standardtag=Val{true}(), concrete_jac=nothing, precs=DEFAULT_PRECS,
      linsolve=nothing, nlsolve=NLNewton(), extrapolant=:constant,
      theta=1, symplectic=false, new_jac_conv_bound=1e-3, 
      controller=:Predictive)

ISSEM: Implicit Split-Step Euler-Maruyama Method (Stiff)

Fully implicit split-step method for handling stiffness in both drift and diffusion terms.

Method Properties

  • Strong Order: 0.5 (Itô sense)
  • Weak Order: 1.0
  • Time stepping: Adaptive (1.0/1.5 heuristic)
  • Noise types: All forms (non-diagonal, scalar, colored noise)
  • SDE interpretation: Itô
  • Implicit treatment: Both drift and diffusion terms (fully implicit)

Parameters

  • theta::Real = 1: Implicitness parameter for drift term
  • symplectic::Bool = false: When true and theta=0.5, uses symplectic implicit midpoint
  • Linear/nonlinear solver options via linsolve and nlsolve

When to Use

  • Recommended for stiff Itô problems with large noise terms
  • When both drift and diffusion cause stability issues
  • Problems where ImplicitEM and ImplicitRKMil are insufficient
  • Fully stiff SDEs requiring implicit treatment of everything

Algorithm Description

Applies implicit treatment to both drift and diffusion using split-step approach:

Step 1: Handle drift implicitly
Step 2: Handle diffusion implicitly

Fully Implicit Features

  • Can handle stiffness in both drift and diffusion
  • More expensive than drift-only implicit methods
  • Most robust for extremely stiff problems
  • Requires solving nonlinear systems for both terms

References

  • Split-step implicit methods for fully stiff SDEs

ISSEulerHeun - Implicit Split-Step Euler-Heun

StochasticDiffEq.ISSEulerHeun — Type
ISSEulerHeun(;chunk_size=0, autodiff=true, diff_type=Val{:central},
             standardtag=Val{true}(), concrete_jac=nothing, precs=DEFAULT_PRECS,
             linsolve=nothing, nlsolve=NLNewton(), extrapolant=:constant,
             theta=1, symplectic=false, new_jac_conv_bound=1e-3, 
             controller=:Predictive)

ISSEulerHeun: Implicit Split-Step Euler-Heun Method (Stiff)

Fully implicit split-step method for Stratonovich SDEs with stiffness in both drift and diffusion terms.

Method Properties

  • Strong Order: 0.5 (Stratonovich sense)
  • Weak Order: 1.0
  • Time stepping: Adaptive (1.0/1.5 heuristic)
  • Noise types: All forms (non-diagonal, scalar, colored noise)
  • SDE interpretation: Stratonovich
  • Implicit treatment: Both drift and diffusion terms (fully implicit)

Parameters

  • theta::Real = 1: Implicitness parameter for drift term
  • symplectic::Bool = false: When true and theta=0.5, uses symplectic implicit midpoint
  • Linear/nonlinear solver options via linsolve and nlsolve

When to Use

  • Recommended for stiff Stratonovich problems with large noise terms
  • When both drift and diffusion cause stability issues in Stratonovich form
  • Stratonovich problems where ImplicitEulerHeun is insufficient
  • Fully stiff Stratonovich SDEs

Algorithm Description

Stratonovich analog of ISSEM with fully implicit treatment of both drift and diffusion terms using split-step approach.

Fully Implicit Features

  • Handles stiffness in both drift and diffusion for Stratonovich SDEs
  • Most expensive but most robust for Stratonovich stiff problems
  • Requires solving nonlinear systems for both drift and diffusion

References

  • Split-step implicit methods for fully stiff Stratonovich SDEs

Method Selection Guide

Problem Classification:

  1. Mildly stiff drift: ImplicitEM, ImplicitRKMil
  2. Stiff additive noise: SKenCarp (highly recommended)
  3. Fully stiff (including diffusion): ISSEM, ISSEulerHeun
  4. Stratonovich problems: ImplicitEulerHeun, ISSEulerHeun

Performance Ranking for Stiff Problems:

  1. SKenCarp - Best for stiff problems with additive noise
  2. ISSEM/ISSEulerHeun - For fully stiff problems
  3. ImplicitRKMil - Higher order for mildly stiff problems
  4. ImplicitEM - Robust fallback option

Understanding Stiffness in SDEs

Drift Stiffness: Large negative eigenvalues in the drift term f(u,t) requiring small time steps for explicit stability.

Diffusion Stiffness: Large coefficients in the diffusion term g(u,t) causing stability issues.

Detection: If explicit methods require very small dt or produce unstable solutions, try implicit methods.

Configuration Options

All implicit methods share common configuration parameters:

# Linear solver options
ImplicitEM(linsolve = KrylovJL_GMRES())

# Nonlinear solver options  
ImplicitEM(nlsolve = NLNewton(max_iter=20))

# Jacobian computation
ImplicitEM(autodiff = true, chunk_size = 4)

# Theta method parameter
ImplicitEM(theta = 0.5)  # Trapezoidal rule

Symplectic Methods

For Hamiltonian SDEs, use symplectic variants:

SImplicitMidpoint()  # Symplectic implicit midpoint
STrapezoid()        # Symplectic trapezoidal rule

Performance Tips

  1. Jacobian: Provide analytical Jacobian when possible
  2. Linear solver: Choose appropriate solver for problem structure
  3. Preconditioning: Use preconditioners for large systems
  4. Theta parameter: θ=0.5 often provides good accuracy/stability balance

References

  • Standard implicit ODE methods adapted to SDEs
  • Milstein, G.N., "Numerical Integration of Stochastic Differential Equations"