home *** CD-ROM | disk | FTP | other *** search
/ The C Users' Group Library 1994 August / wc-cdrom-cusersgrouplibrary-1994-08.iso / listings / v_09_05 / 9n05104a < prev    next >
Text File  |  1991-03-20  |  1KB  |  50 lines
  1.  
  2. // This function solves Kepler's equation by the Newton
  3. //   method.
  4. // Inputs:  m --  mean anomaly
  5. //          e0 -- starting estimate for eccen anom
  6. //          e  -- orbit eccentricity
  7. //          mi -- maximum number of iterations
  8. //          t  -- stopping tolerance
  9. // Outputs: sk returns the final estimate for eccen anom
  10. //   * Input m and e0 and output sk are in degrees. *    
  11.  
  12. double sk(double m,double e0,double e,int mi,double t)
  13. {
  14.   /* Solve Kepler's equation by the Newton method. */
  15.   double en,ep;
  16.   int i;
  17.  
  18.   e0 = e0 * DTR;
  19.   m  = m  * DTR;
  20.  
  21.   i  = 0;
  22.   en = e0;
  23.   do {
  24.     i++;
  25.     ep = en;
  26.  
  27.     // The next line of code is the Newton iteration:
  28.     //   x(n+1) = x(n) - (1/(dy/dx)) * y(n). 
  29.  
  30.     en = en - (en - e*sin(en) - m)/(1 - e*cos(en))   ;
  31.     en = fmod(en,TWOPI);
  32.     if ( fabs(en-ep) < t) break;
  33.   } while (i <= mi );
  34.   return(en*RTD);
  35. }
  36.  
  37.  
  38.                        SAMPLE INPUT/OUTPUT
  39. Initial eccentric anom = 0 
  40. Input mean anom        = 286.879298 
  41. Input maximum its      = 15 
  42. Input stop tolerance   = 0.000001 
  43. Input eccentricity     = 0.100000 
  44. Final eccentric anom   = 281.260008 
  45. Mean anom computed from solution = 286.879298 
  46.  
  47. This solution took five iterations; the solution E was
  48. checked by plugging it back into Kepler's equation
  49. and computing M.
  50.