1361.- Flame (3)

Volvamos al flame. Habíamos dicho que el cociente incremental [y(x+h)-y(x)]/h aproximaba bastante bien y'(x), y que el método de Euler para resolver ecuaciones diferenciales nos daba un error O(h), pensado como un O(h² de Taylor multiplicado por 1/h, el número de pasos.

Si queremos resolver una ecuación tan simple como
y' = ay, y(0)=1
(tan simple que la sabemos resolver exactamente y da e^ax) reemplazamos la derivada por el cociente, llamamos x_k=0 + kh, y tenemos
y(x_k+1) = (1+ah)y(x_k)
Si repetimos esto k veces hacia atrás, queda
y(x_k+1) = (1+ah)^k
porque y(x₀)=y(0)=1. Les dejo que vean que esto converge a la exponencial y=e^ax cuando h tiende a cero...

* * *
El método -Euler me perdone- es malo. Ok., funcionó bien durante años, se lo utilizó como la base de métodos posteriores, y toda la teoría de diferencias finitas arranca ahí. La filosofía es simple, cambiar derivadas por diferenciales (o, como se lo presenta a veces, utilizar el polinomio de Taylor para desarrollar la función incógnita).

Lo malo, para nada obvio, es que la aproximación que hacemos es demasiado general: el método, la aproximación, la filosofía que mencionaba es la misma para cualquier ecuación diferencial. Y por eso el error es O(h): es lo óptimo que se puede conseguir si estamos aproximando la derivada de una función cualquiera.

En otras palabras, dada una ecuación diferencial específica, concreta, es posible que haya métodos mucho mejores que el de Euler -conservando su simplicidad- para esa ecuación. Por ejemplo, para la ecuación anterior, el siguiente:
y(x_k+1)=e^a(k+1)hy(x_k)/e^akh=e^ahy(x_k)
es perfecto! coincide exactamente con la solución en los puntos de la grilla 0, ah, a2h,..., akh,...

* * *
Y esa es la 1ra característica de flame: no hay un método fijo, premoldeado para cualquier ecuación, sino que dependerá mucho de la ecuación que uno quiera resolver.

1357.- Flame (2)

De vez en cuando la vida nos ofrece una ecuación diferencial:
y'(x) = f(x,y), y(x₀)=y₀
irresoluble, que no sale por ninguno de los métodos analíticos.

La 1ra pregunta es ¿qué hacer?
* * *
Hace 250 años Euler tiró una idea [bue, 250... Euler hizo todo en el año 1750 +/-30]:

se puede aproximar la derivada y'(x₀) por el cociente incremental [y(x₀+h)-y(x₀)]/h, y despejando, queda
y(x₀+h) = y(x₀) + h f(x₀,y(x₀))
Parece fácil: conocemos la función en el punto x₀ porque es dato, y ahora con esta cuenta, la conocemos en el punto x₀+h. Llamando a h el paso, hay que ir paso a paso, como dijo Mostaza, y uno llega así a calcular la solución en algún valor x₁ que nos interese.

La 2da pregunta es ¿en cuánto le estamos pifiando?
* * *
Miremos de vuelta el pase de magia: y'(x₀) no es el cociente incremental [y(x₀+h)-y(x₀)]/h, hay un error que podemos llamar T, así que cuando despejamos, la cuenta queda
y(x₀+h) = y(x₀) + h f(x₀,y(x₀))+h T,
y decimos que el error es de orden h.

Uno puede quedarse tranquilo, pensando que cuando h es chico también el error es chico, y que uno entonces puede calcular con precisión casi arbitraria los valores de la ecuación sin muchos líos.

La 3ra pregunta es ¿en serio?
* * *
La respuesta es no. Para entender por qué, uno podría pensar en que x₀=0, y queremos calcular la solución en 1 con cierta precisión. Tomamos como h = 0.1, y hacemos diez pasos, pero resulta que ya en el primero el error es grande. Entonces, tomamos h=0.05, pero ahora tenemos que hacer el doble de pasos, y si bien el error de cada paso es más chico, como tenemos que hacer 20 pasos podría amplificarse. Y si hacemos h=0.01, ahora hay que hacer 100 pasos... Si se quiere, y simplificando mucho, cada vez le pifiamos por menos, pero como tenemos que hacer muchos pasos, la acumulación de esos pequeños errores puede ser terrible.

No voy a entrar en más detalles por ahora, pero que esto no pase es la estabilidad del método, y en el fondo nos lleva al Teorema Fundamental del Análisis numérico: dada la consistencia de un método, es equivalente la convergencia a la estabilidad (la consistencia, para no complicarnos más la vida, es que nuestro método se aproxime a la ecuación diferencial cuando h tiende a cero).

La 4ta pregunta es: ¿y qué error se comete con el método de Euler al calcular valores?
* * *
Un O(h), es de orden h. Si bien h dijimos que era el error de truncamiento, lo calculamos muy a lo bestia. Si hacemos Taylor,
y(x₀+h) - y(x₀)- h f(x₀,y(x₀))=y''(c)h²,
el error real es de orden h². Pero si vamos, por ejemplo, de 0 a 1 con pasitos h, tenemos 1/h pasos, y el error se puede acumular tanto como
y''(c)h² . 1/h = Cte . h

La 5ta pregunta es: ¿es normal no haber entendido un carajo?
* * *
La respuesta es sí. Para entender mejor estos temas conviene agarrar algún libro como la gente, alguno de los apuntes que andan por la red, o aunque más no sea, la wikipedia.

1356.- Flame (1)

Esta semana aprendí qué es flame (también es otra cosa, frecuente en los comentarios de los blogs y los foros, pero eso se hace difícil no saberlo a poco que uno navegue por estos submundos de la web).

Me interesa hablar del otro.

* * *
Si, ya sé que todavía no hablé del cloaking, pero me lo crucé esta semana en el mismo lugar. La idea es tan simple y obvia -y encima, funciona- que resulta genial.

Vamos a ver si le dedico unas líneas la próxima semana.