ODE Solvers

To solve an ODE directly without creating a solve block, use one of the ODE solvers, which solve systems of ODEs of the form

Form of ODE system for solvers

where y is vector of unknown functions of the independent variable x. To solve a higher-order ODE, rewrite it as a system of first order ODEs.

The ODE solvers are divided into two types: solvers for stiff systems and solvers for non-stiff systems. A system of ODEs written in matrix form as y = Ax, is called stiff if the matrix A is nearly singular. Otherwise, the system is non-stiff.

Solvers for Non-Stiff Systems

Adams(init, x1, x2, npoints, D, [tol]) Uses Adams methods.

rkfixed(init, x1, x2, npoints, D) Uses the fourth-order Runge-Kutta fixed-step method.

Rkadapt(init, x1, x2, npoints, D) Uses the fourth-order Runge-Kutta with adaptive step-size.

Bulstoer(init, x1, x2, npoints, D) Uses the Bulirsch-Stoer method, which is slightly more accurate than Runge-Kutta, but requires a smoothly varying system.

Solvers for Stiff Systems

BDF(init, x1, x2, npoints, D, [J], [tol]) Uses backward differentiation formula methods

Radau(init, x1, x2, npoints, D, [J], [M], [tol]) Uses the Radau5 method.

Stiffb(init, x1, x2, npoints, D, AJ) Uses the Bulirsch-Stoer method.

Stiffr(init, x1, x2, npoints, D, AJ) Uses Rosenbrock method.

Hybrid Solver

AdamsBDF(init, x1, x2, npoints, D, [J], [tol]) Determines whether a system is stiff or non-stiff and calls Adams or BDF accordingly.

Note: AdamsBDF is the default solver for Odesolve.

Each of these functions returns an (npoints + 1) × (n + 1) solution matrix, where n is the number of unknowns. The first column of the matrix contains the values of x at which the solutions are evaluated. These values are npoints + 1 evenly spaced numbers between x1 and x2. The remaining columns contain the values of the solutions y0, y1, ..., yn - 1 corresponding to the values of x in the first column.

Arguments:

Note: For the solvers rkfixed, RkAdapt, Bulstoer, Stiffb, and Stiffr, you can specify a scalar tolerance for the solution by defining the variable TOL before calling the solver. For example,

TOL:= 10-6

sets the tolerance for all variables as 10-6.

Note that TOL does not affect the tolerances for the solvers Adams, BDF, or Radau. For these solvers, you must use the optional argument tol, described above, to specify tolerances.

Rewriting a Higher-Order ODE as a System of First-Order ODEs

To rewrite the nth-order ODE

A general nth order ODE

as a system of first-order ODEs, define variables y0, y1, ... yn as follows:

Define new variables for nth order ODE

In terms of the new variables, the original equation becomes

nth order ODE in terms of new variables

This is equivalent to the following system of first-order ODEs:

System of first-order ODEs   where   System of first-order ODEs

If you solve the original equation for yn and substitute the result in the vector, you get a vector function that you can apply the ODE solvers to:

Final form for the system of first-order ODEs

Example:

For example, you can rewrite the second-order ODE

single- function ODE

using vector subscripts, as

Rewritten ODE with implied left-hand-side implied left-hand side

Computing the Jacobian Matrix

You can compute the Jacobian matrix, used by the solvers BDF and Radau using the function Jacob. For the preceding example, the Jacobian matrix J is

Jacobian matrix of partial derivatives of the system of equations in D

To compute the augmented Jacobian, used by the solvers Stiffb and Stiffr, start by computing the first column of the matrix, using the function Jacob, as follows:

First column of augmented Jacobian matrix.

Then concatenate this column with the with the matrix J, using the augment function, as follows:

Augmented Jacobian matrix

Notes:

QuickSheet - First-Order ODE Initial Value Problem
QuickSheet - First-Order System of ODEs

Related Topics