Introduction to Ordinary Differential Equations#

Generalities#

ODE stands for ordinary differential equation.

The ODE module integrates the solution of first-order ordinary differential equations in GEMSEO by computing the solution of an Initial Value Problem (IVP).

An Initial Value Problem (IVP) is composed of the following entities:

  • an initial state \(y_0 \in \mathbb{R}^n\);

  • a time interval \([t_0, t_f]\);

  • a function \(f: [t_0, t_f] \times \mathbb{R}^n \rightarrow \mathbb{R}^n\).

The solution of an IVP consists in identifying the function \(y: [t_0, t_f]\rightarrow \mathbb{R}^n\) solving the following ODE:

\[\begin{split}\begin{cases} y(t_0) &= y_0 \\ \frac{\mathrm{d} y}{\mathrm{d} t}(t) &= f(t, y(t)) \qquad \mbox{for all }t \in [t_0, t_f]. \end{cases}\end{split}\]

Algorithms for the numerical solution#

Multiple algorithms are available in the literature. The algorithms available in GEMSEO, developed in the method solve_ivp of the library scipy.integrate, are:

  • the explicit Runge-Kutta algorithms (RK45, RK23, DOP853),

  • an implicit Runge-Kutta algorithm (Radau),

  • and two algorithms based on a backwards differentiation formula (BDF and LSODA).

The algorithms Radau, BDF, and LSODA require the knowledge of the Jacobian of the function \(f\) with respect to the state \(y\): \(J f = \frac{\partial f}{\partial y}\). The Jacobian can be either passed to the algorithm, or computed by finite differences.

Further algorithms for the solution of IVPs are available in the plugin gemseo-petsc. The plugin gemseo-petsc provides also an adjoint mode to perform a sensitivity analysis on the solution of the ODE with respect to its initial values and the design parameters.

../_images/ODEProblem_ODEResult_attributes_description.png

Terminating events#

For some problems, it might be interesting not to integrate the dynamic for the entire time interval \([t_0, t_f]\), but only up to the realization of a terminating condition. Such conditions are encoded by event functions: real-valued continuous functions \(g_1, \ldots, g_m: [t_0, t_f] \times \mathbb{R}^n \rightarrow \mathbb{R}\). The terminating condition is realized when any of the event function crosses the threshold \(0\).