6.1 EDO con Taylor

Referencia: ejemplo Rodriguez 9.1.1 p335. Chapra 25.1.3 p731 pdf755

En los métodos con Taylor para Ecuaciones Diferenciales Ordinarias (EDO), se aproxima el resultado a n términos de la serie, para lo cual se ajusta la expresión del problema a la derivada correspondiente.

Ejemplo:

Se pide encontrar puntos de la solución en la ecuación diferencial usando los tres primeros términos de la serie de Taylor con h=0.1 y punto inicial x₀=0, y₀=1

y'-y -x +x^2 -1 = 0

---

La solución empieza usando la Serie de Taylor para tres términos:

y_{i+1} = y_{i} + h y'_i + \frac{h^2}{2!} y''_i x_{i+1} = x_{i} + h E = \frac{h^3}{3!} y'''(z) = O(h^3)

Luego, al despejar de la fórmula planteada el valor de  y’, se puede obtener y" a combinar las ecuaciones.

y' = y -x^2 +x +1 y'' = y' -2x + 1 = (y -x^2 +x +1) -2x + 1 y'' = y -x^2 -x +2

ecuaciones que permiten estimar nuevos valores y_i+1 para nuevos puntos  muestra distanciados en i*h desde el punto inicial.

Se empieza evaluando el nuevo punto a una distancia x₁= x₀+h del punto de origen con lo que se obtiene y₁ , repitiendo el proceso para el siguiente punto en forma sucesiva

---

Algoritmo con Python

Para simplificar los calculos se crea una función edo_taylor3t() para encontrar  los valores para una cantidad de muestras distanciadas entre si h veces del punto inicial [x₀,y₀]

# EDO. Método de Taylor 3 términos 
# estima solucion para muestras separadas h en eje x
# valores iniciales x0,y0
# entrega arreglo [[x,y]]
import numpy as np

def edo_taylor3t(d1y,d2y,x0,y0,h,muestras):
    tamano = muestras + 1
    estimado = np.zeros(shape=(tamano,4),dtype=float)
    # incluye el punto [x0,y0]
    estimado[0] = [x0,y0,0,0]
    x = x0
    y = y0
    for i in range(1,tamano,1):
        estimado[i-1,2:]= [d1y(x,y),d2y(x,y)]
        y = y + h*d1y(x,y) + ((h**2)/2)*d2y(x,y)
        x = x+h
        estimado[i,0:2] = [x,y]
    return(estimado)

La función se puede usar en un programa que plantee resolver el problema propuesto:

# PROGRAMA PRUEBA
# Ref Rodriguez 9.1.1 p335 ejemplo.
# prueba y'-y-x+(x**2)-1 =0, y(0)=1

# INGRESO.
# d1y = y', d2y = y''
d1y = lambda x,y: y -x**2 + x + 1
d2y = lambda x,y: y -x**2 - x + 2
x0 = 0
y0 = 1
h = 0.1
muestras = 5

# PROCEDIMIENTO
puntos = edo_taylor3t(d1y,d2y,x0,y0,h,muestras)
xi = puntos[:,0]
yi = puntos[:,1]

# SALIDA
print('estimado[xi, yi, d1yi, d2yi]')
print(puntos)

obteniendo como resultados:

estimado[xi,yi]
[[ 0.          1.        ]
 [ 0.1         1.215     ]
 [ 0.2         1.461025  ]
 [ 0.3         1.73923262]
 [ 0.4         2.05090205]
 [ 0.5         2.39744677]]

Cálculo de Error

La ecuación diferencial ordinaria del ejercicio tiene una solución conocida, lo que permite encontrar el error real en cada punto respecto a la aproximación estimada.

y = e^x + x + x^2

Note que el error crece al distanciarse del punto inicial

# ERROR vs solución conocida
y_sol = lambda x: ((np.e)**x) + x + x**2

yi_psol = y_sol(xi)
errores = yi_psol - yi
errormax = np.max(np.abs(errores))

# SALIDA
print('Error máximo estimado: ',errormax)
print('entre puntos: ')
print(errores)

con los siguientes resultados:

Error máximo estimado:  0.0012745047595
entre puntos: 
[ 0.  0.000170  0.000377  0.000626  0.000922  0.00127 ]

Gráfica

Para visualizar los resultados se usan los puntos estimados y muchos más puntos para la solución conocida en el intervalo de estudio. Con la gráfica se podrán observar las diferencias entre las soluciones encontradas.

# GRAFICA [a,b+2*h]
a = x0
b = h*muestras+2*h
muestreo = 10*muestras+2
xis = np.linspace(a,b,muestreo)
yis = y_sol(xis)

# Gráfica
import matplotlib.pyplot as plt
plt.plot(xis,yis, label='y conocida')
plt.plot(xi[0],yi[0],'o', color='r', label ='[x0,y0]')
plt.plot(xi[1:],yi[1:],'o', color='g', label ='y estimada')
plt.title('EDO: Solución con Taylor 3 términos')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.grid()
plt.show()

Tarea: Realizar el ejercicio con más puntos muestra, donde se visualice que el error aumenta al aumentar la distancia del punto inicial [x₀,y₀]