home *** CD-ROM | disk | FTP | other *** search
/ QBasic & Borland Pascal & C / Delphi5.iso / C / Samples / C-SSP.ARJ / DESOMA.C < prev    next >
Encoding:
Text File  |  1984-08-07  |  1.1 KB  |  38 lines
  1.    desoma(iter,h,x,y,dy,iflag,d)
  2.  
  3.       /* subroutine to solve second-order ordinary differential  */
  4.       /* equations y" = f(x,y,y') by Modified adams' method.     */
  5.  
  6.       int iflag,iter;
  7.       float d[],dy[],h,x[],y[];
  8.  
  9.     {
  10.       int i,lda;
  11.       float dyi,u1,u2,u3,u4,yi;
  12.       double ydopri();
  13.  
  14.       d[0] = 0.;
  15.       d[1] = 0.;
  16.       d[2] = 0.;
  17.       lda = 4;
  18.       desork(lda,h,x,y,dy,iflag);
  19.       u1 = ydopri(x[0],y[0],dy[0],iflag);
  20.       u2 = ydopri(x[1],y[1],dy[1],iflag);
  21.       u3 = ydopri(x[2],y[2],dy[2],iflag);
  22.       u4 = ydopri(x[3],y[3],dy[3],iflag);
  23.  
  24.       for(i = 3; i <= iter-2; i++)
  25.       {
  26.        dyi = dy[i] +h * (55.*u4-59. * u3+37. * u2-9. * u1)/24.;
  27.        yi = y[i] +h*(55.*dy[i]-59.*dy[i-1]+37.*dy[i-2]-9.*dy[i-3])/24.;
  28.        u1 = u2;
  29.        u2 = u3;
  30.        u3 = u4;
  31.        u4 = ydopri(x[i+1],yi,dyi,iflag);
  32.        dy[i+1] = dy[i] +h*(9.*u4+19.*u3-5.*u2+u1)/24.;
  33.        y[i+1] = y[i]+h*(9.*dy[i+1] +19.*dy[i]-5.*dy[i-1]+dy[i-2])/24.;
  34.        u4 = ydopri(x[i+1],y[i+1],dy[i+1],iflag);
  35.        d[i+1] = -19. * (y[i+1]-yi) / 270.;
  36.       }
  37.     }
  38.