5.1 Regla del Trapecio para Integración numérica con Python

Referencia: Chapra 21.1 p621 pdf645, Burden 4.2 p192 pdf202, Rodríguez 7.1.1 p273

La integración numérica con la regla del trapecio usa un polinomio de primer grado como aproximación de f(x) entre los extremos a y b del intervalo.

Es la primera de las formulas cerradas de Newton-Cotes.

I = \int_a^b f(x) dx \cong \int_a^b f_1 (x) dx

usando el rango entre [a,b] el polinomio se aproxima con una línea recta:

f_1 (x) = f(a) + \frac{f(b)-f(a)}{b-a} (x-a)

el área bajo esta línea recta es una aproximación de la integral de f(x) entre los límites a y b.

El resultado del integral es la regla del trapecio:

I = (b-a) \frac{f(a)+f(b)}{2}

que se interpreta como la multiplicación entre la base y altura promedio de un trapecio. También se llega al resultando sumando las áreas de sus componentes: rectángulo y triángulo.

El error de truncamiento se encuentra como el integral del término que le sigue al polinomio de Taylor en la aproximación, es decir el de grado 2, que al integrarlo tiene un orden de h³. (Rodríguez p275)

error_{truncar} = -\frac{h^3}{12}f''(z)

a < z < b

Algoritmo en Python

Ejemplo de Integración con el método del trapecio

Para integrar la función en el intervalo [1,3] con 4, 16, 32 ,64 y 128 tramos,

f(x)= \sqrt {(x)} \sin(x)

Tramos = 4

La base de los trapecios (x_i-x_i+1), si esta base tiene valores equidistantes se sustituye por la variable h .

base =\frac{b-a}{tramos}=\frac{3-1}{4}= 0.5 = h

para cada tramo, el área de un trapecio es:

A_i = base\frac{f(x_i)+f(x_{i+1})}{2}

por lo que cada trapecio en cada subintervalo tendrá un área de:

A_1 = (0.5)\frac{f(1)+f(1.5)}{2}

A_2 = (0.5)\frac{f(1.5)+f(2)}{2}

A_3 = (0.5)\frac{f(2)+f(2.5)}{2}

A_4 = (0.5)\frac{f(2.5)+f(3)}{2}

y el integral de todo el intervalo es la suma de todos los trapecios.

Integral = A_1+ A_2+ A_3+A_4

si se quiere tomar el factor común entre las sumas de áreas de cada tramo, también se tiene la forma de la ecuación con un h equidistante entre cada tramo del intervalo.

Integral = \frac{0.5}{2}\Big(f(1)+2f(1.5) +2f(2)+2f(2.5)+f(3)\Big)

Interpretando la fórmula, el integral es la suma de:

una vez la función evaluada en cada extremo,
mas dos veces cada valor de la función evaluada en los puntos intermedios.
El resultado de la suma se multiplica por h/2.

En caso que los puntos no sean equidistantes la simplificación NO procede, y se sumaran las áreas de los trapecios de cada tramo.

Siguiendo el procedimiento, para cada cantidad de tramos en el intervalo, los resultados serán:

tramos: 4
Integral: 1.99841708623

tramos:  16
Integral:  2.05019783717

tramos:  32
Integral:  2.05277225085

tramos:  64
Integral:  2.05341563776

tramos: 128
Integral:  2.05357647096

Algoritmo en Python – Regla o Método del trapecio

Algoritmo con fórmula simplificada para varios tramos, con puntos de muestras equidistantes h.

# Integración: Regla de los trapecios
# Usando una función fx()
import numpy as np
import matplotlib.pyplot as plt

# INGRESO
fx = lambda x: np.sqrt(x)*np.sin(x)

# intervalo de integración
a = 1
b = 3
tramos = 4

# PROCEDIMIENTO
# Regla del Trapecio
# Usando tramos equidistantes en intervalo
h = (b-a)/tramos
xi = a
suma = fx(xi)
for i in range(0,tramos-1,1):
    xi = xi + h
    suma = suma + 2*fx(xi)
suma = suma + fx(b)
area = h*(suma/2)

# SALIDA
print('tramos: ', tramos)
print('Integral: ', area)

La gráfica como ayuda didáctica, si la cantidad de tramos sea «poca», muestra los trapecios, si aumenta la cantidad de tramos, no es necesario poner líneas verticales en blanco para separar los trapecios:

# GRAFICA
# Puntos de muestra
muestras = tramos + 1
xi = np.linspace(a,b,muestras)
fi = fx(xi)
# Linea suave
muestraslinea = tramos*10 + 1
xk = np.linspace(a,b,muestraslinea)
fk = fx(xk)

# Graficando
plt.plot(xk,fk, label ='f(x)')
plt.plot(xi,fi, marker='o',
         color='orange', label ='muestras')

plt.xlabel('x')
plt.ylabel('f(x)')
plt.title('Integral: Regla de Trapecios')
plt.legend()

# Trapecios
plt.fill_between(xi,0,fi, color='g')
for i in range(0,muestras,1):
    plt.axvline(xi[i], color='w')

plt.show()

Para rellenar el área bajo la curva se usa la instrucción plt.fill_between(xi,0,fi, color=’g’). La división de los trapecios con una línea blanca (white) se realiza con la instrucción plt.axvline(xi[i], color=’w’).

Algoritmo con muestras

Algoritmo usando las muestras y para tramos que pueden tener distancias arbitrarias. En este caso se usa la formula del área del trapecio para cada tramo.

Usado para ejercicios donde los datos proporcionados son muestras. Se adapta el bloque de ingreso y el procedimiento inicia desde «Regla de Trapecio». La cantidad de tramos será el numero de muestras menos uno tramos = len(xi) - 1

# Integración: Regla de los trapecios
# Usando una muestras xi,fi
import numpy as np
import matplotlib.pyplot as plt

def integratrapecio_fi(xi,fi):
    ''' sobre muestras de fi para cada xi
        integral con método de trapecio
    '''
    n = len(xi)
    suma = 0
    for i in range(0,n-1,1):
        dx = xi[i+1]-xi[i]
        untrapecio = dx*(fi[i+1]+fi[i])/2
        suma = suma + untrapecio
    return(suma)

# PROGRAMA -----------------
# INGRESO
xi = [1. , 1.5, 2. , 2.5, 3. ]
fi = [0.84147098, 1.22167687, 1.28594075,
      0.94626755, 0.24442702]

# PROCEDIMIENTO
Area = integratrapecio_fi(xi,fi)

# SALIDA
print('tramos: ',len(xi)-1)
print('Integral con trapecio: ',Area)