home *** CD-ROM | disk | FTP | other *** search
/ QBasic & Borland Pascal & C / Delphi5.iso / C / Samples / C-SSP.ARJ / DESORK.C < prev    next >
Encoding:
Text File  |  1984-08-07  |  1.1 KB  |  39 lines
  1.    desork(iter,h,x,y,dy,iflag)
  2.  
  3.       /* subroutine to solve second-order ordinary differential */
  4.       /* equations y" = f(x,y,y') by Runge-Kutta method.        */
  5.  
  6.       int iflag,iter;
  7.       float dy[],h,x[],y[];
  8.  
  9.     {
  10.       int i;
  11.       float dyn,tk0,tm0,tk1,tm1,tk2,tm2,tk3,tm3,xn,yn;
  12.       double ydopri();
  13.  
  14.       for(i = 0; i <= iter-2; i++)
  15.       {
  16.        tk0 = h * dy[i];
  17.        xn = x[i];
  18.        yn = y[i];
  19.        dyn = dy[i];
  20.        tm0 = h * ydopri(xn,yn,dyn,iflag);
  21.        xn = x[i] + (h/2.);
  22.        yn = y[i] + (tk0/2.);
  23.        dyn = dy[i] + (tm0/2.);
  24.        tk1 = h * (dy[i] + (tm0/2.));
  25.        tm1 = h * ydopri(xn,yn,dyn,iflag);
  26.        yn = y[i] + (tk1/2.);
  27.        dyn = dy[i] + (tm1/2.);
  28.        tk2 = h * (dy[i] + (tm1/2.));
  29.        tm2 = h * ydopri(xn,yn,dyn,iflag);
  30.        xn = x[i] + h;
  31.        yn = y[i] + tk2;
  32.        dyn = dy[i] + tm2;
  33.        tk3 = h * (dy[i] + tm2);
  34.        tm3 = h * ydopri(xn,yn,dyn,iflag);
  35.        y[i+1] = y[i] + (tk0 + 2.*tk1+ 2.*tk2 + tk3)/6.;
  36.        dy[i+1] = dy[i] + (tm0+ 2.*tm1 + 2.*tm2 + tm3)/6.;
  37.       }
  38.     }
  39.