on From ODE45 to others

Overview

This post discusses the practical usage of MATLAB’s ODE solvers, specifically ode45, ode23, and ode113. It covers performance tips regarding output control, mass matrix handling, and solver selection based on stiffness and accuracy requirements.

🏷️ MATLAB ODE Solvers

ode45 is a popular function within MATLAB. It uses an explicit Runge-Kutta 4(5) pair as an integrator for first-order equations:

$$\dot{\xi }=f(\xi ,t)$$

The standard function definition is:

[t,x] = ode45(@fname, tspan, xinit, options)

🏷️ Optimization of tspan

The tspan parameter controls the output size. The default setting with only two elements enforces integration over the entire domain and outputs all intermediate time steps. To prevent rapid memory consumption in large systems, it is more efficient to specify only the time steps of interest in tspan.

🏷️ Handling the Mass Matrix

While the standard case is a simple first-order ODE, some problems involve a mass matrix:

$$M(t,\xi )\dot{\xi }=f(\xi ,t)$$

To solve this, the Mass property must be configured via odeset. For invertible and fixed mass matrices $M$, it is significantly more efficient to solve the system as:

$$\dot{\xi }=M^{-1}f(\xi )$$

🏷️ Solver Configuration

🏷️ Error Control

The default RelTol and AbsTol are $10^{-3}$ and $10^{-6}$ respectively. Setting NormControl to 'on' can improve performance by switching from element-wise error checking to a single integral error norm.

🏷️ Output Properties

• NonNegative: Specifies indices where the solution must remain non-negative.
• OutputFcn: A callback handle for monitoring the integration progress.
• Refine: Controls the density of the output points; often redundant when specific timestamps are provided in tspan.

🏷️ Performance Benchmarks

In wave equation simulations, we encountered systems of the form:

$$M[\dot{\mathbf{u}},\dot{\mathbf{v}},\dot{\mathbf{w}},\dot{\mathbf{z}}]=\mathbf{f}(\mathbf{u},\mathbf{v},\mathbf{w},\mathbf{z})$$

Using a full block-diagonal matrix for $M$ is inefficient. In MATLAB, calculating the components via a single mldivide is significantly faster:

$$[\dot{\mathbf{u}},\dot{\mathbf{v}},\dot{\mathbf{w}},\dot{\mathbf{z}}]=M^{-1}\mathbf{f}(\mathbf{u},\mathbf{v},\mathbf{w},\mathbf{z})$$

Benchmarks show that ode113 can be more economic than ode45 when stringent local accuracy is not the primary goal.

Timing comparison for PML (Minimum outputs)

• ode45: 234 seconds
• ode23: 524 seconds
• ode113: 116 seconds

🏷️ Solver Reference

The following table summarizes the primary non-stiff solvers in MATLAB:

Solver	Method	Type	Best Usage
ode45	Runge-Kutta (4,5)	One-step	First try for most problems.
ode23	Runge-Kutta (2,3)	One-step	Crude tolerances or moderate stiffness.
ode113	Adams-Bashforth-Moulton	Multistep	Stringent tolerances; expensive evaluations.

For stiff systems, consider using ode15s, ode23s, ode23t, or ode23tb.

📚 References