Write functions to calculate the definite integral of a function  ƒ(x)  using all five of the following methods:

 rectangular
 left
 right
 midpoint
 trapezium
 Simpson's
 composite
Your functions should take in the upper and lower bounds (a and b), and the number of approximations to make in that range (n).

Assume that your example already has a function that gives values for  ƒ(x) .

Simpson's method is defined by the following pseudo-code:


 Pseudocode: Simpson's method, composite



procedure quad_simpson_composite(f, a, b, n)
    h := (b - a) / n
    sum1 := f(a + h/2)
    sum2 := 0

    loop on i from 1 to (n - 1)
        sum1 := sum1 + f(a + h * i + h/2)
        sum2 := sum2 + f(a + h * i)
                                            
    answer := (h / 6) * (f(a) + f(b) + 4*sum1 + 2*sum2)




Demonstrate your function by showing the results for:

   ƒ(x) = x3,       where   x   is     [0,1],       with           100 approximations.   The exact result is     0.25               (or 1/4)
   ƒ(x) = 1/x,     where   x   is   [1,100],     with        1,000 approximations.   The exact result is     4.605170+     (natural log of 100)
   ƒ(x) = x,         where   x   is   [0,5000],   with 5,000,000 approximations.   The exact result is   12,500,000
   ƒ(x) = x,         where   x   is   [0,6000],   with 6,000,000 approximations.   The exact result is   18,000,000


See also

   Active object for integrating a function of real time.
   Special:PrefixIndex/Numerical integration for other integration methods.


def leftrect(f, left, right)
  f.call(left)
end
 
def midrect(f, left, right)
  f.call((left+right)/2.0)
end
 
def rightrect(f, left, right)
  f.call(right)
end
 
def trapezium(f, left, right)
  (f.call(left) + f.call(right)) / 2.0
end
 
def simpson(f, left, right)
  (f.call(left) + 4*f.call((left+right)/2.0) + f.call(right)) / 6.0
end
 
def integrate(f, a, b, steps, method)
  delta = 1.0 * (b - a) / steps
  total = 0.0
  steps.times do |i|
    left = a + i*delta
    right = left + delta
    total += delta * send(method, f, left, right)
  end
  total
end
 
def square(x)
  x**2
end
 
def def_int(f, a, b)
  l = case f.to_s
      when /sin>/
        lambda {|x| -Math.cos(x)}
      when /square>/
        lambda {|x| (x**3)/3.0}
      end
  l.call(b) - l.call(a)
end
 
a = 0
b = Math::PI
steps = 10
 
for func in [method(:square), Math.method(:sin)]
  puts "integral of #{func} from #{a} to #{b} in #{steps} steps"
  actual = def_int(func, a, b)
  for method in [:leftrect, :midrect, :rightrect, :trapezium, :simpson]
    int = integrate(func, a, b, steps, method)
    diff = (int - actual) * 100.0 / actual
    printf "   %-10s  %s\t(%.1f%%)\n", method, int, diff
  end
end