Home > src > rk4.m

rk4

PURPOSE ^

RK4 4th order Runge-Kutta integration

SYNOPSIS ^

function [x,Y] = rk4(f,dt,x,P1,P2,P3,Y)

DESCRIPTION ^

RK4  4th order Runge-Kutta integration

 Syntax:
   [x,Y] = rk4(f,dt,x,[P1,P2,P3,Y])

 In:
    f - Name of function in form f(x,P(:)) or
        inline function taking the same parameters.
        In chained case the function should be f(x,y,P(:)).
   dt - Delta time as scalar.
    x - Value of x from the previous time step.
   P1 - Values of parameters of the function at initial time t
        as a cell array (or single plain value). Defaults to empty
        array (no parameters).
   P2 - Values of parameters of the function at time t+dt/2 as
        a cell array (or single plain value). Defaults to P1 and
        each empty (or missing) value in the cell array is replaced
        with the corresponding value in P1.
   P3 - Values of parameters of the function at time t+dt.
        Defaults to P2 similarly to above.
    Y - Cell array of partial results y1,y2,y3,y4 in the RK algorithm
        of the second parameter in the interated function. This can be
        used for chaining the integrators. Defaults to empty.
 
 Out:
    x - Next value of X
    Y - Cell array of partial results in Runge-Kutta algorithm.

 Description:
   Perform one fourth order Runge-Kutta iteration step
   for differential equation

       dx/dt = f(x(t),P{:})

   or in the chained case

       dx/dt = f(x(t),y(t),P{:})
       dy/dt = g(y(t),P{:})

   - Example 1. Simple integration of model

       dx/dt = tanh(x), x(0) = 1

     can be done as follows:

       X = [];
       x = 1;
       f = inline('tanh(x)','x');
       for i=1:100
         x = rk4(f,0.1,x);
         X = [X x];
       end

   - Example 2. Chaining of integrators. Consider a
     model of the form

       dx/dt = x+y,     x(0)=1
       dy/dt = tanh(y), y(0)=2

     The equations can be now integrated as follows:

       XY = [];
       x = 1;
       y = 2;
       fx = inline('x+y','x','y');
       fy = inline('tanh(y)','y');
       for i=1:100
         [y,YY] = rk4(fy,0.1,y);
         x = rk4(fx,0.1,x,{},{},{},YY);
         XY = [XY [x;y]];
       end

   which produces exactly the same result as

      XY = [];
      xy = [1;2];
      fxy = inline('[xy(1)+xy(2);tanh(xy(2))]','xy');
      for i=1:100
         xy = rk4(fxy,0.1,xy);
         XY = [XY xy];
      end

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:
Generated on Fri 12-Aug-2011 15:08:47 by m2html © 2005