Next: Solving the pendulum equation Up: tutorial6 Previous: Example: the nonlinear pendulum

MATLAB code for the second-order Runge-Kutta method (RK2) for two or more first-order equations

First we will solve the linearized pendulum equation (3) using RK2.

We can use a script that is very similar to rk2.m that we wrote last week to solve a single first-order ODE using the RK2 method. So now copy this file, or your version from last week, to rk2_many.m. The many' equations in this case are just two, but the code we write will work for any number of first-order ODEs.

In order to deal with two equations we need to replace every occurrence of the one-dimensional row vector  y(1:n) , where is the number of time-points, say, with a two-dimensional vector  y(1:2,1:n)  where the two rows correspond to the variables and . (We can in fact use a similar approach for any number of equations).

In rk2_many.m, at every place where  y(1), y(n), y(n+1) and so on appear, they should be replaced by y(:,1), y(:,n), y(:,n+1) respectively, e.g. the input statement becomes:

 y(:,1) = input(' enter y1(0) and y2(0) as a column vector: '); 

and the basic RK2 iteration:

 y(:,n+1) = y(:,n)+k2; 

How many such changes did you make in rk2_many.m?

You must alter all calls to rhs.m to rhs2.m in the script rk2_many.m.

Add a line of code to plot i.e.  y(1,:)  as a function of , i.e.  t (see plotting ).

The MATLAB function rhs2.m is very similar to rhs.m from last week. However, the variable  y  is now a column vector, and the function must return a column vector for the slopes  dydt. Thus for our linearized pendulum problem (3) we need:

    function [dydt]=rhs2(t,y)
%
% function returns the values of the rhs's
% for the system:
%    dy1/dt = y2
%    dy2/dt = -y1
%
dydt=[y(2); -y(1)];
`

Next: Solving the pendulum equation Up: tutorial6 Previous: Example: the nonlinear pendulum