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

Referencia: Chapra 21.1 p621 pdf645, Burden 4.2 p192 pdf202, Rodriguez 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.

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 fiene un orden de h3. (Rodriguez 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 (xi-xi+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 conceptual en Python

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 incluso muestras arbitrariamente espaciadas
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
# Puntos de muestra
muestras = tramos + 1
xi = np.linspace(a,b,muestras)
fi = fx(xi)

# Regla del Trapecio
# Usando puntos muestreados
# incluso arbitrariamente espaciados
suma = 0
for i in range(0,tramos,1):
    dx = xi[i+1]-xi[i]
    Atrapecio = dx*(fi[i]+fi[i+1])/2
    suma = suma + Atrapecio
integral = suma

# SALIDA
print('tramos: ', tramos)
print('integral: ', integral)