Métodos numéricos Competencia 4 Integración y derivación

4. DIFERENCIACIÓN E INTEGRACIÓN NUMÉRICA

En un primer acercamiento a un curso de cálculo, se aprende a calcular derivadas e integrales de

forma analítica, cuando las expresiones aparecen explícitamente con una regla de correspondencia
de la forma 𝑦 = 𝑓(𝑥). Sin embargo, en la gran mayoría de los problemas de aplicación esto no pasa,
más bien se presenta la función como una tabla de valores.

Es en estos casos cuando el uso de los métodos clásicos de integración y derivación pierden sentido
y cobran importancia los métodos numéricos.

4.1 DERIVACIÓN NUMÉRICA

La forma clásica de representar una derivada numérica es mediante el método de diferencias finitas,
una diferencia finita se obtiene de quitar el límite al cociente de Newton para el cálculo de derivadas.
Si presentamos la derivada de una función mediante definición podemos escribir:
𝑑𝑦 𝑦(𝑥 + ℎ) − 𝑦(𝑥)
= (1. 𝑎)
𝑑𝑥 ℎ
𝑑𝑦 𝑦(𝑥) − 𝑦(𝑥 − ℎ)
= (1. 𝑏)
𝑑𝑥 ℎ

Ambas, (1.a) y (1.b) nos ayudan a obtener la misma derivada, la única diferencia, radica en cómo se
toma el incremento h. En (1.a) se toma el incremento h hacia adelante, mientras que en (1.b) se
toma hacia atrás. Si quitamos el límite de esas expresiones obtenemos nuestra primera
aproximación a la derivada, esto queda:
𝑑𝑦 𝑦(𝑥 + ℎ) − 𝑦(𝑥)
≈ (2. 𝑎)
𝑑𝑥 ℎ
𝑑𝑦 𝑦(𝑥) − 𝑦(𝑥 − ℎ)
≈ (2. 𝑏)
𝑑𝑥 ℎ

A la expresión (2.a) se le llama diferencia finita hacia delante, mientras que a (2.b) se le llama
diferencia finita hacia atrás. Ambas ecuaciones nos ayudan a calcular la tasa de variación de “y” con
respecto de “x”, la única diferencia radica en que la primera considera los datos hacia adelante
mientras que la segunda los datos hacia atrás. Para reducir el error, se suele tomar el promedio
entre ambas, calculándose, de esa manera, una “diferencia centra”. Promediando (2.a) y (2.b)
obtenemos:

𝑦(𝑥 + ℎ) − 𝑦(𝑥) 𝑦(𝑥) − 𝑦(𝑥 − ℎ)

( )+( ) 𝑦(𝑥 + ℎ) − 𝑦(𝑥 − ℎ)
ℎ ℎ
=
2 2ℎ
𝑦(𝑥 + ℎ) − 𝑦(𝑥 − ℎ)
=
2ℎ
La expresión final queda:
𝑑𝑦 𝑦(𝑥 + ℎ) − 𝑦(𝑥 − ℎ)
≈ (3)
𝑑𝑥 2ℎ

1
Métodos numéricos Competencia 4 Integración y derivación

Para calcular la segunda derivada sacamos la derivada de la primera derivada:

𝑑2 𝑦 𝑑 𝑑𝑦 (4)
= ( )
𝑑𝑥 2 𝑑𝑥 𝑑𝑥

Podemos calcularlas igual hacia adelante y hacia atrás de la siguiente forma:

𝑑𝑦 𝑑𝑦
𝑑2 𝑦 𝑑𝑥 |𝑥+ℎ − 𝑑𝑥 |𝑥
= (5. 𝑎)
𝑑𝑥 2 ℎ
𝑑𝑦 𝑑𝑦
𝑑2 𝑦 𝑑𝑥 |𝑥 − 𝑑𝑥 |𝑥−ℎ (5. 𝑏)
=
𝑑𝑥 2 ℎ

Para aproximar, quitamos el límite y nos quedamos con una diferencia finita:
𝑑𝑦 𝑑𝑦
𝑑2 𝑦 𝑑𝑥 |𝑥+ℎ − 𝑑𝑥 |𝑥
≈ (6. 𝑎)
𝑑𝑥 2 ℎ
𝑑𝑦 𝑑𝑦
𝑑2 𝑦 𝑑𝑥 |𝑥 − 𝑑𝑥 |𝑥−ℎ (6. 𝑏)
≈
𝑑𝑥 2 ℎ

Reemplazamos ahora (2.a) y (2.b) y obtenemos las derivadas segundas hacia delante y hacia atrás.
Para las derivadas hacia delante obtenemos:

𝑦((𝑥 + ℎ) + ℎ) − 𝑦(𝑥 + ℎ) 𝑦(𝑥 + ℎ) − 𝑦(𝑥)

𝑑𝑦 𝑑𝑦 ( )−( )
𝑑2 𝑦 𝑑𝑥 |𝑥+ℎ − 𝑑𝑥 |𝑥 ℎ ℎ
≈ =
𝑑𝑥 2 ℎ ℎ
𝑦(𝑥 + 2ℎ) − 𝑦(𝑥 + ℎ) − 𝑦(𝑥 + ℎ) + 𝑦(𝑥) 𝑦(𝑥 + 2ℎ) − 2𝑦(𝑥 + ℎ) + 𝑦(𝑥)
= =
ℎ2 ℎ2
Por otro lado, la segunda derivada hacia atrás queda:

𝑑𝑦 𝑑𝑦 𝑦(𝑥) − 𝑦(𝑥 − ℎ) 𝑦(𝑥 − ℎ) − 𝑦((𝑥 − ℎ) − ℎ)

𝑑2 𝑦 𝑑𝑥 |𝑥 − 𝑑𝑥 |𝑥−ℎ ( ℎ
)−(
ℎ
)
≈ =
𝑑𝑥 2 ℎ ℎ
𝑦(𝑥) − 𝑦(𝑥 − ℎ) − 𝑦(𝑥 − ℎ) + 𝑦(𝑥 − 2ℎ) 𝑦(𝑥) − 2𝑦(𝑥 − ℎ) + 𝑦(𝑥 − 2ℎ)
= =
ℎ2 ℎ2

2
Métodos numéricos Competencia 4 Integración y derivación

Finalmente obtenemos lo siguiente:

𝑑2 𝑦 𝑦(𝑥 + 2ℎ) − 2𝑦(𝑥 + ℎ) + 𝑦(𝑥)

≈ (7. 𝑎)
𝑑𝑥 2 ℎ2
2
𝑑 𝑦 𝑦(𝑥) − 2𝑦(𝑥 − ℎ) + 𝑦(𝑥 − 2ℎ)
≈ (7. 𝑏)
𝑑𝑥 2 ℎ2

La expresión (7.a) sirve para calcular las derivadas segundas hacia delante, mientras que (7.b) las
calcula hacia atrás.

Para la derivada segunda central, se calcula mediante el cociente de la primera derivada hacia
adelante y la primera derivada hacia atrás:

𝑑𝑦 𝑑𝑦 𝑦(𝑥 + ℎ) − 𝑦(𝑥) 𝑦(𝑥) − 𝑦(𝑥 − ℎ)

𝑑2 𝑦 𝑑𝑥 |𝑥+ℎ − 𝑑𝑥 |𝑥−ℎ ( ℎ
)−(
ℎ
)
≈ =
𝑑𝑥 2 ℎ ℎ
𝑦(𝑥 + ℎ) − 𝑦(𝑥) − 𝑦(𝑥) + 𝑦(𝑥 − ℎ)
=
ℎ2
Finalmente obtenemos la expresión de la segunda derivada central:

𝑑2 𝑦 𝑦(𝑥 + ℎ) − 2𝑦(𝑥) + 𝑦(𝑥 − ℎ) (8)

≈
𝑑𝑥 2 ℎ2

Siguiendo este procedimiento, se puede calcular cualquier derivada de cualquier orden natural.

En resumen, las derivadas primera y segunda son:

Primera derivada Segunda derivada

𝑑2 𝑦 𝑦(𝑥 + 2ℎ) − 2𝑦(𝑥 + ℎ) + 𝑦(𝑥)
𝑑𝑦 𝑦(𝑥 + ℎ) − 𝑦(𝑥) ≈
Hacia adelante ≈ 𝑑𝑥 2 ℎ2
𝑑𝑥 ℎ
𝑑2 𝑦 𝑦(𝑥) − 2𝑦(𝑥 − ℎ) + 𝑦(𝑥 − 2ℎ)
𝑑𝑦 𝑦(𝑥) − 𝑦(𝑥 − ℎ) ≈
Hacia atrás ≈ 𝑑𝑥 2 ℎ2
𝑑𝑥 ℎ
𝑑2 𝑦 𝑦(𝑥 + ℎ) − 2𝑦(𝑥) + 𝑦(𝑥 − ℎ)
𝑑𝑦 𝑦(𝑥 + ℎ) − 𝑦(𝑥 − ℎ) ≈
Central ≈ 𝑑𝑥 2 ℎ2
𝑑𝑥 2ℎ
Si los datos a analizar son los que estén por delante del punto 𝑥0 resultará más conveniente trabajar
diferencias hacia delante. Por otro lado, si los datos están por detrás de 𝑥0 , nos convendrá mejor
trabajar con diferencias hacía atrás. En el caso de que existan datos hacia delante y hacia atrás,
resultará mejor trabajar diferencias centrales.

En el libro de métodos numéricos aplicados con Software de Shoichiro Nakamura, podrá encontrar
más expresiones para calcular derivadas numéricamente. Más concretamente en las páginas 161 y
162 del libro (181 y 182 del PDF) hay fórmulas que fueron deducidas mediante polinomios de Taylor
que utilizan más de un par de puntos. https://evflores.files.wordpress.com/2014/02/metodos-
numericos-aplicados-con-software-by-sholchlro-nakamura1.pdf

3
Métodos numéricos Competencia 4 Integración y derivación

Ejemplo 4.1. Calcular la segunda derivada de la función 𝑦 =ln 𝑥 en el punto 𝑥 = 1, considere un

valor ℎ = 0.25 y derivada central

𝑑2 𝑦 𝑦(𝑥 + ℎ) − 2𝑦(𝑥) + 𝑦(𝑥 − ℎ)

≈
𝑑𝑥 2 ℎ2
𝑑2 𝑙𝑛 𝑙𝑛 𝑥 𝑙𝑛 𝑙𝑛 (1 + 0.25) − 2 𝑙𝑛 (1) + 𝑙𝑛 (1 − 0.25)
|𝑥=1 ≈ = −1.03261
𝑑𝑥 0.252
Ejemplo 4.2. Calcular numéricamente la derivada parcial dada en el punto (𝑥, 𝑦) = (1,2), utilice
diferencias hacia delante, hacia atrás y centrales, suponga un valor h=0.5:

Hacia Delante

𝜕(3𝑥 2 𝑦) 𝑑(6𝑥 2 ) 6(1 + 0.5)2 − 6(1)2

|(1,2) = |𝑥=1 ≈ = 15
𝜕𝑥 𝑑𝑥 0.5
Hacia Atrás

𝜕(3𝑥 2 𝑦) 𝑑(6𝑥 2 ) 6(1)2 − 6(1 − 0.5)2

|(1,2) = |𝑥=1 ≈ =9
𝜕𝑥 𝑑𝑥 0.5
Central

𝜕(3𝑥 2 𝑦) 𝑑(6𝑥 2 ) 6(1 + 0.5)2 − 6(1 − 0.5)2

|(1,2) = |𝑥=1 ≈ = 12
𝜕𝑥 𝑑𝑥 2(0.5)

Ejemplo 4.3. Calcular la primera derivada de 𝑦 = 𝑥𝑒 𝑥 , usando diferencias centrales, hacia delante
y hacia atrás, evaluada en 𝑥 = 2 con un valor de ℎ = 0.01:
Centrales

𝑑𝑥𝑒 𝑥 (2 + 0.01)𝑒 (2+0.01) − (2 − 0.01)𝑒 2−0.01

| ≈ = 22.1678
𝑑𝑥 𝑥=2 2 ∗ (0.01)
Hacia delante

𝑑𝑥𝑒 𝑥 (2 + 0.01)𝑒 (2+0.01) − 2𝑒 2

|𝑥=2 ≈ = 22.3155
𝑑𝑥 0.01
Hacia atrás

𝑑𝑥𝑒 𝑥 2𝑒 2 − (2 − 0.01)𝑒 2−0.01

|𝑥=2 ≈ = 22.020
𝑑𝑥 0.01
Valor real

𝑑2 𝑥𝑒 𝑥
| = 22.1671683..
𝑑𝑥 𝑥=2
**La mejor aproximación para este ejercicio es la central

4
Métodos numéricos Competencia 4 Integración y derivación

Ejemplo 4.4. Dada la siguiente tabla de valores calcular la primera derivada en x=0.5 usando
derivada central, hacía adelante y hacia atrás:

x y
0 1
0.1 0.99500417
0.2 0.98006658
0.3 0.95533649
0.4 0.92106099
0.5 0.87758256
0.6 0.82533561
0.7 0.76484219
0.8 0.69670671
0.9 0.62160997
1 0.54030231
1.1 0.45359612
1.2 0.36235775
Calculamos la central en 0.5
𝑦(𝑥 + ℎ) − 𝑦(𝑥 − ℎ)
𝑦´ =
2ℎ
𝑥 = 0.5, ℎ = 0.1
𝑦(0.5 + 0.1) − 𝑦(0.5 − 0.1) 𝑦(0.6) − 𝑦(0.4) 0.82533561 − 0.92106099
𝑦´(0.5) = = =
0.2 0.2 0.2
= −0.4786269
Calculamos la derivada hacia delante:
𝑦(𝑥 + ℎ) − 𝑦(𝑥)
𝑦´ =
ℎ
𝑥 = 0.5, ℎ = 0.1
𝑦(0.5 + 0.1) − 𝑦(0.5) 𝑦(0.6) − 𝑦(0.5) 0.82533561 − 0.87758256
𝑦´ = = = = −0.5224695
0.1 0.1 0.1
Calculamos la derivada hacia atrás:
𝑦(𝑥) − 𝑦(𝑥 − ℎ)
𝑦´(𝑥) =
ℎ
𝑥 = 0.5, ℎ = 0.1
𝑦(0.5) − 𝑦(0.4) 0.87758256 − 0.92106099
𝑦´(𝑥) = = = −0.4347843
0.1 0.1

5
Métodos numéricos Competencia 4 Integración y derivación

Ejemplo 4.5. Deducir una expresión para calcular la derivada de tercer orden hacia delante
𝑦´´(𝑥 + ℎ) − 𝑦´´(𝑥)
𝑦´´´(𝑥) ≈
ℎ
Recordamos la fórmula de segunda derivada:

𝑑2 𝑦 𝑦(𝑧 + 2ℎ) − 2𝑦(𝑧 + ℎ) + 𝑦(𝑧)

≈
𝑑𝑧 2 ℎ2
Reemplazamos:
𝑦(𝑥 + ℎ + 2ℎ) − 2𝑦(𝑥 + ℎ + ℎ) + 𝑦(𝑥 + ℎ) 𝑦(𝑥 + 2ℎ) − 2𝑦(𝑥 + ℎ) + 𝑦(𝑥)
( )−( )
ℎ2 ℎ2
𝑦´´´(𝑥) ≈
ℎ
𝑦(𝑥 + ℎ + 2ℎ) − 2𝑦(𝑥 + ℎ + ℎ) + 𝑦(𝑥 + ℎ) − 𝑦(𝑥 + 2ℎ) + 2𝑦(𝑥 + ℎ) − 𝑦(𝑥)
𝑦´´´(𝑥) ≈
ℎ3
𝑦(𝑥 + 3ℎ) − 2𝑦(𝑥 + 2ℎ) + 𝑦(𝑥 + ℎ) − 𝑦(𝑥 + 2ℎ) + 2𝑦(𝑥 + ℎ) − 𝑦(𝑥)
𝑦´´´(𝑥) ≈
ℎ3
𝑦(𝑥 + 3ℎ) − 3𝑦(𝑥 + 2ℎ) + 3𝑦(𝑥 + ℎ) − 𝑦(𝑥)
𝑦´´´(𝑥) ≈
ℎ3
Ejemplo 4.6. De la tabla del ejemplo 4.4 calcular la tercera derivada hacia delante, utilizando la
fórmula obtenida en el problema 4.5
𝑦(𝑥 + 3ℎ) − 3𝑦(𝑥 + 2ℎ) + 3𝑦(𝑥 + ℎ) − 𝑦(𝑥)
𝑦´´´(𝑥) ≈
ℎ3
ℎ = 0.1, 𝑥 = 0.5
𝑦(0.8) − 3𝑦(0.7) + 3𝑦(0.6) − 𝑦(0.5)
𝑦´´´(0.5) ≈
(0.1)3
0.69670671 − 3 ∗ 0.76484219 + 3 ∗ 0.82533561 − 0.87758256
𝑦´´´(0.5) ≈ = 0.60441
(0.1)3

6
Métodos numéricos Competencia 4 Integración y derivación

Ejemplo 4.7. Calcular numéricamente la derivada direccional de: 𝑓(𝑥, 𝑦, 𝑧) = 𝑥 2 𝑦 2 𝑧 2 , en

𝑖+𝑗+𝑘
dirección del vector y en el punto (𝑥0 , 𝑦0 , 𝑧0 ) = (1, 1, 1). Pista: Recuerda que la derivada
√3
direccional se puede calcular de la siguiente manera: Dfu = ∇𝑓(𝑥0 , 𝑦0 , 𝑧0 ) ∙ 𝒖. Las derivadas
calcúlalas numéricamente con un paso h=0.001. Utilice derivada central

𝜕(𝑥 2 𝑦 2 𝑧 2 ) 𝜕(𝑥 2 𝑦 2 𝑧 2 ) 𝜕(𝑥 2 𝑦 2 𝑧 2 )

∇𝑓(1,1,1) = 𝑓𝑥 𝑖 + 𝑓𝑦 𝑗 + 𝑓𝑧 𝑘 = | 𝑖+ | 𝑗+ | 𝑘
𝜕𝑥 (1,1,1)
𝜕𝑦 (1,1,1)
𝜕𝑧 (1,1,1)

𝑑(𝑥 2 ) 𝑑(𝑦 2 ) 𝑑(𝑧 2 )

= | 𝑖+ | 𝑗+ | 𝑘
𝑑𝑥 𝑥=1 𝑑𝑦 𝑦=1 𝑑𝑧 𝑧=1

(1 + 0.001)2 − (1 − 0.001)2 (1 + 0.001)2 − (1 − 0.001)2 (1 + 0.001)2 − (1 − 0.001)2

= 𝑖+ 𝑗+ 𝑘
2 ∗ (0.001) 2 ∗ (0.001) 2 ∗ (0.001)

= 2𝑖 + 2𝑗 + 2𝑘

Se calcula el producto punto:

𝑖+𝑗+𝑘 2+2+2 6 6 √3 6√3

= ∇𝑓(𝑥0 , 𝑦0 , 𝑧0 ) ∙ 𝒖 = (2𝑖 + 2𝑗 + 2𝑘) ∙ ( )= = = =
√3 √3 √3 √3 √3 3

= 2√3

7
Métodos numéricos Competencia 4 Integración y derivación

4.2 INTEGRACIÓN NUMÉRICA

En esta sección, se presentan métodos para calcular integrales de manera numérica, se hará
principal énfasis en los métodos de la regla trapezoidal, las fórmulas de Simpson y técnicas para
calcular numéricamente integrales múltiples.

4.2.1 Regla del trapecio, integración lineal dado 2 puntos

Consiste en dividir la función en muchos trapecios pequeños y calcular su área de cada uno
de ellos para después sumarla:

Ilustración 1. Representación gráfica del método del trapecio

La forma de calcular el área de un trapecio es

𝑏+𝐵
𝐴=( )ℎ
2
Realizando la analogía:

𝑏 = 𝑓(𝑥1 )
𝐵 = 𝑓(𝑥2 )
ℎ = 𝑥2 − 𝑥1

8
Métodos numéricos Competencia 4 Integración y derivación

Entonces, el área de un trapecio es:

𝑓(𝑥1 ) + 𝑓(𝑥2 ) 𝑓(𝑥1 ) + 𝑓(𝑥2 )
𝐴12 = (𝑥2 − 𝑥1 ) = ℎ
2 2
Si dividimos la función, en “n” trapecios con longitud ℎ iguales, entonces, el área total se
calcula como sigue:

Primero Calculamos la longitud de h, para esto hay que considerar que para partir en 𝑛
intervalos iguales, necesitamos un total de 𝑛 + 1 puntos, si partimos desde el punto 𝑥𝑜 la
numeración iría: 𝑥𝑜 , 𝑥1 , 𝑥2 , … , 𝑥𝑛
𝑥𝑛 − 𝑥0
ℎ=
𝑛
𝑛

𝐴 𝑇 = ∑ 𝐴𝑘 = 𝐴1 + 𝐴2 + ⋯ + 𝐴𝑛
𝑘=1

𝑓(𝑥0 ) + 𝑓(𝑥1 ) 𝑓(𝑥1 ) + 𝑓(𝑥2 ) 𝑓(𝑥𝑛−1 ) + 𝑓(𝑥𝑛 )

= ℎ+ ℎ + ⋯+ ℎ
2 2 2
ℎ
= [𝑓(𝑥0 ) + 𝑓(𝑥1 ) + 𝑓(𝑥1 ) + 𝑓(𝑥2 ) + 𝑓(𝑥2 ) + ⋯ + 𝑓(𝑥𝑛−1 ) + 𝑓(𝑥𝑛−1 ) + 𝑓(𝑥𝑛 )]
2
ℎ
= [𝑓(𝑥0 ) + 2𝑓(𝑥1 ) + 2𝑓(𝑥2 ) + 2𝑓(𝑥3 ) + ⋯ + 2𝑓(𝑥𝑛−1 ) + 𝑓(𝑥𝑛 )]
2
𝑛−1
ℎ
= [𝑓(𝑥0 ) + 2 ∑ 𝑓(𝑥𝑘 ) + 𝑓(𝑥𝑛 )]
2
𝑘=1

En conclusión, si usted tiene la integral de la función 𝑓 definida desde 𝑥0 hasta 𝑥𝑛 se puede

aproximar como:

𝑥𝑛 𝑛−1
ℎ
∫ 𝑓(𝑥)𝑑𝑥 ≈ [𝑓(𝑥0 ) + 2 ∑ 𝑓(𝑥𝑘 ) + 𝑓(𝑥𝑛 )] (9)
𝑥0 2
𝑘=1
La expresión 9 se le conoce como regla trapezoidal o de ½

Ejemplo 4.8. Utilizando la regla trapezoidal (1/2) evaluar el valor erf(17)

𝑥
2 2
erf(𝑥) = ∫ 𝑒 −𝑢 𝑑𝑢
√𝜋
0

Evaluamos:
17 17
2 −𝑢2
2 2
erf(17) = ∫𝑒 𝑑𝑢 = (∫ 𝑒 −𝑢 𝑑𝑢)
√𝜋 √𝜋
0 0

2 17−0
Con 22 trapecios, 𝑓(𝑢) = 𝑒 −𝑢 , ℎ = 22

9
Métodos numéricos Competencia 4 Integración y derivación

x0 0
xn 17
n 22
h 0.77272727

puntos f(x)
x0 0 1 1
x1 0.77272727 0.55040141 1.100802812
x2 1.54545455 0.09177368 0.183547356
x3 2.31818182 0.00463571 0.00927141
x4 3.09090909 7.0937E-05 0.000141874
x5 3.86363636 3.2884E-07 6.57686E-07
x6 4.63636364 4.6181E-10 9.2362E-10
x7 5.40909091 1.9647E-13 3.92941E-13
x8 6.18181818 2.5322E-17 5.06431E-17
x9 6.95454545 9.8865E-22 1.9773E-21
x10 7.72727273 1.1694E-26 2.33875E-26
x11 8.5 4.1901E-32 8.38019E-32
x12 9.27272727 4.5483E-38 9.09668E-38
x13 10.0454545 1.4957E-44 2.99138E-44
x14 10.8181818 1.49E-51 2.98002E-51
x15 11.5909091 4.4967E-59 8.99343E-59
x16 12.3636364 4.1111E-67 8.22225E-67
x17 13.1363636 1.1386E-75 2.27728E-75
x18 13.9090909 9.5536E-85 1.91073E-84
x19 14.6818182 2.4284E-94 4.85671E-94
x20 15.4545455 1.87E-104 3.7398E-104
x21 16.2272727 4.362E-115 8.7238E-115
x22 17 3.082E-126 3.0824E-126
suma 2.293764111

Rfinal 1.000000133
Valor real 1

10
Métodos numéricos Competencia 4 Integración y derivación

Ejemplo 4.9. Resolver mediante la regla del trapecio la integral:

50
sin 𝑢
∫ 𝑑𝑢
𝑢
0.5

sin 𝑢
𝑛 = 30, 𝑓(𝑢) =
𝑢
Resultado con 800 pasos

x0 0.5
xn 50
n 800
h 0.061875
Algunas iteraciones
puntos f(x)
x0 0.5 0.958851077 0.95885108
x1 0.561875 0.948207104 1.89641421
x2 0.62375 0.936405791 1.87281158
x3 0.685625 0.923474056 1.84694811
x4 0.7475 0.909441351 1.8188827

Resultado

Vfinal 1.058567704
Vreal 1.058509654
Error 0.0055%

4.2.1 Reglas de Simpson (1/3) integración parabólica dado 3 puntos

Este método consiste en considerar que la curva es una función cuadrática, es decir, si
consideramos que la función la podemos escribir como un polinomio de Taylor de segundo
grado:

𝑓(𝑥) = 𝑓(𝑎) + 𝑓´(𝑎)(𝑥 − 𝑎) + 𝑓´´(𝑎)(𝑥 − 𝑎)2 + 𝑅𝑛+1 (10)

Podemos reescribir como:

𝑓(𝑥) ≈ 𝑎𝑥 2 + 𝑏𝑥 + 𝑐 (10)

Supongamos que tenemos el intervalo de integración [𝑥0 , 𝑥2 ] se divide en 2 partes, para lo

cual se consideran 3 puntos 𝑥𝑜 , 𝑥1 , 𝑥2 los cuales están a la misma distancia h, es decir:
𝑥2 − 𝑥0 𝑥2 + 𝑥0
ℎ= ⇒ 𝑥1 =
2 2
11
Métodos numéricos Competencia 4 Integración y derivación

Si sustituimos (10) en la integral definida de 𝑓(𝑥) desde 𝑥0 hasta 𝑥2 obtenemos:

𝑥2 𝑥2 𝑥
𝑎𝑥 3 𝑏𝑥 2 2

∫ 𝑓(𝑥) 𝑑𝑥 ≈ ∫ (𝑎𝑥 2 + 𝑏𝑥 + 𝑐) 𝑑𝑥 = + + 𝑐𝑥|

𝑥0 𝑥0 3 2 𝑥
(11)
0

Al evaluar los límites queda:

𝑥2
𝑎(𝑥23 − 𝑥03 ) 𝑏(𝑥22 − 𝑥02 )
∫ 𝑓(𝑥) 𝑑𝑥 ≈ + + 𝑐(𝑥2 − 𝑥0 )
𝑥0 3 2 (12)

Nótese que se pueden factorizar los binomios que aparecen en la integral interior, además
𝑥2 −𝑥0 𝑥2 +𝑥0
que si consideramos el hecho de sabemos que ℎ = 2
y que 𝑥1 = 2
obtenemos:

𝑥23 − 𝑥03 = (𝑥2 − 𝑥0 )(𝑥22 + 𝑥2 𝑥0 + 𝑥02 ) = 2ℎ(𝑥22 + 𝑥2 𝑥0 + 𝑥02 )

𝑥22 − 𝑥02 = (𝑥2 − 𝑥0 )(𝑥2 + 𝑥0 ) = 2ℎ(𝑥2 + 𝑥0 ) = (2ℎ)(2𝑥1 )

𝑥2 − 𝑥0 = 2ℎ
Reemplazando eso en (11) se obtiene:
𝑥2
2ℎ𝑎(𝑥22 + 𝑥2 𝑥0 + 𝑥02 ) 𝑏(2ℎ)(2𝑥1 )
∫ 𝑓(𝑥) 𝑑𝑥 ≈ + + 2𝑐ℎ
𝑥0 3 2

Se factoriza 2h
𝑥2
𝑎(𝑥22 + 𝑥2 𝑥0 + 𝑥02 ) 𝑏(2𝑥1 )
∫ 𝑓(𝑥) 𝑑𝑥 ≈ 2ℎ ( + + 𝑐)
𝑥0 3 2

Se busca común denominador

𝑥2
2𝑎(𝑥22 + 𝑥2 𝑥0 + 𝑥02 ) 3𝑏(2𝑥1 ) 6𝑐
∫ 𝑓(𝑥) 𝑑𝑥 ≈ 2ℎ ( + + )
𝑥0 6 6 6

Se factoriza 1/6:
𝑥2
2ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈ (2𝑎(𝑥22 + 𝑥2 𝑥0 + 𝑥02 ) + 6𝑏𝑥1 + 6𝑐)
𝑥0 6

Se simplifica la fracción y se separa 6𝑏𝑥1 :

𝑥2
ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈ (2𝑎(𝑥22 + 𝑥2 𝑥0 + 𝑥02 ) + 4𝑏𝑥1 + 2𝑏𝑥1 + 6𝑐)
𝑥0 3

Se desarrolla el producto 2𝑎(𝑥22 + 𝑥2 𝑥0 + 𝑥02 ) y se escribe 2𝑏𝑥1 = 𝑏(𝑥2 + 𝑥0 ):

𝑥2
ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 2𝑎𝑥22 + 2𝑎𝑥2 𝑥0 + 2𝑎𝑥02 + 4𝑏𝑥1 + 𝑏(𝑥2 + 𝑥0 ) + 6𝑐 )
𝑥0 3

Separamos los términos 2𝑎𝑥22 y 2𝑎𝑥02 y desarrollamos 𝑏(𝑥2 + 𝑥0 ) = 𝑏𝑥2 + 𝑏𝑥0 :

12
Métodos numéricos Competencia 4 Integración y derivación

𝑥2
ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑎𝑥22 + 𝑎𝑥22 + 2𝑎𝑥2 𝑥0 + 𝑎𝑥02 + 𝑎𝑥02 + 4𝑏𝑥1 + 𝑏𝑥2 + 𝑏𝑥0 + 6𝑐 )
𝑥0 3

Reagrupamos en un trinomio cuadrado perfecto

𝑥2
ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑎𝑥22 + 𝑎(𝑥2 + 𝑥0 )2 + 𝑎𝑥02 + 4𝑏𝑥1 + 𝑏𝑥2 + 𝑏𝑥0 + 6𝑐 )
𝑥0 3

Reemplazamos 2𝑥1 = 𝑥2 + 𝑥0 :
𝑥2
ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑎𝑥22 + 𝑎(2𝑥1 )2 + 𝑎𝑥02 + 4𝑏𝑥1 + 𝑏𝑥2 + 𝑏𝑥0 + 6𝑐 )
𝑥0 3

Elevamos al cuadrado (2𝑥1 )2 y separamos 6c=4c+c+c:

𝑥2
ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑎𝑥22 + 4𝑎𝑥12 + 𝑎𝑥02 + 4𝑏𝑥1 + 𝑏𝑥2 + 𝑏𝑥0 + 𝑐 + 4𝑐 + 𝑐 )
𝑥0 3

Reorganizamos:
𝑥2
ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑎𝑥22 + 𝑏𝑥2 + 𝑐 + 4(𝑎𝑥12 + 𝑏𝑥1 + 𝑐) + 𝑎𝑥02 + 𝑏𝑥0 + 𝑐 )
𝑥0 3

Finalmente se obtiene:
𝑥2
ℎ (13)
∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑓(𝑥2 ) + 4𝑓(𝑥1 ) + 𝑓(𝑥0 ) )
𝑥0 3

A la expresión 13 se le llama, método de 1/3 o integración cuadrática, y se utiliza cuando se

tiene un intervalo dividido en 2 o bien 3 puntos separados a la misma distancia.

13
Métodos numéricos Competencia 4 Integración y derivación

Método compuesto de (1/3): Consiste en dividir todo el intervalo de integración [𝑥0 , 𝑥2𝑛 ]
en 2n intervalos iguales. Nótese que debe ser una cantidad par de intervalos, porque el
método inicial con el que se hizo la prueba, fue dividido en 2. Ya hecho eso, sumamos todas
las áreas obtenidas.

Supóngase, que se tiene el intervalo [𝑥0 , 𝑥2𝑛 ], obtendríamos los siguientes intervalos:

[𝑥0 , 𝑥1 ], [𝑥1 , 𝑥2 ], [𝑥2 , 𝑥3 ], … , [𝑥2𝑛−1 , 𝑥2𝑛 ]

𝑥2𝑛 − 𝑥0
Cada uno de ellos, al ser iguales, tiene una longitud: ℎ=
2𝑛
Se escriben las integrales:

Intervalo Integral
𝑥2
ℎ
[𝑥0 , 𝑥1 ] ∪ [𝑥1 , 𝑥2 ] = [𝑥0 , 𝑥2 ] ∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑓(𝑥0 ) + 4𝑓(𝑥1 ) + 𝑓(𝑥2 ) )
𝑥0 3
𝑥4
ℎ
[𝑥2 , 𝑥3 ] ∪ [𝑥3 , 𝑥4 ] = [𝑥2 , 𝑥4 ] ∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑓(𝑥2 ) + 4𝑓(𝑥3 ) + 𝑓(𝑥4 ) )
𝑥2 3
𝑥6
ℎ
[𝑥4 , 𝑥5 ] ∪ [𝑥5 , 𝑥6 ] = [𝑥4 , 𝑥6 ] ∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑓(𝑥4 ) + 4𝑓(𝑥5 ) + 𝑓(𝑥6 ) )
𝑥4 3
⋮ ⋮
𝑥2𝑛
ℎ
[𝑥2𝑛−2 , 𝑥2𝑛−1 ] ∪ [𝑥2𝑛−1 , 𝑥2𝑛 ] = [𝑥2𝑛−2 , 𝑥2𝑛 ] ∫ 𝑓(𝑥) 𝑑𝑥 ≈ ( 𝑓(𝑥2𝑛−2 ) + 4𝑓(𝑥2𝑛−1 ) + 𝑓(𝑥2𝑛 ) )
𝑥2𝑛−2 3
La integral sobre intervalo [𝑥0 , 𝑥2𝑛 ] la obtenemos sumando todas las integrales anteriores.
Entonces realizamos la suma de ellas, posteriormente reemplazamos las integrales por el
método de 1/3 y simplificamos en notación sigma:
𝑥2𝑛 𝑥2 𝑥4 𝑥6 𝑥2𝑛
∫ 𝑓(𝑥) 𝑑𝑥 = ∫ 𝑓(𝑥) 𝑑𝑥 + ∫ 𝑓(𝑥) 𝑑𝑥 + ∫ 𝑓(𝑥) 𝑑𝑥 + ⋯ + ∫ 𝑓(𝑥) 𝑑𝑥
𝑥0 𝑥0 𝑥2 𝑥4 𝑥2𝑛−2

ℎ ℎ ℎ
≈ ( 𝑓(𝑥0 ) + 4𝑓(𝑥1 ) + 𝑓(𝑥2 ) ) + ( 𝑓(𝑥2 ) + 4𝑓(𝑥3 ) + 𝑓(𝑥4 ) ) + ⋯ + ( 𝑓(𝑥2𝑛−2 ) + 4𝑓(𝑥2𝑛−1 ) + 𝑓(𝑥2𝑛 ) )
3 3 3
ℎ
= [𝑓(𝑥0 ) + 4𝑓(𝑥1 ) + 𝑓(𝑥2 ) + 𝑓(𝑥2 ) + 4𝑓(𝑥3 ) + 𝑓(𝑥4 ) + ⋯ + 𝑓(𝑥2𝑛−2 ) + 4𝑓(𝑥2𝑛−1 ) + 𝑓(𝑥2𝑛 )]
3
ℎ
= [𝑓(𝑥0 ) + 4𝑓(𝑥1 ) + 2𝑓(𝑥2 ) + 4𝑓(𝑥3 ) + 2𝑓(𝑥4 ) + ⋯ + 2𝑓(𝑥2𝑛−2 ) + 4𝑓(𝑥2𝑛−1 ) + 𝑓(𝑥2𝑛 )]
3
𝑖𝑚𝑝𝑎𝑟𝑒𝑠 𝑝𝑎𝑟𝑒𝑠
2𝑛−1 2𝑛−2
ℎ
= 𝑓(𝑥0 ) + 4 ( ∑ 𝑓(𝑥𝑘 )) + 2 (∑ 𝑓(𝑥𝑗 )) + 𝑓(𝑥2𝑛 )
3
𝑘=1 𝑗
[ ]
Es decir:
𝑖𝑚𝑝𝑎𝑟𝑒𝑠 𝑝𝑎𝑟𝑒𝑠
2𝑛−1 2𝑛−2
𝑥2𝑛
ℎ
∫ 𝑓(𝑥) 𝑑𝑥 ≈
3
𝑓(𝑥0 ) + 4 ∑ 𝑓(𝑥𝑘 ) + 2 ∑ 𝑓(𝑥𝑗 ) + 𝑓(𝑥2𝑛 ) (14)
𝑥0 𝑘=1 𝑗
[ ( ) ]
( )
La expresión (14) se llama método de integración de (1/3) compuesto, se usa cuando
tenemos 2n+1 puntos o bien 2n intervalos.

14