[Implementação] Integração aproximada pelo Método do Trapézio - Python

Este artigo é uma continuação da discussão teórica feita aqui. Se você ainda não leu, por favor, leia. Aqui vamos simplesmente implementar, em Python, tudo o que foi discutido lá.

Começando do começo

Conforme já discutido aqui (se você não leu, é a última chance), a equação para a integração aproximada pelo Método do Trapézio é a seguinte:

$\int\limits_a^b {f(x)dx} \approx \frac{\Delta x}{2}\sum\limits_{i=1}^n {[f(x_i) + f(x_{i-1})]}$

E a integral que queremos calcular é a seguinte:

$\int\limits_0^1 {\left(x^5+27x^4-\frac{13}{4}x^3+\frac12 x^2 + 37x + 15\right) dx}$

Primeiramente, vamos definir essa função (para facilitar nosso trabalho):

def f(x):
    return x**5 + 27*(x**4) - ((13/4) * (x**3)) + ((0.5)*(x**2)) + (37*x) + 15

Calculando a integral

Vamos, agora, implementar a função que, de fato, calcula a integral. Primeiro, vamos definir $x_i$. Vimos que $x_0 = a$, então $x_1 = a + \Delta x$ e $x_i = a + i\Delta x$.
No código, usaremos as seguintes variáveis:

A = área
Dx = $\Delta x$ = $\frac{b-a}{n}$
xi = $x_i$
xi1 = $x_{i-1}$
a = $a$ = limite inferior de integração
b = $b$ = limite superior de integração

Nos basta, então, seguir a equação que temos. Primeiro vou definir $\Delta x$ e algumas variáveis e, então, calcularemos o somatório (veja que, primeiro, defini $x_i$ e $x_{i-1}$ e depois somei chamando a função $f$).

def integral(a, b, n):
    A = 0
    Dx = (b-a)/n
    for i in range(1, n+1):
        xi  = a + Dx*i
        xi1 = a + Dx*(i-1)
        A += f(xi) + f(xi1)
    A *= (Dx/2)
    return A

Por fim, basta multiplicar o somatório por $\frac{\Delta x}{2}$.
Agora que temos a função de integração, vamos receber $a$, $b$ e $n$ do usuário e exibir o resultado:

a = float(input("Por favor, insira o limite inferior a = "))
b = float(input("Por favor, insira o limite superior b = "))
n = int(input("Por favor, insira o a quantidade de iterações desejada: n = "))
print(integral(a,b,n))

E, pronto, temos o programa finalizado (código completo no final). Apenas para efeito de teste, integrei no intervalo de 0 a 1 com $n=1$, $n=10$, $n=100$ e $n=1000$. A integral exata (obtida através da primitiva) também é mostrada abaixo.

$\frac{1}{2}\sum\limits_{i=1}^1 {[f(x_i) + f(x_{i-1})]} = 46.125$

$\frac{1}{20}\sum\limits_{i=1}^{10} {[f(x_i) + f(x_{i-1})]} = 38.50761$

$\frac{1}{200}\sum\limits_{i=1}^{100} {[f(x_i) + f(x_{i-1})]} = 38.42170207350001$

$\frac{1}{2000}\sum\limits_{i=1}^{1000} {[f(x_i) + f(x_{i-1})]} = 38.42084202083233$

$\int\limits_0^1 {f(x)dx} = \frac{9221}{240} \approx 38.42083333333333$

Concluímos que conseguimos uma aproximação excelente com esse método. Além disso, o cálculo em Python se mostrou muito mais preciso que no Portugol Studio.

Veja também a implementação em Portugol Studio.

Código final:

raw

def f(x):
    return x**5 + 27*(x**4) - ((13/4) * (x**3)) + ((0.5)*(x**2)) + (37*x) + 15

def integral(a, b, n):
    A = 0
    Dx = (b-a)/n
    for i in range(1, n+1):
        xi  = a + Dx*i
        xi1 = a + Dx*(i-1)
        A += f(xi) + f(xi1)
    A *= (Dx/2)
    return A

a = float(input("Por favor, insira o limite inferior a = "))
b = float(input("Por favor, insira o limite superior b = "))
n = int(input("Por favor, insira o a quantidade de iterações desejada: n = "))
print(integral(a,b,n))