Trapezoidal Rule in F#

Numerical is approximate; analytical is exact.

A numerical method is a tool used to find a numerical solution to a problem which can be difficult to solve analytically.

Mostly done by computers, this approach is a complete set of instructions that are executed in order to find a solution. With this being said,
the Trapezoidal Rule is a numerical method used to find definite integrals. The most common application for this method is to find the
area under a function in a given set of intervals. Finding this area is typically done using Riemann sums, but this method proves to be more exact because it uses trapezoids instead of rectangles.

• Reimann sum
• Trapezoidal Rule

In this case, I implemented this numerical method using F# which is a functional programming language. First, I chose the function to integrate; in this case sin(x):

let f x = System.Math.Sin(x)

The algorithm has to receive the intervals (a,b) through which it will find the area and a number (N) of trapezoids to form.
With these parameters, it can now calculate the length of each trapezoid:

let h = (b - a)/N

• By Mkwadee at English Wikipedia, CC BY-SA 3.0, Link

The formula for this method is the following:

I created two functions in order to calculate the formula above:

let rec sum_recursive (xi:float) (h:float) f (N:float) (iter:float) =
    if iter = N then
      xi
    else
      let xi_1 = xi + h*f(h*iter)
      sum_recursive xi_1 h f N (iter+1.0)

let trapezoidal (a:float) (b:float) f (N:float) =
      let h = (b - a)/N
      let left_side = h/2.0*(f(a)+f(b))
      left_side + sum_recursive 0.0 h f N 0.0

This method can be used for things like: finding how much fluid is flowing through a specific area, calculating the work needed to move an object at a certain speed, or even finding the center of mass of an object.

Overall, learning to program a numerical method has enabled me to relate mathematical concepts, that are normally just seen in class; with programming, which is something I actually use in my day to day life. Furthermore, this was also my first time working with a functional programming language such as F#; which I ended up enjoying.

Here's my complete code
Hope you find it useful!

References

Understanding the trapezoidal rule
Trapezoidal Rule
How to approximate area with the trapezoid rule

Comments

[Bohdan Stupak]

Great job but I'm not sure why using functional-first language to use it OO-style. To me using mutable in F# is a flag that something is done in a non-functional fashion.
How about considering the code below?

let trapezoidal2 a b N = 
    let h = (b - a)/N
    let initialAccValue = h/2.0*(f(a)+f(b))
    [|1 .. System.Convert.ToInt32(N)|]
    |> Array.mapi(fun i x -> i+1,x)
    |> Array.fold (fun acc (i,_) -> acc + h*f(float(i)*h+a)) initialAccValue     

The key tricks here are Array.mapi which allows us to work with tuple (element,index) instead of just an element; and substitution of iterative increment of xi by h to multiplication index*h+a

[Maria Ines Montenegro]

Since this was my first time with a functional programming language I think I hadn't fully grasped the concepts. I ended up using recursion in order to remove the mutable variables; which I updated in my post. Thank you very much for your feedback!

[Ashley M Campbell]

This is so true. The more trapezoids you use the more exact the integral tends to be. Also, it is very simple to use compared to other methods. Great article!

[samantha]

Really cool approach