[Paper Review] ODE Test Problems: a MATLAB suite of initial value problems

Steven Roberts, Andrey A. Popov+1 more|arXiv (Cornell University)|Jan 14, 2019

Meteorological Phenomena and Simulations18 references17 citations

TL;DR

This paper introduces ODE Test Problems (OTP), a MATLAB-based, object-oriented framework offering a comprehensive suite of initial value problems for testing numerical methods in scientific computing. It supports diverse problem types—linear, nonlinear, stiff, chaotic, and PDE-discretized systems—via an intuitive API with presets, exact solutions, and built-in visualization, enabling seamless integration with ODE solvers and data assimilation tools.

ABSTRACT

ODE Test Problems (OTP) is an object-oriented MATLAB package offering a broad range of initial value problems which can be used to test numerical methods such as time integration methods and data assimilation (DA) methods. It includes problems that are linear and nonlinear, homogeneous and nonhomogeneous, autonomous and nonautonomous, scalar and high-dimensional, stiff and nonstiff, and chaotic and nonchaotic. Many are real-world problems from fields such as chemistry, astrophysics, meteorology, and electrical engineering. OTP also supports partitioned ODEs for testing IMEX methods, multirate methods, and other multimethods. Functions for plotting solutions and creating movies are available for all problems, and exact solutions are provided when available. OTP is desgined for ease of use-meaning that working with and modifying problems is simple and intuitive.

Motivation & Objective

• To provide a standardized, easy-to-use framework for testing numerical methods on a wide range of initial value problems in scientific computing.
• To unify disparate, ad-hoc ODE test problems into a single, extensible MATLAB package with minimal boilerplate code.
• To support advanced testing needs such as IMEX, multirate, and data assimilation methods through partitioned and structured problem designs.
• To enable efficient benchmarking of ODE solvers by exposing Jacobians, adjoint vector products, and linear solvers.
• To facilitate reproducibility and experimentation by providing exact solutions, default presets, and visualization capabilities.

Proposed method

• The framework is implemented entirely in MATLAB with no external dependencies or MEX files, ensuring portability and ease of use.
• Each problem is encapsulated as an object with access to its right-hand side (RHS) function, time span, initial condition, and additional properties like Jacobians and adjoint operations.
• Presets are used to standardize problem configurations (e.g., 'Canonical' or 'GC' for specific initial conditions and parameters), accessible via tab-completion and named-value argument support.
• The package supports both explicit and implicit solvers by exposing the RHS function and its Jacobian, enabling integration with MATLAB’s ODE suite like ode45 and ode15s.
• For PDE-based problems like QGSO and Gray–Scott, spatial discretization is handled via finite differences, with operators like the Laplacian and Jacobian approximated using standard numerical methods.
• Built-in support for linear solvers (Cholesky, multigrid, GMRES) and distance functions enables efficient implementation of data assimilation and adjoint-based methods.

Experimental results

Research questions

• RQ1 How can a unified, extensible, and user-friendly framework be designed to support diverse ODE test problems for numerical method evaluation?
• RQ2 To what extent can a MATLAB-based ODE test suite support advanced numerical workflows, including IMEX, multirate, and data assimilation methods?
• RQ3 How effectively can exact solutions and adjoint operations be integrated into a high-level, object-oriented MATLAB interface for scientific benchmarking?
• RQ4 What is the performance and usability impact of providing built-in visualization and movie generation for complex dynamical systems?
• RQ5 Can a single framework efficiently handle both small-scale chaotic systems (e.g., Lorenz-96) and large-scale PDE-discretized models (e.g., QGSO) with minimal configuration?

Key findings

• The OTP framework enables users to define and solve complex ODEs with minimal code, such as running a 40-variable Lorenz-96 system with a single line: `model = csl.otp.lorenz96.presets.Canonical;`
• The framework successfully supports both explicit and implicit solvers, with Jacobian-based integration (e.g., via `ode15s`) achieving stable and accurate results on stiff problems like the 40-variable Lorenz-96 system.
• The QGSO model, with 16,129 state variables, was successfully implemented and tested using multiple linear solvers (Cholesky, multigrid, GMRES), demonstrating scalability and configurability for large-scale problems.
• The Gray–Scott model with 5,000 variables and stiff dynamics was effectively tested using the 'Random' preset, confirming the framework's utility for stiff integrator validation.
• The framework enables direct comparison of solution accuracy by modifying solver parameters—e.g., changing the linear solver tolerance in QGSO resulted in a relative error of 0.822 at final time, demonstrating sensitivity and configurability.
• Built-in movie generation and visualization functions allow real-time analysis of complex dynamics, such as the chaotic behavior in Lorenz-63 and the vorticity patterns in QGSO.