complexsets
complexsets

Reputation: 15

How to compute an integral of a given function without SymPy?

I'm having difficulties creating a simple for loop code to integrate a given function without using SymPy. I'm thinking using some sort of Riemann approximation but I don't know how to do this exactly. The code I have so far is:

def xsquared(x):
    n = 2
    return x**n

def integral(fun, xmin, xmax):
    total = 0
    for a in range(xmin, xmax):
        x = a
        total += fun(x*1.235)
    return total 

print(integral(xsquared, 0, 4))

The output gives 21.3 but how do I do this without inputting a number referring to the "fun(x*1.235)" part?

Any help would be appreciated.

Upvotes: 1

Views: 123

Answers (1)

Cheche
Cheche

Reputation: 1516

To compute integral by Riemann means that you are computing the limit of Riemann's sum by making partitions innner (Wikipedia).

For that, your function could have a new parameter interval: dx or maybe you could just guess the interval by splitting the full range into N equal sized intervals. Here is an example asking dx as argument.

Then, your function should be:

def riemann(fun, xmin, xmax, dx):
    total = 0
    a = xmin
    while a < xmax:
        total += fun(a + dx/2)*dx
        a += dx
    return total

Example Outputs

print(riemann(xsquared, 0, 4, 0.1))
> 21.330000000000013

print(riemann(xsquared, 0, 4, 0.25))
> 21.3125

print(riemann(xsquared, 0, 4, 0.5))
> 21.25

Analythic resolution gives: 64/3 ~ 21.33333

You're aproximating then the integral by computing the area of the rectangle having:

  • height: the function value at interval middle point fun(a + dx/2)
  • width: the interval length (dx)

Note: if xmax < xmin, you should verify that dx < 0.

Upvotes: 1

Related Questions