NotasTema5_0 - Axel

Vista previa del material en texto

Tema 5: Solución numérica de ecuaciones y sistemas de ecuaciones diferenciales
M. I. Luis Ángel Santamaŕıa Padilla
Facultad de Ingenieŕıa, UNAM
1 Solución numérica de ecuaciones diferenciales de
primer orden
Método de Euler
Ejemplo 1
Método de Taylor de orden superior
Ejemplo 2
Métodos de Runge-Kutta de segundo orden
Ejemplo 3
Métodos de Runge-Kutta de tercer orden (RK3)
Ejemplo 4
Métodos de Runge-Kutta de cuarto orden (RK4)
Ejemplo 5
2 Solución numérica de sistemas de ecuaciones de
primer orden
Ecuaciones diferenciales en variables de estado
Ejemplo 6
Ejemplo 7
Método de Euler para sistemas
Ejemplo 8
Métodos numéricos de Heun y RK4 clásico para
sistemas de EDOs
Ejemplo 9
3 Solución numérica de problemas con valores en la
frontera (PVF)
Método de disparo
Ejemplo 10
Método de diferencias finitas
Ejemplo 11
Objetivo tema 5
El estudiante comparará algunos métodos de aproximación para la solución de ecuaciones y sistemas de
ecuaciones diferenciales, sujetas a condiciones iniciales o de frontera
Se comenzará con métodos para una ecuación de primer orden de una sola variable con problema de valor inicial
de la forma
y′ = f(x, y), y(a) = y0, x0 = a ≤ x ≤ b = xn
donde y0 es la condición inicial, x es la variable independiente que toma valores de [a, b] y se asume que existe
solución única en el intervalo [a, b].
El intervalo se divide en n segmentos de igual tamaño, tal que
x1 = x0 + h, x2 = x0 + 2h, . . . , xn = x0 + nh
La solución en x0 está disponible por la condición inicial. El objetivo es encontrar estimados de la solución en los
puntos subsecuentes x1, x2, . . . , xn
Solución numérica de ecuaciones diferenciales de primer orden
Métodos de un paso
Estiman la solución yi+1 en la ubicación xi+1 extrapolando la solución estimada yi en la ubicación xi.
Dependiendo del método, será la manera en que se calcule el nuevo punto estimado.
La figura describe un método simple de un paso
La pendiente ϕ se usa para extrapolar, a partir de
yi al nuevo estimado yi+1
yi+1 = yi + hϕ i = 0, 1, . . . , n− 1
Cada método utiliza una forma particular para
calcular la pendiente ϕ
Solución numérica de ecuaciones diferenciales de primer orden Método de Euler
Método de Euler
La expansión en series de Taylor de y(x1) desde x0, está dada
por
y(x1) = y(x0 + h) = y(x0) + hy
′(x0) +
1
2!
h2y′′(x0) + . . .
Se conservan los términos lineales
y(x1) = y(x0) + hy
′(x0) +
1
2!
h2y′′(ξi)
Para algún ξi entre xi y xi+1.
Se introduce la notación yi = y(xi), yi+1 = y(xi+1).
Se sabe que y′(xi) = f(xi, yi)
La solución estimada se puede encontrar
como
yi+1 = yi + hf(xi, yi)
i = 0, 1, 2, . . . , n− 1
Este método es conocido como
método de Euler
Comparando con yi+1 = yi + hϕ, se
observa que la pendiente en xi se
estima por f(xi, yi), el cual es la
primera derivada en xi (y
′(xi))
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 1
Ejemplo 1
Considere el problema de valor inicial
2y′ + y = e−x, y(0) =
1
2
, 0 ≤ x ≤ 1
Encuentre y(0.1) e y(0.2), utilizando el método de Euler con un paso de h = 0.1.
Solución.
2y′ = e−x − y
y′ =
1
2
(
e−x − y
)
︸ ︷︷ ︸
f(x,y)
Se busca y(0.1), y(0.2)
Recordar
yi+1 = yi + hf(xi, yi)
Para calcular y(0.1) se tiene que y(0.1) = y1
y1 = y0 + hf(x0, y0) =
1
2
+ 0.1f(0, 1/2)
y1 = 0.5250
Para calcular y(0.2) se tiene que y(0.2) = y2
y2 = y1 + hf(x1, y1) = 0.5250 + 0.1f(0.1, 0.5250)
y2 = 0.5440
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 1
xi Euler Real % Error
0 0.5000 0.5000 0
0.1 0.5250 0.5220 0.5734
0.2 0.5440 0.5385 1.0151
0.3 0.5577 0.5502 1.3603
0.4 0.5669 0.5578 1.6328
0.5 0.5721 0.5617 1.8489
0.6 0.5738 0.5624 2.0204
0.7 0.5725 0.5604 2.1561
0.8 0.5687 0.5562 2.2624
0.9 0.5628 0.5499 2.3443
1.0 0.5550 0.5419 2.4057
0 0.2 0.4 0.6 0.8 1
0.5
0.51
0.52
0.53
0.54
0.55
0.56
0.57
0.58
Solución numérica de ecuaciones diferenciales de primer orden Método de Taylor de orden superior
Métodos de Taylor de orden superior
Se retienen más términos en la serie de Taylor
y(xi+1) = y(xi) + hy
′(xi) +
1
2!
h2y′′(xi) + . . .+
1
k!
hky(k)(xi) +
1
(k + 1)!
hk+1y(k+1)(ξi)
donde ξi se encuentra entre xi y xi+1
El método de Taylor de orden k está dado por
yi+1 = yi + hpk(xi, yi) i = 0, 1, 2, . . . , n− 1
donde
pk(xi, yi) = f(xi, yi) +
1
2!
hf ′(xi, yi) + . . .+
1
k!
hk−1f (k−1)(xi, yi)
El método de Euler es un método de Taylor de primer orden
Obervaciones:
Entre mayor sea el orden del
método de Taylor, más exactas
serán los estimados
La reducción del error requiere del
cálculo de las derivadas de f(x, y)
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 2
Ejemplo 2
Considere el problema de valor inicial
2y′ + y = e−x, y(0) =
1
2
, 0 ≤ x ≤ 1
Encuentre y(0.1), y(0.2), utilizando el método de Taylor de segundo orden con un paso de h = 0.1.
Solución.
f(x, y) =
1
2
(e−x − y)
⇓
f ′(x, y) =
1
2
−e−x − y′︸︷︷︸
y′=f(x,y)

=
1
2
(
−e−x − 1
2
(e−x − y)
)
f ′(x, y) =
1
4
(
−3e−x + y
)
Para el método de Taylor de segundo orden
yi+1 = yi + hp2(xi, yi), i = 0, 1, 2, . . . , n− 1
donde
p2(xi, yi) = f(xi, yi) +
1
2!
hf ′(xi, yi)
Entonces
yi+1 = yi + h
(
f(xi, yi) +
1
2!
hf ′(xi, yi)
)
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 2
Ejemplo 2
Para calcular y(0.1) se tiene que y(0.1) = y1
y1 = y0 + h
(
f(x0, y0) +
1
2
hf ′(x0, y0)
)
y1 = 0.5 + 0.1
(
f(0, 0.5) +
1
2
(0.1)f ′(0, 0.5)
)
= 0.5220
Para calcular y(0.2) se tiene que y(0.2) = y2
y2 = y1 + h
(
f(x1, y1) +
1
2
hf ′(x1, y1)
)
y2 = 0.5220 + 0.1
(
f(0.1, 0.5220) +
1
2
(0.1)f ′(0.1, 0.5220)
)
y2 = 0.5385
0 0.2 0.4 0.6 0.8 1
0.5
0.51
0.52
0.53
0.54
0.55
0.56
0.57
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de segundo orden
Métodos de Runge-Kutta
Estos métodos generan soluciones estimadas con la misma exactitud de los métodos de Taylor,
pero sin tener que calcular derivadas.
Recordar que
y′ = f(x, y), y(a) = ya, a = x0 ≤ x ≤ xn = b
Se expresa como
yi+1 = y1 + hϕ(xi, yi)
donde ϕ(xi, yi):
Es una función de incremento
Es una pendiente conveniente en el intervalo [xi, xi+1] que se utiliza para extrapolar yi+1 desde yi
El número de puntos utilizados en [xi, xi+1] define el orden del método de Runge-Kutta
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de segundo orden
Métodos de Runge-Kutta de segundo orden (RK2)
La función de incremento se
expresa como
ϕ(xi, yi) = a1k1 + a2k2 tal que
yi+1 = yi + h(a1k1 + a2k2)
con
k1 = f(xi, yi)
k2 = f(xi + b1h, yi + c11k1h)
donde a1, a2, b1 y c11 son
constantes
Estas constantes se calculan con los primeros tres términos de la
serie de Taylor
yi+1 = yi + hy
′(xi) +
1
2!
h2y′′(xi) +O(h
3)
el término y′(xi) es f(xi, yi), mientras
y′′(xi) = f
′(xi, yi) =
∂f
∂x
∣∣∣∣
(xi,yi)
+
∂f
∂y
∣∣∣∣
(xi,yi)
y′
=
∂f
∂x
∣∣∣∣
(xi,yi)
+
∂f
∂y
∣∣∣∣
(xi,yi)
f(xi, yi)
Sustituyendo y′ e y′′
yi+1 = yi + hf(xi, yi) +
1
2
h2
∂f
∂x
∣∣∣∣
(xi,yi)
+
1
2
h2
∂f
∂y
∣∣∣∣
(xi,yi)
f(xi, yi)
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de segundo orden
Métodos de Runge-Kutta de segundo orden (RK2)
Se calculará yi+1 utilizando otro enfoque.
Se sabe que k2 = f(xi + b1h, yi + c11k1h) es una función de dos variables y se puede expandir alrededor de
(xi, yi) como
k2 = f(xi + b1h, yi + c11k1h) = f(xi, yi) + b1h
∂f
∂x
∣∣∣∣
(xi,yi)
+ c11h
∂f
∂y
∣∣∣∣
(xi,yi)
Sustituyendo en yi+1 = yi + h(a1k1 + a2k2) se encuentra
yi+1 = yi + h
(
a1f(xi, yi) + a2
(
f(xi, yi) + b1h
∂f
∂x
∣∣∣∣
(xi,yi)
+ c11h
∂f
∂y
∣∣∣∣
(xi,yi)
))
yi+1 = yi + (a1 + a2)hf(xi, yi) + a2b1h
2 ∂f
∂x
∣∣∣∣
(xi,yi)
+ a2c11h
2 ∂f
∂y
∣∣∣∣
(xi,yi)
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de segundo orden
Métodos de Runge-Kutta de segundo orden (RK2)Retomando las ecuaciones
yi+1 = yi + hf(xi, yi) +
1
2
h2
∂f
∂x
∣∣∣∣
(xi,yi)
+
1
2
h2
∂f
∂y
∣∣∣∣
(xi,yi)
f(xi, yi)
yi+1 = yi + (a1 + a2)hf(xi, yi) + a2b1h
2 ∂f
∂x
∣∣∣∣
(xi,yi)
+ a2c11h
2 ∂f
∂y
∣∣∣∣
(xi,yi)
Se igualan las ecuaciones y se observa
a1 + a2 = 1; a2b1 =
1
2
; a2c11 =
1
2
Observaciones
Se tienen cuatro incógnitas y tres ecuaciones
No existe un conjunto únido de solución
Si se asigna un valor a alguna de las constantes, las restantes se pueden calcular
Existen varias versiones de los métodos RK2
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de segundo orden
Métodos de Runge-Kutta de segundo orden (RK2)
Recordar
a1 + a2 = 1; a2b1 =
1
2
; a2c11 =
1
2
yi+1 = yi + h(a1k1 + a2k2)
k1 = f(xi, yi)
k2 = f(xi + b1h, yi + c11k1h)
Método mejorado de Euler
Se supone a2 = 1
a1 = 0; b1 =
1
2
; c11 =
1
2
Se obtiene
yi+1 = y1 + hk2
k1 = f(xi, yi)
k2 = f
(
xi +
1
2
h, yi +
1
2
k1h
)
Método de Heun
Se supone a2 =
1
2
a1 =
1
2
; b1 = 1; c11 = 1
Se obtiene
yi+1 = y1 +
1
2
h(k1 + k2)
k1 = f(xi, yi)
k2 = f (xi + h, yi + k1h)
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de segundo orden
Métodos de Runge-Kutta de segundo orden (RK2)
Método de Ralston
Se supone a2 =
2
3
a1 =
1
3
; b1 =
3
4
; c11 =
3
4
Se obtiene
yi+1 = y1 +
1
3
h(k1 + 2k2)
k1 = f(xi, yi)
k2 = f
(
xi +
3
4
h, yi +
3
4
k1h
)
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 3
Ejemplo 3
Considere el PVI
y′ − x2y = 2x2; y(0) = 1, 0 ≤ x ≤ 1
Calcule el valor estimado de y1 = y(0.1)
utilizando los métodos RK2 y un paso h = 0.1.
Solución.
La solución anaĺıtica es y = 3ex
3/3 − 2
y′ = 2x2 + x2y = x2(2 + y)
f(x, y) = x2(2 + y)
Euler mejorado
k1 = f(x0, y0) = f(0, 1) = 0
k2 = f
(
x0 +
1
2
h, y0 +
1
2
k1h
)
= f(0.05, 1)
k2 = 0.0075
y1 = y0 + hk2 = 1 + (0.1)(0.0075)
y1 = 1.0008
Heun
k1 = f(x0, y0) = f(0, 1) = 0
k2 = f (x0 + h, y0 + k1h) = f(0.1, 1)
k2 = 0.0300
y1 = y0 +
1
2
h(k1 + k2) = 1 + (0.5)(0.1)(0.0300)
y1 = 1.0015
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 3
Ejemplo 3
Ralston
k1 = f(x0, y0) = f(0, 1) = 0
k2 = f
(
x0 +
3
4
h, y0 +
3
4
k1h
)
= f(0.075, 1)
k2 = 0.0169
y1 = y0 +
1
3
h(k1 + 2k2) = 1 +
1
3
(0.1)(2(0.0169))
y1 = 1.0011
0 0.2 0.4 0.6 0.8 1
1
1.2
1.4
1.6
1.8
2
2.2
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de tercer orden (RK3)
Métodos de Runge-Kutta de tercer orden (RK3)
La función de incremento
ϕ(xi, yi) = a1k1 + a2k2 + a3k3, tal que
yi+1 = yi + h(a1k1 + a2k2 + a3k3)
donde
k1 = f(xi, yi)
k2 = f(xi + b1h, yi + c11k1h)
k3 = f(xi + b2h, yi + c21k1h+ c22k2h)
donde a1, a2, a3, b1, b2, c11, c21 y c22 son constantes.
Cada conjunto de estas variables determina un método
RK3
Observaciones
Las constantes se encuentran utilizando una
expansión en series de Taylor utilizando los
primeros cuatro términos
Se tendrán seis ecuaciones y ocho incógnitas
Se asignan valores a dos variables y se calculan las
seis restantes
Se presentarán dos métodos RK3
Método RK3 clásico
Método de Heun
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de tercer orden (RK3)
Métodos de Runge-Kutta de tercer orden (RK3)
Método RK3 clásico
yi+1 = yi +
1
6
h (k1 + 4k2 + k3)
donde
k1 = f (xi, yi)
k2 = f
(
xi +
1
2
h, yi +
1
2
k1h
)
k3 = f (xi + h, yi − k1h+ 2k2h)
Método de Heun
yi+1 = yi +
1
4
h (k1 + 3k3)
donde
k1 = f (xi, yi)
k2 = f
(
xi +
1
3
h, yi + k1h
)
k3 = f
(
xi +
2
3
h, yi +
2
3
k2h
)
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 4
Ejemplo 4
Considere el PVI
y′ − x2y = 2x2, y(0) = 1, 0 ≤ x ≤ 1
Calcule el valor estimado de y1 = y(0.1) con los
métodos RK3 y un paso h = 0.1.
Solución.
y′ = 2x2 + x2y
y′ = x2(2 + y)
f(x, y) = x2(2 + y)
RK3 clásico
k1 = f (x0, y0)
k1 = f(0, 1) = 0
k2 = f
(
x0 +
1
2
h, y0 +
1
2
k1h
)
k2 = f (0.05, 1) = 0.0075
k3 = f (x0 + h, y0 − k1h+ 2k2h)
k3 = f (0.1, 1.0015) = 0.0300
y1 = y0 +
1
6
h (k1 + 4k2 + k3)
y1 = 1 +
1
6
(0.1) ((0) + 4(0.0075) + (0.0300))
= 1.0010
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 4
Ejemplo 4
Heun RK3
k1 = f (x0, y0)
k1 = f(0, 1) = 0
k2 = f
(
x0 +
1
3
h, y0 + k1h
)
k2 = f (0.0333, 1) = 0.0033
k3 = f
(
x0 +
2
3
h, y0 +
2
3
k2h
)
k3 = f (0.0667, 1.0002) = 0.0133
y1 = y0 +
1
4
h (k1 + 3k3)
y1 = 1 +
1
4
(0.1) ((0) + 3(0.0133))
= 1.0010
0 0.2 0.4 0.6 0.8 1
1
1.2
1.4
1.6
1.8
2
2.2
Solución numérica de ecuaciones diferenciales de primer orden Métodos de Runge-Kutta de cuarto orden (RK4)
Métodos de Runge-Kutta de cuarto orden (RK4)
La función de incremento
ϕ(xi, yi) = a1k1 + a2k2 + a3k3 + a4k4, tal que
yi+1 = yi + h(a1k1 + a2k2 + a3k3 + a4k4)
donde
k1 = f(xi, yi)
k2 = f(xi + b1h, yi + c11k1h)
k3 = f(xi + b2h, yi + c21k1h+ c22k2h)
k4 = f(xi + b3h, yi + c31k1h+ c32k2h+ c33k3h)
Método RK4 clásico
yi+1 = yi +
1
6
h (k1 + 2k2 + 2k3 + k4)
donde
k1 = f (xi, yi)
k2 = f
(
xi +
1
2
h, yi +
1
2
k1h
)
k3 = f
(
xi +
1
2
h, yi +
1
2
k2h
)
k4 = f (xi + h, yi + k3h)
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 5
Ejemplo 5
Considere el PVI
y′ − x2y = 2x2, y(0) = 1, 0 ≤ x ≤ 1
Calcule el valor estimado de y1 = y(0.1) con
RK4 clásico y un paso h = 0.1.
Solución.
y′ = 2x2 + x2y
y′ = x2(2 + y)
f(x, y) = x2(2 + y)
RK4 clásico
k1 = f (x0, y0)
k1 = f(0, 1) = 0
k2 = f
(
x0 +
1
2
h, y0 +
1
2
k1h
)
k2 = f (0.05, 1) = 0.0075
k3 = f
(
x0 +
1
2
h, y0 +
1
2
k2h
)
k3 = f (0.05, 1.0004) = 0.0075
k4 = f (x0 + h, y0 + k3h)
k4 = f (0.1, 1.0008) = 0.0300
y1 = y0 +
1
6
h (k1 + 2k2 + 2k3 + k4)
y1 = 1 +
1
6
h ((0) + 2(0.0075) + 2(0.0075) + (0.0300))
y1 = 1.0010
Solución numérica de ecuaciones diferenciales de primer orden Ejemplo 5
0 0.5 1 1.5 2
0
5
10
15
20
25
30
35
40
45
Solución numérica de sistemas de ecuaciones de primer orden
Sistemas de ecuaciones diferenciales
Observaciones
Se transformarán modelos en sistemas de ecuaciones diferenciales de primer orden
Se resolverán los sistemas de ecuaciones diferenciales de primer orden numéricamente.
title
Transformación en un sistema de EDO de primer orden
Variables de estado. Forman el conjunto ḿınimo de variables linealmente independientes que describen
totalmente el estado del sistema.
Se determinan de la siguiente manera:
Se elige el número de variables de estado como tantas condiciones iniciales se dispongan
Son aquellas para las cuales se tienen condiciones iniciales
Las variables de estado se representarán por ui. Si fueran tres variables se denotaŕıan por u1, u2, u3
Solución numérica de sistemas de ecuaciones de primer orden Ecuaciones diferenciales en variables de estado
Ecuaciones diferenciales en variables de estado
Si un sistema tiene m variables de estado. u1, u2, . . . , um
Hay m ecuaciones en variables de estado
Cada una de ellas es una EDO de primer orden de la forma
u̇i = f(t, u1, u2, . . . , um) i = 1, 2, . . . ,m
El lado izquierdo es la primer derivada respecto al tiempo t de la variable de estado ui
El lado derecho, es una función algebraica de las variables de estado y posiblemente de t
El sistema de EDO de primer orden se puede expresar en forma matricial como
u̇ = f(t,u); u(t) =

u1
u2
...
um
 ; f(t,u) =

f1(t,u)
f2(t,u)
...
fm(t,u)

donde u es el vector de estados
Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 6
Ejemplo 6
Considere el PVI descrito por
3y′′′ − y′′ − 5y′ − 3y = e−x/2
y(0) = 0; y′(0) = −1; y′′(0) = 1
Obtenga una representación en variables de estado.
Solución.
Se tienen tres condiciones inicales
Se tendrán tres variables de estado u1, u2, u3
Las variables de estado son aquellas variables para
las cuales se tienen condiciones iniciales. Entonces
u1 =y, u2 = y
′, u3 = y
′′
Hay tres ecuaciones de estado, que se forman como
u′1 = y
′ = u2
u′2 = y
′′ = u3
u′3 = y
′′′ =
1
3
[
y′′ + 5y′ + 3y + e−x/2
]
u′3 = y
′′′ =
1
3
[
u3 + 5u2 + 3u1 + e
−x/2
]
Sujeto a
u1(0) = y(0) = 0
u2(0) = y
′(0) = −1
u3(0) = y
′′(0) = 1
Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 6
Ejemplo 6
En forma vectorial
u̇ = f(t,u)
u =
u1u2
u3
 , f =
 u2u3
1
3
[
u3 + 5u2 + 3u1 + e
−x/2
]
 , u0 =
 0−1
1

Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 7
Ejemplo 7
Sea el siguiente modelo matemático de un sistema
dinámico
ẍ1 + 2ẋ1 + 2(x1 − x2) = e−t sin(t)
ẋ2 − 2(x1 − x2) = 0
Sujeto a x1(0) = 0, x2(0) = 0, ẋ1(0) = 1. Obtenga un
representación en variables de estado.
Solución.
Se tienen tres condiciones iniciales
Habrá tres variables de estado u1, u2, u3
La forma de seleccionar las variables de estado es la
siguiente: primero se asignan a las variables sin
derivar y después a las que tienen derivadas.
Entonces
u1 = x1, u2 = x2, u3 = ẋ1
Hay tres ecuaciones de estado, que se forman como
u̇1 = ẋ1 = u3
u̇2 = ẋ2 = 2(x1 − x2)
u̇2 = 2(u1 − u2)
u̇3 = ẍ1 = −2ẋ1 − 2(x1 − x2) + e−t sin(t)
u̇3 = −2u3 − 2(u1 − u2) + e−t sin(t)
Sujeto a
u1(0) = x1(0) = 0
u2(0) = x2(0) = 0
u3(0) = ẋ1(0) = 1
Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 7
Ejemplo 7
En forma vectorial
u̇ = f(t,u)
u =
u1u2
u3
 , f =
 u32(u1 − u2)
−2u3 − 2(u1 − u2) + e−t sin(t)
 , u0 =
00
1

Solución numérica de sistemas de ecuaciones de primer orden Método de Euler para sistemas
Solución numérica de sistemas de ecuaciones diferenciales de primer orden
Los métodos de solución de
EDO de primer orden se
pueden extender a la
solución de sistemas de EDO
de primer orden
Se presentarán los métodos
de
Euler
Heun
RK4 clásico
Método de Euler para sistemas
El intervalo [a, b] se divide en subintervalos de la misma longitud h, tal que
x1 = x0 + h, x2 = x0 + 2h, , . . . , xn = x0 + nh
El método de Euler para sistemas es de la forma
ui+1 = ui + hf(xi,ui); i = 0, 1, 2, . . . , n− 1
Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 8
Ejemplo 8
Considere la EDO 3y′′′ − y′′ − 5y′ − 3y = e−x/2, sujeta a y(0) = 0, y′(0) = 1, 0 ≤ x ≤ 1.
Utilice el método de Euler para sistemas, con paso h = 0.1 y encuentre un estimado para y2 = y(0.2).
Solución.
Del ejemplo anterior se tiene
u̇ = f(t,u)
u =
u1u2
u3
 , f =
 u2u3
1
3
[
u3 + 5u2 + 3u1 + e
−x/2
]
 , u0 =
 0−1
1

Para encontrar y(0.2) se necesita encontrar el vector solución u2 = u(0.2) y tomar su primer componente,
que corresponde con y(0.2)
f(x0,u0) = f(0,u0) =
 u2(0)u3(0)
1
3
[
u3(0) + 5u2(0) + 3u1(0) + e
−0/2
]
 =
 −11
1
3
[1 + 5(−1) + 3(0) + 1]

Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 8
Ejemplo 8
f(0,u0) =
−11
−1

u1 = u0 + hf(x0,u0) =
 0−1
1
+ 0.1
−11
−1
 =
−0.1−0.9
0.9

Para el siguiente paso
u2 = u1 + hf(x1,u1)
f(0.1,u1) =
 u2(0.1)u3(0.1)
1
3
[
u3(0.1) + 5u2(0.1) + 3u1(0.1) + e
−0.1/2
]
 =
 −0.90.9
1
3
[
0.9 + 5(−0.9) + 3(0.1) + e−0.1/2
]

Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 8
Ejemplo 8
f(0,u0) =
 −0.90.9
−0.9829

u2 = u1 + hf(x1,u1) =
−0.1−0.9
0.9
+ 0.1
 −0.90.9
−0.9829
 =
−0.19−0.81
0.8017

Recordar que se busca y(0.2) el cual corresponde con el primer elemento de u2
y(0.2) = −0.19
Solución numérica de sistemas de ecuaciones de primer orden Métodos numéricos de Heun y RK4 clásico para sistemas de EDOs
Método de Heun para sistemas (RK2)
ui+1 = ui +
1
2
h (k1 + k2)
i = 0, 1, 2, . . . , n− 1
donde
k1 = f (xi,ui)
k2 = f (xi + h,ui + hk1)
Método RK4 clásico para sistemas
ui+1 = ui +
1
6
h (k1 + 2k2 + 2k3 + k4)
i = 0, 1, 2, . . . , n− 1
donde
k1 = f (xi,ui)
k2 = f
(
xi +
1
2
h,ui +
1
2
hk1
)
k3 = f
(
xi +
1
2
h,ui +
1
2
hk2
)
k4 = f (xi + h,ui + hk3)
Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 9
Ejemplo 9
Considere la EDO 3y′′′ − y′′ − 5y′ − 3y = e−x/2, sujeta a y(0) = 0, y′(0) = 1, 0 ≤ x ≤ 1.
Utilice el método RK4 para sistemas, con paso h = 0.1 y encuentre un estimado para y1 = y(0.1).
Solución.
ui+1 = ui +
1
6
h (k1 + 2k2 + 2k3 + k4)
k1 = f(x0,u0) = f(0,u0)
k1 =
 u2(0)u3(0)
1
3
[
u3(0) + 5u2(0) + 3u1(0) + e
−0.1/2
]
 =
 −11
1
3
[1 + 5(−1) + 3(0) + 1]

k1 =
−11
−1

Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 9
Ejemplo 9
k2 = f
(
x0 +
1
2
h,u0 +
1
2
hk1
)
= f
(
0.05,u0 +
1
2
hk1
)
u0 +
1
2
hk1 =
 0−1
1
+ 0.5(0.1)
−11
−1
 =
−0.05−0.95
0.95

k2 =
 −0.950.95
1
3
[
0.95 + 5(−0.95) + 3(−0.05) + e−0.05/2
]
 =
 −0.950.95
−0.9916

Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 9
Ejemplo 9
k3 = f
(
x0 +
1
2
h,u0 +
1
2
hk2
)
= f
(
0.05,u0 +
1
2
hk2
)
u0 +
1
2
hk2 =
 0−1
1
+ 0.5(0.1)
 −0.950.95
−0.9916
 =
−0.0475−0.9525
0.9504

k3 =
−0.95250.9504
−0.9931

k4 = f (x0 + h,u0 + hk3) = f (0.05,u0 + hk3)
u0 + hk3 =
 0−1
1
+ (0.1)
−0.95250.9504
−0.9931
 =
−0.0953−0.9050
0.9007

k4 =
−0.90500.9007
−0.9863

Solución numérica de sistemas de ecuaciones de primer orden Ejemplo 9
Ejemplo 9
Finalmente
u1 = u0 +
1
6
h (k1 + 2k2 + 2k3 + k4)
u1 =
 0−1
1
+ 1
6
(0.1)
−11
−1
+ 2
 −0.950.95
−0.9916
+ 2
−0.95250.9504
−0.9931
+
−0.90500.9007
−0.9863

u1 =
−0.0952−0.9050
0.9007

La solución y(0.1)
y(0.1) = −0.0952
Solución numérica de problemas con valores en la frontera (PVF)
Solución numérica de problemas con valores en la frontera (PVF)
Un PVI se refiere a cuando una EDO de orden n es acompañada de n condiciones iniciales especificadas en
el mismo valor de la variable independiente
Los PVF son aquellos en que los valores de la variable independiete, en las condiciones iniciales, se pueden
especificar en distintos valores (normalmente en los extremos)
Se presentarán dos métodos de solución
Método de disparo.
El PVF se transforma en un PVI al adivinar en las condiciones iniciales
El PVI se resuelve y se prueba la solución para verificar que las condiciones de fronte originales se satisfagan.
Este método utiliza las técnicas RK4
Método de diferencias finitas
Se divide el intervalo del sistema en varios subintervalos
Se sustituyen las derivadas por aproximaciones en diferencias finitas
Se obtiene un sistema de ecuaciones algebraicas y se resuelve
Solución numérica de problemas con valores en la frontera (PVF)
PVF de segundo orden
Considere una ecuación diferencial de segundo orden en su forma general
y′′ = f(x, y, y′), a ≤ x ≤ b
Sujeto a dos condiciones de frontera, normalmente especificads en los extremos a y b.
Se pueden proporcionar distintas condiciones de frontera a los valores de y o y′
Condiciones de frontera
Dirichlet. Se dan valores de y en los extremos
y(a) = ya; y(b) = yb
Neumann. Se dan valores de y′ en los extremos
y′(a) = y′a; y
′(b) = y′b
Mixtas
c1y
′
1(a) + c2y(a) = Ba; c3y
′
1(b) + c4y(b) = Bb
Solución numérica de problemas con valores en la frontera (PVF) Método de disparo
Método de disparo
Procedimiento:
Un PVF que involucra EDO de orden n viene con n condiciones de frontera
Usando variables de estado, la EDO de orden n se transforma en un sistema de EDO de primer orden
Resolver este nuevo sistema requiere n condiciones iniciales
Las condiciones de frontera dadas en el extremo izquierdo del intervalo, sirven como condiciones iniciales
ah́ı, mientras las restantes se deberán adivinar
El sistema EDO se resuelve numéricamente con RK4 y el valor de la solución resultante en el otro punto
extremo se compara con las condiciones de frontera ah́ı.
Si la exactitud no es aceptable, las condiciones iniciales que se adivinaron, se vuelven a adivinar y se resuelve
de nuevo elsistema
El proceso se repite hasta que la solución en el extremo derecho coincida con la condición de frontera
preescirta ah́ı.
Solución numérica de problemas con valores en la frontera (PVF) Ejemplo 10
Ejemplo 10
Resuelva el siguiente PVF utilizando el método de
disparo
ü = 0.02u+ 1; u(0) = 10; u(10) = 100
Solución.
Es una EDO de segundo orden, se transforma a un
sistema de EDO con varibales de estado x1 = u, x2 = u̇.
Se obtiene
ẋ1 = x2
ẋ2 = 0.02x1 + 1
Con condiciones iniciales
x1(0) = 10
x2(0) = ?
Estrategia para hallar x2(0)
Se adivina el valor para x2(0) y se resuelve el
sistema
De la solución se extrae la solución de la primer
variable de estado x1
El valor de esta variable del lado derecho [u(10)] se
compara con la condición de frontera en dicho
extremo [u(10) = 100]
Se adivina otro valor para x2(0) y se repite el
proceso
Para cada adivinanza de x2(0) se encuentra un
valor para u(10)
Como la EDO es lineal, estos últimos valores
también son lineales
Una interpolación lineal de dichos valores proveerá
un valor único para x2(0) que resultará en la
solución correcta
Solución numérica de problemas con valores en la frontera (PVF) Ejemplo 10
Ejemplo 10
Como primer aproximación se elige x2(0) = 10, tal
que x1(0) = 10, x2(0) = 10
Se resuelve numéricamente y se obtiene
u(10)a1 = 217.5208
El valor es mayor a u(10) = 100
Se elige x2(0) = 0 y se ejecuta el método numérico
u(10)a2 = 80.6910
Este punto está por debajo de u(10) = 100
Se utilizan los puntos obtenidos para interpolar de
la siguiente manera
x = [80.6910 217.5208]
y = [0 10]
Se interpolará para xi = 100
Se encuentra que yi = 1.4112
Lo anterior implica
u2(0) = 1.4112
Se vuelve a ejecutar el método con
x1(0) = 10, x2(0) = 1.4112
Se verifica que la solución satisface x2(10) = 100
Solución numérica de problemas con valores en la frontera (PVF) Ejemplo 10
Observaciones:
Un PVI lineal es relativamente sencillo de resolver usando el método de disparo
Los datos obtenidos con dos valores inciales se pueden interpolar linealmente para obtener el estimado
correcto
Lo anterior no es suficiente para un PVI no lineal
Solución numérica de problemas con valores en la frontera (PVF) Método de diferencias finitas
Método de diferencias finitas
Es comunmente utilizado como alternativa al método de disparo
El intervalo [a, b] se divide en N subintervalos de paso h = b−a
N
Se obtendrán N + 1 puntos, con x1 = a, xn+1 = b como puntos extremos y x2, . . . , xN los puntos internos
En cada punto interno, las derivadas involucradas en la EDO se transformarán en un sistema de N − 1
ecuaciones algebraicas que pueden ser resueltas por cualquier método de solución de ecuaciones lineales
Normalmente se utilizan fórmulas de diferencias centradas para aproximar las derivadas, ya que son más
exactas
Para la primer y segunda derivada
dy
dx
∼=
yi+1 − yi−1
2h
d2y
dx2
∼=
yi−1 − 2yi + yi+1
h2
Solución numérica de problemas con valores en la frontera (PVF) Ejemplo 11
Ejemplo 11
Sea ü = 0.02u+ 1, sujeto a u(0) = 10, u(10) = 100
Resuelva usando el método de diferencias finitas, utilizando diferencias centradas para aproximar las derivadas y
considere h = ∆t = 2.
Solución.
En cada punto interior de la malla, la segunda derivada se reemplaza con una fórmula de diferencias
centradas para obtener
ui−1 − 2ui + ui+1
∆t2
= 0.02ui + 1
Se tendrán N = 10/2 = 5 intervalos, tal que se tienen N − 1 = 4 puntos internos
Simplificando la ecuación anterior
ui−1 −
(
2 + 0.02∆t2
)
ui + ui+1 = ∆t
2; i = 2, 3, 4, 5
Solución numérica de problemas con valores en la frontera (PVF) Ejemplo 11
Ejemplo 11
¡Recuerde que u1 = 10, u6 = 100 se conocen de las condiciones de frontera!
Generando el sistema de ecuaciones a partir de la ecuación previa se tiene
ui−1 −
(
2 + 0.02∆t2
)
ui + ui+1 = ∆t
2
u1 − 2.08u2 + u3 = 4
u2 − 2.08u3 + u4 = 4
u3 − 2.08u4 + u5 = 4
u4 − 2.08u5 + u6 = 4
⇒
10− 2.08u2 + u3 = 4
u2 − 2.08u3 + u4 = 4
u3 − 2.08u4 + u5 = 4
u4 − 2.08u5 + 100 = 4
⇒
−2.08u2 + u3 = −6
u2 − 2.08u3 + u4 = 4
u3 − 2.08u4 + u5 = 4
u4 − 2.08u5 = −96
En forma matricial se tiene
−2.08 1 0 0
1 −2.08 1 0
0 1 −2.08 1
0 0 1 −2.08


u2
u3
u4
u5
 =

−6
4
4
−96

Resolviendo el sistema de
ecuaciones lineales
u2 = 15.0489
u3 = 25.3018
u4 = 41.5788
u5 = 65.1821
	Solución numérica de ecuaciones diferenciales de primer orden
	Método de Euler
	Ejemplo 1
	Método de Taylor de orden superior
	Ejemplo 2
	Métodos de Runge-Kutta de segundo orden
	Ejemplo 3
	Métodos de Runge-Kutta de tercer orden (RK3)
	Ejemplo 4
	Métodos de Runge-Kutta de cuarto orden (RK4)
	Ejemplo 5
	Solución numérica de sistemas de ecuaciones de primer orden
	Ecuaciones diferenciales en variables de estado
	Ejemplo 6
	Ejemplo 7
	Método de Euler para sistemas
	Ejemplo 8
	Métodos numéricos de Heun y RK4 clásico para sistemas de EDOs
	Ejemplo 9
	Solución numérica de problemas con valores en la frontera (PVF)
	Método de disparo
	Ejemplo 10
	Método de diferencias finitas
	Ejemplo 11