Integral by the Trapezoidal rule algorithm I wrote is giving a 10 times larger answer

254 views Asked by At

I made a function in Python which calculates a definite integral according to the Trapezoidal rule: Trapezoidal rule formula

That's the code:

from math import ceil

def Trapez_rule(f, a, b, n):
    '''Calculates an estimation of a definite integral of a function f(x), between the boundries a, b, by dividing the area to n equal areas'''
    sum = (f(a) + f(b)) / 2
    for i in range(ceil((b * n))):
        sum += f(a + i / n)
    sum *= (b - a) / n
    return sum

The answer it gives is 10 times higher that it should have returned. I can't find the source of the problem.

2

There are 2 answers

0
minseong On

I went ahead and fixed your code, and also renamed the function to fit with the official style guide PEP-8.

def trapezium_rule_integral(f, a, b, n):
    '''Calculates an estimate of the definite integral of a function f, between
       the boundaries a and b, by dividing the area to n equal areas'''
    height = (b - a) / n
    x = a
    ys = []
    while x <= b:
        ys.append(f(x))
        x += height
    estimate = 0.5 * height * ( (ys[0] + ys[-1]) + 2 * (sum(ys[1:-1])) )
    return estimate
3
Boendal On

Assume:

a=10
b=20
n=5

These lines are the problem:

for i in range(ceil((b * n))):
    sum += f(a + i / n)

i go from 0 to 99
when i = 99 then:

f(a + i / n) => f(10 + 99/5) => f(29)
  1. You divide two ints 99/5 => 29 and not 29.8.
  2. But you only want to have it in range from 10 to 20.
  3. You use n false look at the post solution below, so this should work:

    def Trapez_rule(f, a, b, n): h = (b-a) / float(n) sum = (f(a) + f(b)) / 2.w for i in range(1,n-1): sum += f(a + i * h) sum *= h return sum