Linda Petzold

The Finite State Projection (FSP) method approximates the Chemical Master Equation (CME) by restricting the dynamics to a finite subset of the (typically infinite) state space, enabling direct numerical solution with computable error bounds. Adaptive variants update this subset in time, but multiscale systems with widely separated reaction rates remain challenging, as low-probability bottleneck states can carry essential probability flux and the dynamics alternate between fast transients and slowly evolving stiff regimes. We propose a flux-based adaptive FSP method that uses probability flux to drive both state-space pruning and time-step selection. The pruning rule protects low-probability states with large outgoing flux, preserving connectivity in bottleneck systems, while the time-step rule adapts to the instantaneous total flux to handle rate constants spanning several orders of magnitude. Numerical experiments on stiff, oscillatory, and bottleneck reaction networks show that the method maintains accuracy while using substantially smaller state spaces.

Gradient-based optimization of neural differential equations and other parameterized dynamical systems fundamentally relies on the ability to differentiate numerical solutions with respect to model parameters. In many stiff systems, it has been observed that sensitivities to parameters controlling fast-decaying modes become vanishingly small during training, leading to optimization difficulties. In this paper, we show that this vanishing gradient phenomenon is not an artifact of any particular method but a universal feature of all A-stable and L-stable stiff numerical integration schemes. We analyze the rational stability function for general stiff integration schemes and demonstrate that the relevant parameter sensitivities, governed by the derivative of the stability function, decay to zero for large stiffness. Explicit formulas for common stiff integration schemes are provided, which illustrate the mechanism in detail. Finally, we rigorously prove that the slowest possible rate of decay for the derivative of the stability function is O(|z|−1), revealing a fundamental limitation: all A-stable and L-stable time-stepping methods inevitably suppress parameter gradients in stiff regimes, posing a significant barrier for training and parameter identification in stiff neural ordinary differential equations.

We introduce a new approach for directly training stiff neural ordinary differential equations using explicit rational Taylor series methods. These explicit schemes, of second and third order, are A-stable and require only a single linear solve or matrix inversion per time step, providing both high efficiency and strong numerical stability. In contrast to traditional implicit methods, which are often computationally intensive due to the need to solve nonlinear systems at each step, the explicit rational Taylor methods deliver accurate and stable learning of stiff dynamics at a significantly reduced computational cost. Our results show that these methods effectively learn stiff systems, such as the van der Pol oscillator, at large step sizes without encountering stability issues that challenge implicit schemes. This work broadens the capabilities of data-driven simulation, supporting efficient and accurate modeling of complex dynamics in stiff neural differential equations, mesh-based simulation, and physics-informed neural networks.

Stiff ordinary differential equations (ODEs) are common in many science and engineering fields, but standard neural ODE approaches struggle to accurately learn these stiff systems, posing a significant barrier to the widespread adoption of neural ODEs. In our earlier work, we addressed this challenge by utilizing single-step implicit methods for solving stiff neural ODEs. While effective, these implicit methods are computationally costly and can be complex to implement. This paper expands on our earlier work by exploring explicit exponential integration methods as a more efficient alternative. We evaluate the potential of these explicit methods to handle stiff dynamics in neural ODEs, aiming to enhance their applicability to a broader range of scientific and engineering problems. We found the integrating factor Euler (IF Euler) method to excel in stability and efficiency. While implicit schemes failed to train the stiff van der Pol oscillator, the IF Euler method succeeded, even with large step sizes. However, IF Euler’s first-order accuracy limits its use, leaving the development of higher-order methods for stiff neural ODEs an open research problem.

Xianjun Yang, Wei Cheng, Xujiang Zhao, Wenchao Yu, Linda Ruth Petzold, Haifeng Chen. Findings of the Association for Computational Linguistics: NAACL 2025. 2025.

Stiff ordinary differential equations (ODEs) are common in many science and engineering fields, but standard neural ODE approaches struggle to accurately learn these stiff systems, posing a significant barrier to widespread adoption of neural ODEs. In our earlier work, we addressed this challenge by utilizing single-step implicit methods for solving stiff neural ODEs. While effective, these implicit methods are computationally costly and can be complex to implement. This paper expands on our earlier work by exploring explicit exponential integration methods as a more efficient alternative. We evaluate the potential of these explicit methods to handle stiff dynamics in neural ODEs, aiming to enhance their applicability to a broader range of scientific and engineering problems. We found the integrating factor Euler (IF Euler) method to excel in stability and efficiency. While implicit schemes failed to train the stiff Van der Pol oscillator, the IF Euler method succeeded, even with large step sizes. However, IF Euler's first-order accuracy limits its use, leaving the development of higher-order methods for stiff neural ODEs an open research problem.

Stiff systems of ordinary differential equations (ODEs) are pervasive in many science and engineering fields, yet standard neural ODE approaches struggle to learn them. This limitation is the main barrier to the widespread adoption of neural ODEs. In this paper, we propose an approach based on single-step implicit schemes to enable neural ODEs to handle stiffness and demonstrate that our implicit neural ODE method can learn stiff dynamics. This work addresses a key limitation in current neural ODE methods, paving the way for their use in a wider range of scientific problems.

Stiff systems of ordinary differential equations (ODEs) arise in a wide range of scientific and engineering disciplines and are traditionally solved using implicit integration methods due to their stability and efficiency. However, these methods are computationally expensive, particularly for applications requiring repeated integration, such as parameter estimation, Bayesian inference, neural ODEs, physics-informed neural networks, and MeshGraphNets. Explicit exponential integration methods have been proposed as a potential alternative, leveraging the matrix exponential to address stiffness without requiring nonlinear solvers. This study evaluates several state-of-the-art explicit single-step exponential schemes against classical implicit methods on benchmark stiff ODE problems, analyzing their accuracy, stability, and scalability with step size. Despite their initial appeal, our results reveal that explicit exponential methods significantly lag behind implicit schemes in accuracy and scalability for stiff ODEs. The backward Euler method consistently outperformed higher-order exponential methods in accuracy at small step sizes, with none surpassing the accuracy of the first-order integrating factor Euler method. Exponential methods fail to improve upon first-order accuracy, revealing the integrating factor Euler method as the only reliable choice for repeated, inexpensive integration in applications such as neural ODEs and parameter estimation. This study exposes the limitations of explicit exponential methods and calls for the development of improved algorithms.

Stiff systems of ordinary differential equations (ODEs) are pervasive in many science and engineering fields, yet standard neural ODE approaches struggle to learn them. This limitation is the main barrier to the widespread adoption of neural ODEs. In this paper, we propose an approach based on single-step implicit schemes to enable neural ODEs to handle stiffness and demonstrate that our implicit neural ODE method can learn stiff dynamics. This work addresses a key limitation in current neural ODE methods, paving the way for their use in a wider range of scientific problems.