Cadabra manual: Ndsolve

ndsolve

Numerically solve ordinary differential equations

The ndsolve function determines numerical solutions to ordinary differential equations. These should be provided in the form of a system of first order equations. The result is one or more interpolating functions, which can be plotted or evaluated.

In order to use ndsolve, one has to set up the dependent and independent variables, and also determine how one wants to write the derivatives. So the starting point is something similar to:

x::Coordinate; \dot{#}::Derivative(x); y::Depends(x);

$\displaystyle{}\text{Property Coordinate attached to }x.$

$\displaystyle{}\text{Property Derivative attached to }\backslash\texttt{dot}\{\#\}.$

$\displaystyle{}\text{Property Depends attached to }y.$

The differential equation we want to solve is entered as an equation with the derivative on the left-hand side.

eq1 := \dot{y} = y \cos(x + y);

$\displaystyle{}\dot{y} = y \cos\left(x +y\right)$

This can now be solved with ndsolve, by passing it the differential equation, the initial conditions (one in our case, for $y$) and the range of the independent variable $x$,

sol1 = ndsolve(eq1, {$y$: 1}, ($x$, 0, 30));

$\displaystyle{}\left[y = \square{} (x)\right]$

The result above is an equation for $y$ in terms of an anonymous function (the 'box'), which can be used directly in the plot function (anonymous functions are similar to Mathematica's InterpolatingFunction).

from cdb.graphics.plot import * plot(sol1, {$x$: (0,30)} );

The logic is similar for systems of ODEs. Declare the dependent coordinate, the functions to solve for, and the differential equation:

t::Coordinate. {x, y}::Depends(1t). \dot{#}::Derivative(t). x:=x. y:=y. t:=t.

eq3 := [ \dot{x} = -y - x**2, \dot{y} = 2 x - y**3 ];

$\displaystyle{}\left[\dot{x} = -\,y -\,{x}^{2\,}, \dot{y} = 2\,x -\,{y}^{3\,}\right]$

Then solve using ndsolve. This now returns a set of equations, one for each function that it has solved for, in terms of two anonymous functions.

res3 = ndsolve(eq3, {x: 1, y: 1}, (t, 0, 20));

$\displaystyle{}\left[x = \square{} (t), y = \square{} (t)\right]$

plot(res3, {$t$: (0,20)});