CUADRATURA GAUSSIANA
El método de cuadratura de Gauss es un excelente método numérico para evaluar integrales
definidas de funciones, por medio de sumatorias simples y fáciles de implementar.
Este método es similar a la regla del trapecio con múltiples intervalos, pero aquí los intervalos no
necesariamente son iguales. Este criterio de integración numérica, lo que en realidad hace es
buscar dos o más puntos de la función, raíces de polinomios ortogonales o de Legendre, de modo
que se reduzca el error que conlleva métdos del trapecio que muchas veces son muy grandes.
Caracteristicas:
- Se puede extender a tres o más puntos
- Emplea valores no equidistantes de f(x)
- El método se puede aplicar a todo el intervalo de interés o puede dividirse en varios
subintervalos y aplicar el método de Gauss a cada uno de ellos.
Algoritmo general de método
𝑏 𝑏−𝑎 𝑏−𝑎 𝑎+𝑏
∫𝑎 𝑓(𝑥) 𝑑𝑥 = ∑𝑛𝑖=1 𝑤𝑖 ∙ 𝑓 ( 𝑧𝑖 + )
2 2 2
donde:
a = límite inferior de integración
b = límite superior de integración
n = número de puntos a emplear
𝑤𝑖 = valor del coeficiente de la cuadratura de Gauss
𝑧𝑖 = valor de las abscisas de la cuadratura de Gauss
Las raíces de los polinomios de Legendre 𝑃𝑛+1 (𝑧) y sus factores de ponderación se pueden
encontrar en forma tabular en la siguiente tabla:
Ejemplo 1)
Calcule el valor de la siguiente integral usando el método de la cuadratura de Gauss con 2 puntos
1.2 2
∫0.2 𝑒 𝑥 𝑑𝑥
Solución
n=2
𝑏 𝑏−𝑎 𝑏−𝑎 𝑎+𝑏
∫𝑎 𝑓(𝑥) 𝑑𝑥 = ∑2𝑖=1 𝑤𝑖 ∙ 𝑓 ( 𝑧𝑖 + )
2 2 2
Al desarrollar la sumatoria tenemos
𝑏 𝑏−𝑎 𝑏−𝑎 𝑎+𝑏 𝑏−𝑎 𝑎+𝑏
∫𝑎 𝑓(𝑥) 𝑑𝑥 = 2
[ 𝑤1 ∙ 𝑓 (
2
𝑧1 +
2
) + 𝑤2 ∙ 𝑓 ( 2
𝑧2 +
2
)]
Para la cuadratura Gasussiana de dos puntos los valores obtenidos de la tabla son:
𝑤1 = 𝑤2 = 1 → 𝑧1 = −0.5773502692 𝑧2 = 0.5773502692
Sustituyendo
1.2 𝑏−𝑎 𝑏−𝑎 𝑎+𝑏 𝑏−𝑎 𝑎+𝑏
∫0.2 = 2
[ 𝑤1 ∙ 𝑓 (
2
𝑧1 +
2
) + 𝑤2 ∙ 𝑓 ( 2
𝑧2 +
2
)]
1.2
2 1.2 − 0.2 1.2 − 0.2 0.2 + 1.2
∫ 𝑒 𝑥 𝑑𝑥 = [1 ∙𝑓 ( (−0.5773502692) + )+
2 2 2
0.2
1.2 − 0.2 0.2 + 1.2
1 ∙𝑓 ( (0.5773502692) + )]
2 2
1.2 2
∫0.2 𝑒 𝑥 𝑑𝑥 = 0.5 [1 ∙ 𝑓 (0.41132) + 1 ∙ 𝑓(0.98868)]
2
𝑓(0.41132) = 𝑒 0.41132 = 1.18434
2
𝑓(0.98868) = 𝑒 0.98868 = 2.65777
Por lo tanto:
1.2 2
∫0.2 𝑒 𝑥 𝑑𝑥 = 0.5 [ 1 ∙ 1.18434 + 1 ∙ 2.65777 ]
1.2 2
∫0.2 𝑒 𝑥 𝑑𝑥 = 1.92106
ECUACIONES DIFERENCIALES
Una ecuación diferencial es una ecuación que involucra a las derivadas de una función con la
propia función y/o las variables de las que depende. En sus aplicaciones, las funciones
generalmente representan cantidades y las derivadas son las tasas de variación de estas
cantidades.
Estas relaciones entre cantidades y sus tasas de cambio son frecuentes en áreas como la Física,
la Biología, la Ingeniería o la Economía. Existe una gran cantidad de fenómenos reales
complejos que pueden ser formulados en términos de ecuaciones diferenciales, lo cual permite
predecir sus comportamientos futuros, entre otras cosas. Las ecuaciones diferenciales también
pueden describir cómo cambian las poblaciones, cómo se transfiere el calor, cómo se desintegra el
material radiactivo, etc.
En la práctica la mayoría de las ecuaciones diferenciales no pueden resolverse mediante métodos
comunes y se recurre a los métodos numéricos, que permiten la utilización de equipo con alto
poder computable para resolver el problema.
METODO DE EULER
El método de Euler es un método numérico que sirve para resolver sistemas de ecuaciones
diferenciales ordinarias de primer orden.
Para poder aplicar este método, primero debe establecerse un intervalo solución y el tamaño del
espaciamiento ‘h’, y conocer las condiciones iniciales de cada ecuación diferencial del sistema.
El método consiste en usar la ecuación de Euler:
𝑦𝑛+1 = 𝑦𝑛 + ℎ 𝑓(𝑥𝑛 , 𝑦𝑛 )
Ejemplo 2)
Usar el método de Euler para aproximar la solución de la ecuación diferencial, dado en los puntos
𝑥 = 0.1, 0.2, 0.3, 0.4 𝑦 0.5, usando tamaño de paso h=0.1
𝑑𝑦 𝑥 𝑥
=− ó 𝑦′ = − , con 𝑦0 = 4
𝑑𝑥 𝑦 𝑦
Solución
𝑥0 = 0 𝑦0 = 4
0
𝑥1 = 0.1 𝑦1 = 𝑦0 + ℎ 𝑓(𝑥0 , 𝑦0 ) = 4 + 0.1 (− ) = 4
4
0.1
𝑥2 = 0.2 𝑦2 = 𝑦1 + ℎ 𝑓(𝑥1 , 𝑦1 ) = 4 + 0.1 (− ) = 3.9975
4
0.2
, 𝑥3 = 0.3 𝑦3 = 𝑦2 + ℎ 𝑓(𝑥2 , 𝑦2 ) = 3.9975 + 0.1 (− ) = 3.9926
3.9975
0.3
𝑥4 = 0.4 𝑦4 = 𝑦3 + ℎ 𝑓(𝑥3 , 𝑦3 ) = 3.9926 + 0.1 (− ) = 3.9850
3.9926
0.5
𝑥5 = 0.5 𝑦5 = 𝑦4 + ℎ 𝑓(𝑥4 , 𝑦4 ) = 3.2411 + 0.1 (− ) = 3.9749
3.2411
h= 0.1
i Xi yi
0 0 4.0000
1 0.1 4.0000
2 0.2 3.9975
3 0.3 3.9925
4 0.4 3.9850
5 0.5 3.9749
Ejemplo 3)
Usar el método de Euler para aproximar la solución de la ecuación diferencial,
𝑑𝑦
= 0.1 √𝑦 + 0.4 𝑥 2 ó 𝑦 ′ = 0.1 √𝑦 + 0.4 𝑥 2
𝑑𝑥
para n = 10 con y(2) = 4
Calcule una aproximación para y(2.5)
Solución
𝑏−𝑎 2.5−2
ℎ= = = 0.05
𝑛 10
𝑥0 = 2 𝑦0 = 4
𝑥1 = 2.05 𝑦1 = 𝑦0 + ℎ 𝑓(𝑥0 , 𝑦0 ) = 4 + 0.05 (0.1 √4 + 0.4 ∗ 22 ) = 4.09
𝑥2 = 0.2 𝑦2 = 𝑦1 + ℎ 𝑓(𝑥1 , 𝑦1 ) = 4.09 + 0.05 (0.1 √4.09 + 0.4 ∗ 2.052 ) = 4.1842
h= 0.05
i Xi yi
0 2 4.0000
1 2.05 4.0900
2 2.1 4.1842
3 2.15 4.2826
4 2.2 4.3854
5 2.25 4.4927
6 2.3 4.6045
7 2.35 4.7210
8 2.4 4.8423
9 2.45 4.9686
10 2.5 5.0997
METODO DE RUNGE-KUTTA
Es uno de los métodos más empleados para la solución de ecuaciones diferenciales y obtienen
una mejor aproximación. Las ecuaciones son una derivación de la serie de Taylor y su ventaja es
su fácil aplicación. Existen 4 o 5 grandes grupos de métodos de solución Runge-Kutta
dependiendo del orden de precisión que se requiera.
RUNGE-KUTTA 2° ORDEN
Solo toma 2 términos de la serie de Taylor, consiste en determinar dos pendientes 𝑘1 al inicio y 𝑘2
al final, en el intervalo 𝑥𝑖 𝑦 𝑥𝑖+1
de manera que consiste en obtener una ecuación considerando la siguiente forma genérica de
aproximación el valor de 𝑦𝑖+1
1
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 𝑘2 ) ∆𝑥
2
ℎ
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 𝑘2 )
2
donde
𝑘1 = 𝑓(𝑥𝑖 , 𝑦𝑖 ) pendiente al inicio del intervalo
𝑘2 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + ∆𝑥) pendiente al final del intervalo
RUNGE-KUTTA 3° ORDEN
1
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 4𝑘2 + 𝑘3 ) ∆𝑥
6
ℎ
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 4𝑘2 + 𝑘3 )
6
donde
𝑘1 = 𝑓(𝑥𝑖 , 𝑦𝑖 ) pendiente al inicio del intervalo
1 1
𝑘2 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + 𝑘1 ∆𝑥) pendiente en punto medio del intervalo
2 2
𝑘3 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 − 𝑘1 ∆𝑥 + 2𝑘2 ∆𝑥 ) pendiente al final del intervalo
RUNGE-KUTTA 4° ORDEN
1
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 2𝑘2 + 2𝑘3 + 𝑘4 ) ∆𝑥
6
ℎ
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 2𝑘2 + 2𝑘3 + 𝑘4 )
6
donde
𝑘1 = 𝑓(𝑥𝑖 , 𝑦𝑖 ) pendiente al inicio del intervalo
1 1
𝑘2 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + 𝑘1 ∆𝑥) pendiente en punto medio del intervalo
2 2
1 1
𝑘3 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + 𝑘2 ∆𝑥 ) pendiente en punto medio del intervalo
2 2
𝑘4 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + 𝑘3 ∆𝑥 ) pendiente al final del intervalo
Mientras mayor sea el orden más preciso será el resultado.
Ejemplo 4)
Resolver EDO aplicando Runge-Kutta de 3° orden
𝑑𝑦
= 4𝑒 0.8𝑥 − 0.5 𝑦 en el intervalo [0,4]
𝑑𝑥
Para las condiciones iniciales 𝑥0 = 0 𝑦 𝑦0 = 2 con ∆𝑥 = 1
Solución:
Intervalo 1
𝑥𝑖 = 𝑥0 = 0 𝑥𝑖+1 = 1
𝑦𝑖 = 𝑦0 = 2 𝑦𝑖+1 = ?
Runge-Kutta 3° orden
ℎ
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 4𝑘2 + 𝑘3 )
6
donde
𝑘1 = 𝑓(𝑥𝑖 , 𝑦𝑖 ) pendiente al inicio del intervalo
1 1
𝑘2 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + 𝑘1 ∆𝑥) pendiente en punto medio del intervalo
2 2
𝑘3 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 − 𝑘1 ∆𝑥 + 2𝑘2 ∆𝑥 ) pendiente al final del intervalo
Cálculo:
𝑘1 = 𝑓(𝑥𝑖 , 𝑦𝑖 ) = 𝑓(0,2)
𝑘1 = 4𝑒 0.8 (0) − 0.5 (2) = 3
1 1
𝑘2 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + 𝑘1 ∆𝑥)
2 2
1 1
= 𝑓 (0 + (1), 2 + (3) 1 )
2 2
= 𝑓 (0.5,3.5)
𝑘2 = 4𝑒 0.8 (0.5) − 0.5 (3.5) = 4.217298
𝑘3 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 − 𝑘1 ∆𝑥 + 2𝑘2 ∆𝑥 )
= 𝑓(0 + 1 , 2 − 3(1) + 2 (4.217298 ∗ 1)
= 𝑓 ( 1, 7.434596)
𝑘3 = 4𝑒 0.8 (1) − 0.5 (7.434596) = 5.1848657
Por lo tanto
ℎ
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 4𝑘2 + 𝑘3 )
6
1
𝑦𝑖+1 = 2 + 6 (3 + 4(4.217298) + 5.1848657)
𝑦𝑖+1 = 6.175676
Intervalo 2
𝑥𝑖 = 1 𝑥𝑖+1 = 2
𝑦𝑖 = 6.175676 𝑦𝑖+1 = ?
𝑘1 = 𝑓(𝑥𝑖 , 𝑦𝑖 ) = 𝑓(1, 6.175676)
𝑘1 = 4𝑒 0.8 (1) − 0.5 (6.175676) = 5.814325
1 1
𝑘2 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + 𝑘1 ∆𝑥)
2 2
1
= 𝑓 (1 + (1),
2
6.175676 + 12 (5.814325) 1 )
= 𝑓 (1.5, 9.08283)
𝑘2 = 4𝑒 0.8 (1.5) − 0.5 ( 9.08283) = 8.739048
𝑘3 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 − 𝑘1 ∆𝑥 + 2𝑘2 ∆𝑥 )
= 𝑓(1 + 1 , 6.175676 − 5.814325 (1) + 2 (8.739048 ∗ 1)
= 𝑓 ( 2, 17.839447)
𝑘3 = 4𝑒 0.8 (2) − 0.5 (17.839447) = 10.8924062
Por lo tanto
ℎ
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 4𝑘2 + 𝑘3 )
6
1
𝑦𝑖+1 = 6.175676 + 6 (5.814325 + 4(8.739048) + 10.8924062)
𝑦𝑖+1 = 14.786162
Intervalo 3
𝑥𝑖 = 2 𝑥𝑖+1 = 3
𝑦𝑖 = 14.786162 𝑦𝑖+1 = ?
𝑘1 = 𝑓(𝑥𝑖 , 𝑦𝑖 ) = 𝑓(2, 14.786162)
𝑘1 = 4𝑒 0.8 (2) − 0.5 (14.786162) = 12.4190481
1 1
𝑘2 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 + 𝑘1 ∆𝑥)
2 2
1
= 𝑓 (2 + (1),
2
14.786162 + 12 (12.4190481) 1 )
= 𝑓 (2.5, 20.995687)
𝑘2 = 4𝑒 0.8 (2.5) − 0.5 (20.995687 ) = 19.05838077
𝑘3 = 𝑓(𝑥𝑖 + ∆𝑥 , 𝑦𝑖 − 𝑘1 ∆𝑥 + 2𝑘2 ∆𝑥 )
= 𝑓(2 + 1 , 14.786162 − 12.4190481 (1) + 2 (19.05838077 ∗ 1)
= 𝑓 ( 3, 40.48387664
𝑘3 = 4𝑒 0.8 (3) − 0.5 (40.48387664) = 23.850672
Por lo tanto
ℎ
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 4𝑘2 + 𝑘3 )
6
1
𝑦𝑖+1 = 14.786162 + 6 (12.4190481 + 4(19.05838077) + 23.850672)
𝑦𝑖+1 = 33.5367196
Intervalo 4
𝑥𝑖 = 3 𝑥𝑖+1 = 4
𝑦𝑖 = 33.5367196 𝑦𝑖+1 = ?
𝑘1 = 𝑓(3, 33.5367196 ) = 27.3243457
𝑘2 = 𝑓(3.5, 47.19889186) = 42.17914115
𝑘3 = 𝑓(4, 90.5706562) = 52.84479269
Por lo tanto
ℎ
𝑦𝑖+1 = 𝑦𝑖 + (𝑘1 + 4𝑘2 + 𝑘3 )
6
1
𝑦𝑖+1 = 33.5367196 + 6 (27.3243457 + 4(42.17914115) + 52.84479269)
𝑦𝑖+1 = 75.0176
Xi yi
0 2.0000
1 6.1756
2 14.7861
3 33.5367
4 75.0176