Next: , Previous: Basic Math, Up: ode


8.2 Simple examples using ode

The following examples should illustrate the procedure of stating an initial value problem and solving it with ode. If these examples are too elementary, see Input Language, for a formal specification of the ode input language. There is also a directory containing examples of ode input, which is distributed along with the GNU plotting utilities. On most systems it is installed as /usr/share/ode or /usr/local/share/ode.

Our first example is a simple one, namely

     y'(t) = y(t)

with the initial condition

     y(0) = 1

The solution to this differential equation is

     y(t) = e^t.

In particular

     y(1) = e^1 = 2.718282

to seven digits of accuracy.

You may obtain this result with the aid of ode by typing on the command line the sequence of commands

     ode
     y' = y
     y = 1
     print t, y
     step 0, 1

Two columns of numbers will appear. Each line will show the value of the independent variable t, and the value of the variable y, as t is `stepped' from 0 to 1. The last line will be

     1 2.718282

as expected. You may use the ‘-p’ option to change the precision. If, for example, you type ‘ode -p 10’ rather than ‘ode’, you will get ten digits of accuracy in the output, rather than seven (the default).

After the above output, ode will wait for further instructions. Entering for example the line

     step 1, 0

should yield two more columns of numbers, containing the values of t and y that are computed when t is stepped back from 1 to 0. You could type instead

     step 1, 2

to increase rather than decrease t. To exit ode, you would type a line containing only ‘.’, i.e. a single period, and tap `return'. ode will also exit if it sees an end-of-file indicator in its input stream, which you may send from your terminal by typing control-D.

Each line of the preceding example should be self-explanatory. A ‘step statement sets the beginning and the end of an interval over which the independent variable (here, t) will range, and causes ode to set the numerical scheme in motion. The initial value appearing in the first ‘step’ statement (i.e., 0) and the assignment statement

     y = 1

are equivalent to the initial condition y(0) = 1. The statements ‘y' = y and ‘y = 1’ are very different: ‘y' = y’ defines a way of computing the derivative of y, while ‘y = 1’ sets the initial value of y. Whenever a ‘step’ statement is encountered, ode tries to step the independent variable through the interval it specifies. Which values are to be printed at each step is specified by the most recent ‘print’ statement. For example,

     print t, y, y'

would cause the current value of the independent variable t, the variable y, and its derivative to be printed at each step.

To illustrate ode's ability to take its input or the initial part of its input from a file, you could prepare a file containing the following lines:

     # an ode to Euler
     y  = 1
     y' = y
     print t, y, y'

Call this file euler. (The ‘#’ line is a comment line, which may appear at any point. Everything from the ‘# to the end of the line on which it appears will be ignored.) To process this file with ode, you could type on your terminal

     ode -f euler
     step 0, 1

These two lines cause ode to read the file euler, and the stepping to take place. You will now get three quantities (t, y, and y') printed at each of the values of t between 0 and 1. At the conclusion of the stepping, ode will wait for any further commands to be input from the terminal. This example illustrates that

     ode -f euler

is not equivalent to

     ode < euler

The latter would cause ode to take all its input from the file euler, while the former allows subsequent input from the terminal. For the latter to produce output, you would need to include a ‘step’ line at the end of the file. You would not need to include a ‘.’ line, however. ‘.’ is used to terminate input only when input is being read from a terminal.

A second simple example involves the numerical solution of a second-order differential equation. Consider the initial value problem

     y''(t) = -y(t)
     y(0) = 0
     y'(0) = 1

Its solution would be

     y(t) = sin(t)

To solve this problem using ode, you must express this second-order equation as two first-order equations. Toward this end you would introduce a new function, called yp say, of the independent variable t. The pair of equations

     y' = yp
     yp' = -y

would be equivalent to the single equation above. This sort of reduction of an n'th order problem to n first order problems is a standard technique.

To plot the variable y as a function of the variable t, you could create a file containing the lines

     # sine : y''(t) = -y(t), y(0) = 0, y'(0) = 1
     sine' = cosine
     cosine' = -sine
     sine = 0
     cosine = 1
     print t, sine

(y and yp have been renamed sine and cosine, since that is what they will be.) Call this file sine. To display the generated data points on an X Window System display as they are generated, you would type

     ode -f sine | graph -T X -x 0 10 -y -1 1
     step 0, 2*PI
     .

After you type the ode line, graph -T X will pop up a window, and after you type the ‘step’ line, the generated dataset will be drawn in it. The ‘-x 0 10’ and ‘-y -1 1’ options, which set the bounds for the two axes, are necessary if you wish to display points in real time: as they are generated. If the axis bounds were not specified on the command line, graph -T X would wait until all points are read from the input before determining the bounds, and drawing the plot.

A slight modification of this example, showing how ode can generate several datasets in succession and plot them on the same graph, would be the following. Suppose that you type on your terminal the following lines.

     ode -f sine | graph -T X -C -x 0 10 -y -1 1
     step 0, PI
     step PI, 2*PI
     step 2*PI, 3*PI
     .

Then the sine curve will be traced out in three stages. Since the output from each ‘step’ statement ends with a blank line, graph -T X will treat each section of the sine curve as a different dataset. If you are using a color display, each of the three sections will be plotted in a different color. This is a feature provided by graph, which normally changes its linemode after each dataset it reads. If you do not like this feature, you may turn it off by using graph -T X -B instead of graph -T X.

In the above examples, you could use any of the other variants of graph instead of graph -T X. For example, you could use graph -T ps to obtain a plot in encapsulated Postscript format, by typing

     ode -f sine | graph -T ps > plot.ps
     step 0, 2*PI
     .

You should note that of the variants of graph, the variants graph -T png, graph -T pnm, graph -T gif, graph -T svg, graph -T ai, graph -T ps, graph -T cgm, graph -T fig, graph -T pcl and graph -T hpgl do not produce output in real time, even when the axis bounds are specified with the ‘-x’ and ‘-y options. So if any of these variants is used, the plot will be produced only when input from ode is terminated, which will occur when you type ‘..

In the preceding examples, the derivatives of the dependent variables were specified by comparatively simple expressions. They are allowed to be arbitrarily complicated functions of the dependent variables and the independent variable. They may also involve any of the functions that are built into ode. ode has a fair number of functions built in, including abs, sqrt, exp, log, log10, sin, cos, tan, asin, acos, atan, sinh, cosh, tanh, asinh, acosh, and atanh. Less familiar functions which are built into it are besj0, besj1, besy0, besy1, erf, erfc, inverf, lgamma, gamma, norm, invnorm, ibeta, and igamma. These have the same definitions as in the plotting program gnuplot. (All functions take a single argument, except for ibeta, which takes three, and igamma, which takes two). ode also knows the meaning of the constant ‘PI’, as the above examples show. The names of the preceding functions are reserved, so, e.g., ‘cos’ and ‘sin’ may not be used as names for variables.

Other than the restriction of avoiding reserved names and keywords, the names of variables may be chosen arbitrarily. Any sequence of alphanumeric characters starting with an alphabetic character may be used; the first 32 characters are significant. It is worth noting that ode identifies the independent variable by the fact that it is (or should be) the only variable that has not appeared on the left side of a differential equation or an initial value assignment. If there is more than than one such variable then no stepping takes place; instead, an error message is printed. If there is no such variable, a dummy independent variable is invented and given the name ‘(indep)’, internally.