next up previous contents
Next: Numerical Integration Up: Calculus Previous: Calculus

Differential Equations

For exact (symbolic) differentiations and solutions, see gif Symbolic Math.

The function ode23 is used to solve ordinary differential equations with 2nd and 3rd order Runge-Kutta formulas. The variant ode23p also plots the results. The function ode45 uses 4th and 5th order Runge-Kutta formulas and so is slower but more accurate; it otherwise acts as ode23. A simple demo is available with odedemo.

The basic format is [t, y] = ode23( tex2html_wrap_inline2867 ). Here t is the single independent variable and y is a vector (or scalar) of the dependent variables. The scalars tex2html_wrap_inline2873 and tex2html_wrap_inline2875 specify the initial and final times; tex2html_wrap_inline2877 is a vector with the initial state of the dependent variables. The system of equations to be solved is in the m-file mysystem.m.

(An optional fifth argument can specify the tolerance, which defaults to tex2html_wrap_inline2881 ; an optional boolean sixth argument specifies whether to give status displays while integrating, defaulting to 0 (False).)

The results are a vector t of times and a matrix y in which the nth column gives the values of the nth dependent variable ( tex2html_wrap_inline2891 ) at those positions. Thus plot(t, y) suffices for basic display. (See ode23p for automatic phase plane plots.)

The m-file defines a function with the same name as the system of equations, with one output (a vector of derivatives) and two inputs (a scalar time, and a vector holding the variables whose derivatives are returned). It is important to realize that y, tex2html_wrap_inline2895 , and tex2html_wrap_inline2897 are independent variables as far as MATLAB is concerned until you define their relationships. This means that the y vector will typically contain some intermediate derivatives of your basic dependent variables. The output vector is defined to be the derivative vector of the input vector; the process by which you specify this defines the system.

As an example, consider the system of equations

eqnarray968

Here x and y are our basic dependent variables, and t is our independent variable (since we used time derivatives). Since the system is second-order in x, we must consider tex2html_wrap_inline2909 as a necessary part of the state of the system and include it in our vectors of dependent variables. Thus we want a function that takes a scalar t and a vector holding (in whatever order) tex2html_wrap_inline2913 , and returns a vector holding in the same order tex2html_wrap_inline2915 , i.e. the derivative of the input vector.

function der = mysystem(t, state)
% MYSYSTEM  to solve the system x'' + y' - 2x' = 3, y' = 3x + 5.
%       mysystem(t, [x(t) y(t) x'(t)])
%       returns     [x'(t) y'(t) x''(t)]
%       Suitable for use with ode23, ode45.

% First assign the input values to variables with names that
% make what they are easier to remember.
x = state(1);
y = state(2);
xdot = state(3);

% Now calculate the output values.
% The first, the derivative of x, is trivial since we require
% it as an input value.  We'll put a line to do it anyway as a reminder.
xdot = xdot;

% We can get the second derivative of x in terms of x' and y' from
% our first equation, but we need y' first.  We'll use the second for that.
ydot = 3*x + 5;

% Now we can get the second derivative of x.
xdoubledot = 3 - 2*xdot - ydot;

% Construct the vector of derivatives in the same order as
% the inputs.  We took them as x, y, x', so we must return x', y', x''.
der = [xdot, ydot, xdoubledot];

With this function in the file mysystem.m, we can solve the system over any stretch of time with any initial values:

.1ex>> t0 = 0;
.1ex>> tf = 10;
.1ex>> x0 = 0;
.1ex>> y0 = 1;
.1ex>> xdot0 = -1;
.1ex>> state0 = [x0 y0 xdot0];
.1ex>> [t state] = ode23('mysystem', t0, tf, state0);
.1ex>> x = state(:, 1); y = state(:, 2); xdot = state(:, 3);
.1ex>> figure(1); plot(t, x); title('x(t)')
.1ex>> figure(2); plot(t, y); title('y(t)')
.1ex>> figure(3); plot(t, xdot); title('xdot(t)')
.1ex>> figure(4); plot(x, y); title('phase plot')

While we have considered the independent variable to be time, it is obvious no such interpretation is required.

MATLAB does not have built-in functionality to solve systems with multiple independent variables.


next up previous contents
Next: Numerical Integration Up: Calculus Previous: Calculus

sepherke
Sat Mar 21 21:42:28 EST 1998