/***************************************************************************/ /* The following is a C subroutine implementing the fourth-order */ /* Runge-Kutta method for integrating the following differential equation: */ /* */ /* dx/dt = f(x,y,input) */ /* dy/dt = g(x,y) */ /* */ /* In the routine, 'x1' and 'y1' are the input values of the current */ /* time step, and 'x2' and 'y2' are the output values of the next */ /* time step. The step length is given by 'h'. */ /***************************************************************************/ runge_kutta(x1,y1,x2,y2,input) double x1, y1, *x2, *y2, input; { double kx1, kx2, kx3, kx4, ky1, ky2, ky3, ky4; kx1 = h*f(x1,y1,input); ky1 = h*g(x1,y1); kx2 = h*f(x1+kx1/2.0, y1+ky1/2.0, input); ky2 = h*g(x1+kx1/2.0, y1+ky1/2.0); kx3 = h*f(x1+kx2/2.0, y1+ky2/2.0, input); ky3 = h*g(x1+kx2/2.0, y1+ky2/2.0); kx4 = h*f(x1+kx3, y1+ky3, input); ky4 = h*g(x1+kx3, y1+ky3); *x2 = x1 + (kx1 + 2.0*kx2 + 2.0*kx3 + kx4)/6.0; *y2 = y1 + (ky1 + 2.0*ky2 + 2.0*ky3 + ky4)/6.0; }