home *** CD-ROM | disk | FTP | other *** search
/ Otherware / Otherware_1_SB_Development.iso / amiga / misc / icalc.lha / ICalc / Scripts / trapzd.ic < prev    next >
Encoding:
Text File  |  1992-08-03  |  1.0 KB  |  40 lines
  1. #
  2. # trapzd
  3. #
  4. # Routine from "Numerical Recipes in C" by Press et. al.
  5. # It is used by qtrap.ic, qsimp.ic and qromb.ic.
  6. #
  7. # This routine computes the n'th stage of refinement of an extended
  8. # trapezoidal rule. f is a function of global variable x (e.g. x*x, sin(x))
  9. # to be integrated between the limits a and b. When called with n=1, the
  10. # routine returns the crudest estimate of the integral. Subsequent calls
  11. # with n=2,3,... (in that sequential order) will improve the accuracy by
  12. # adding 2^(n-2) additional interior points.
  13. # 
  14. trapzd_s = trapzd_it = 0    # create required globals
  15. #
  16. func trapzd(~f,a,b,n) = {
  17.     local tmp,tnm,sum,del,j
  18.  
  19.     if (n == 1) {
  20.         trapzd_it = 1
  21.         x = a        # set x to value a
  22.         tmp = f        # evaluate f at a
  23.         x = b
  24.         tmp += f 
  25.         return trapzd_s = 0.5*(b-a)*tmp
  26.     } else {
  27.         tnm = trapzd_it
  28.         del = (b-a)/tnm
  29.         x = a + 0.5*del
  30.         sum = 0
  31.         for (j = 1; j <= trapzd_it; j += 1) {
  32.             sum += f    # i.e. sum += f evaluated at x
  33.             x += del
  34.         }
  35.         trapzd_it *= 2
  36.         trapzd_s = 0.5*(trapzd_s + (b-a)*sum/tnm)
  37.         return trapzd_s
  38.     }
  39. }
  40.