Método de Euler: Solución Numérica de Ecuaciones Diferenciales
Contenido del curso
Los modelos cuya solución está basada en un algoritmo numérico los denominamos modelos numéricos. El algoritmo numérico más sencillo que existe es el método de Euler el cual discutiremos a continuación.
Método de Euler
Cuando tenemos una ecuación diferencial de primer orden para una función desconocida y(t) lo más general posible, escribimos:
Donde f(y, t) representa una expresión que contiene todos los demás términos de la ecuación diferencial que no acompañen a la derivada. Ahora, como vimos con el problema de la propagación de una epidemia, un sistema dinámico debe tener un punto de inicio o condición inicial. En general siempre debemos definir cual es esta condición inicial para nuestro sistema así:
Con y₀ siendo el valor inicial de la variable que estamos buscando predecir. En este sentido, el método de Euler propone construir la solución de la ecuación diferencial paso a paso suponiendo que la derivada se aproxima por medio de diferencias finitas pero muy pequeñas:
Estas diferencias finitas representan intervalos de tiempo tan pequeños como uno los considere y paso a paso se puede ir calculando el siguiente valor de la variable y con base en el valor del tiempo anterior:
Esto representa un proceso iterativo donde vamos construyendo la solución paso a paso. El precio de esto es que el aproximar una derivada por intervalos finitos implica que en cada paso se genera un error que se acumula en los siguientes pasos. Matemáticamente hablando, la única manera de reducir este error es considerando un intervalo de tiempo 𝚫t ➝ 0, pero computacionalmente no podemos considerar números infinitamente pequeños, entonces procuramos tomarlo tan pequeño como nuestros recursos de cómputo nos lo permitan. Esta es la gran diferencia entre el cálculo de infinitesimales y los métodos numéricos.
Dado que los métodos numéricos son muy sencillos conceptualmente hablando, pero pueden representar procesos que deban repetirse muchas veces, el uso de un computador para su ejecución es lo ideal y por lo tanto hoy en día todos los métodos numéricos se implementan en un lenguaje de programación para poder realizar múltiples pasos en tiempos razonablemente cortos y dejándole todo el trabajo a una máquina. Como dato curioso estos algoritmos eran calculados a mano antes de que aparecieran las grandes computadoras, pero ya te imaginarás lo engorroso que eso era en aquella época.
Solución de Euler en Python
Tomaremos como ejemplo una ecuación diferencial sencilla junto con una condición inicial que usaremos como base para aplicar el método de Euler:
Es decir esto representa un problema donde la solución será una función que es igual a su derivada y que su valor para t = 0 sea 1. El algoritmo de Euler en este caso quedaría de la siguiente manera:
Para inicial el algoritmo debemos definir el intervalo de tiempo, digamos que tomamos 𝚫t = 0.01 y calculamos los pasos así:
Y así sucesivamente hasta el número de pasos que sea necesario. Por ejemplo si necesitamos construir la solución hasta un tiempo t = 10 , requerimos calcular al menos 1000 pasos con un 𝚫t = 0.01. Esto es mejor hacerlo en un computador, que es justamente lo que haremos a continuación en las siguientes líneas de código:
De aquí, la solución numérica correspondiente se vé así:
Y para los que ya saben que la solución exacta (la que obtenemos de los métodos del cálculo) de esta ecuación ecuación diferencial es la función exponencial:
Se darán cuenta que al comparar ambos tipos de soluciones, se observa una pequeña diferencia que va creciendo conforme avanzamos en el tiempo y que es mayor si escogemos intervalos de tiempo más largos:
Y si quieres seguir jugando con este código te darás cuenta de que esta diferencia desaparece a medida que el intervalo de tiempo es cada vez más pequeño, solo que cada vez tendrás que realizar más pasos.
Y es así como vemos una alternativa para obtener soluciones a modelos de sistemas dinámicos a partir del paradigma numérico, con la advertencia que los métodos numéricos tienen el problema que el error es directamente proporcional al intervalo de tiempo que se usa y esto, a su vez, es inversamente proporcional al número de operaciones que deben realizarse. Es decir, que para obtener mejores soluciones es preciso realizar más pasos y por lo tanto consumir más recursos computacionales.
El notebook de esta clase lo encuentras en este link de Google Colab. En nuestra próxima clase discutiremos en detalle las diferencias entre las soluciones exactas que obtenemos con álgebra y cálculo contra las soluciones numéricas.
Cristian Iván Solórzano Maya
Roger Christian Cansaya Olazabal
Lorena Mondragón
Roger Christian Cansaya Olazabal
Ricardo Posada
Eduardo Monzón
Anthony Jean Paul Blaz Lazo
Muy buen aporte respecto a los métodos numéricos, sería excelente ampliar o crear un curso de análisis numérico y optimización, para la ruta de data science. Pues tienen grandes aplicaciones y al parecer solo los que somos matemáticos estamos al tanto del tema
Por si no queda claro como llegamos a la formula.
import numpy as np import matplotlib.pyplot as plt def exact_sol(ts): return np.exp(ts) def num_sol(ts, dt, tf=10, y0=1): ys = [y0] ts = [0] num_steps = int(tf/dt) for _ in range(num_steps): ts.append(ts[-1]+dt) ys.append((1+dt)*ys[-1]) return ts, ys if __name__=='__main__': ys = [1] # creamos un areglo inicial de valores en y ts = [0] # creamos un areglo inicial de valores en t dt = 0.01 num_steps = 1000 for i in range(num_steps): ts.append(ts[-1]+dt) # calculamos el proximo t y lo agregamos al arreglo ys.append((1+dt)*ys[-1]) # calculamos el proximo y y lo agregamos al arreglo #print (ts) #print(ys) plt.figure(figsize=(15,8)) plt.plot(ts,ys) plt.xticks(fontsize=20) plt.yticks(fontsize=20) plt.show() plt.figure(figsize=(12,6)) for dt in [0.1, 0.05, 0.01]: ts, ys = num_sol(ts, dt) plt.plot(ts, ys, '--', label='dt = {}'.format(dt)) plt.plot(ts,exact_sol(ts), label='exact') plt.xticks(fontsize=20) plt.yticks(fontsize=20) plt.legend(fontsize=20) plt.show()
excelente y elegante código
Nunca perder de vista que al usar modelos nunméricos hay que tener siempre en cuenta el error que se va a generar por usarlos. Este va a ser menor mientras más nos acercamos a lo infinitesimal, pero ello requiere que usemos más capacidad de cómputo. Lo ideal será encontrar un equilibrio entre cantidad de error y recursos a usar, esto dependerá del tipo de problema a resolver.
Genial!