Catedra 12

Algoritmos simplecticos

  • Verlet
  • Soermer

Estos algoritmos tienen aplicacion a ecuaciones de segundo orden tales que la primera derivada no aparece a la derecha. Ademas conservan la energia del sistema.

$$ \begin{matrix} \frac{d^2y}{dx^2} = f(x, y) \\ f = ma \end{matrix} $$

Verlet Clasico

$$ \begin{matrix} y'' = f(x,y) = a \\ \\ y_{n+1} = y_n + \frac{dy}{dx}\biggr\vert_n+\frac{1}{2}\frac{d^2y}{dx^2}h^2 + O\left(h^3\right) \qquad\qquad\\ \\ y_{n+1} = y_n + v_nh + \frac{1}{2}a_nh^2 + O\left(h^3\right) \qquad\qquad (1)\\ \\ y_{n-1} = y_n - v_bh + \frac{1}{2}a_nh^2 + O\left(h^3\right) \qquad\qquad (2) \end{matrix} $$

Sumando las ecuaciones (1) y (2) $\implies$ $$ \begin{matrix} y_{n+1} + y_{n-1} = 2y_n + a_nh^2 + O\left(h^4\right) \\ \\ y_{n+1} = 2y_n-y_{n-1} + a_nh^2 + O\left(h^4\right) \end{matrix} $$

Condiciones iniciales:

  • $y(x_0) = y_0$
  • $y'(x_0) = v_0$

Ventajas:

  • No depende de $v_n$
  • Rapido

Desventajas:

  • No es ''Self-starting'' i.e. requiere de otro algoritmo par empezar.

Verlet con velocidad

$$\begin{matrix} v_{n+1} = v_n + \frac{1}{2}\left(a_n + a_{n+1}\right)h \\ \\ y_{n+1} = y_n + v_nh + \frac{1}{2}a_nh^2 \end{matrix}$$

Asi $v_{n+1} = \frac{y_{n+2} - y_n}{2h}$

Algoritmo de Beeman

Es lo mismo que verlet ($O\left(h^4\right)$) pero conserva mejor la energia del sistema. $$\begin{matrix} v_{n+1} = v_n + \frac{1}{6}\left(2a_{n+1} + 5a_n - a_{n-1}\right)h + O\left(h^3\right) \\ \\ y_{n+!} = y_n + v_nh + \frac{1}{6}\left(4a_n - a_{n-1}\right)h^2 + O\left(h^4\right) \end{matrix}$$

Paso adaptativo

Si hay dos o mas escalas de tiempo asociadas al problema, entonces no existe un paso adaptativo para resolver el problema.

Ejemplos

  • 2 resortes conectados y con frecuencia de oscilacion muy distintas
  • $$\begin{matrix} u' = 998u + 1998v \\ \\ v' = -99u - 199v \end{matrix}$$ Con condiciones iniciales $u(0) = 1,\quad v(0) = 0$

Para el segundo ejemplo: Transformacion $$\begin{matrix} u = 2y-z \rightarrow u' = 2y' - z' \\ \\ v = -y+z \rightarrow v' = -y' -z' \end{matrix}$$

Esto da como solucion: $$\begin{matrix} u = 2e^{-x}-e^{-1000x} \\ \\ v = -e^{-x} + e^{-1000x} \end{matrix}$$

Se requiere pasos pequenos para resolver el sistema.

Estrategia

Ejemplo

Consideremos $y' = -cy$ con $c>0$ constante.

$$\begin{matrix} \frac{y_{n+1}-y_n}{h} = -cy_n \\ \\ y_{n+1} = y_n - chy_n \\ \\ y_{n+1} = y_n(1-ch) \end{matrix}$$

Que pasa si $1-ch < -1$? $\implies y_1 \rightarrow\infty$ as $n\rightarrow\infty\implies2<ch$

Cambiando el procedimiento a

$$ y'_{n+1} = \frac{y_{n+1} - y_n}{h} $$

Se tiene

$$\begin{matrix} \frac{y_{n+1} - y_n}{h} = -cy_{n+1} \\ \\ y_{n+1}(1 + ch) = y_n \\ \\ y_{n+1} = \frac{y_n}{1+ch} \\ \\ \rightarrow y_n = \frac{y_0}{\left(1+ch\right)^n} \end{matrix}$$

Esto nunca diverge.

De forma mas general:

$y' = f(t,y)$

Implicito:

$$y_{n+1} = y_n + hf(t_{n+1}, y_{n+1})$$

La ecuacion a resolver es posiblemente no lineal.

En los casos no lineales, la idea es linealizar (al mismo orden que el metodo de discretizacion de la derivada).

Ejemplo:

Dada la ecuacion $$\begin{matrix} y' = f(t,y) \\ y_n = \bar{y}_n \end{matrix}$$ con $y_n = \bar{y}_n$ condicion de borde.

$$ y(t) = \bar{y}_n + \int_{t_n}^t f(t', y)dt' $$

Para $t = t_{n+1} \rightarrow y_{n+1} = \bar{y}_n + \int_{t_n}^{t_{n+1}}f(t', y)dt'$. Esto se puede integrar usando, por ejemplo, la regla trapezoidal. Luego

$$ y_{n+1} = \bar{y}_n + \frac{h}{2}\left[f\left(t_{n+1}, y_{n+1}\right) + f\left(t_n, y_n\right)\right] + O\left(h^2\right) $$
  • En EDP, se llama algoritmo de Crank-Nicolson.

En general la ecuacion anterior es no lineal. El primer termino de la integral se debe expandir hasta orden $O\left(h^2\right)$

$$ f(t_{n+1}, y_{n+1}) = f(t_{n+1}, y_n) + (y_{n+1} - y_n)\frac{\partial f}{\partial y}\biggr\vert_{t_{n+1}, y_n} + O\left(h^2\right) $$

Asi $$\begin{matrix} y_{n+1} = y_n + \frac{h}{2}\left[f(t_{n+1}, y_n) + (y_{n+1} - y_n)\frac{\partial f}{\partial y}\biggr|vert_{t_{n+1}, y_{n}} + O\left(h^2\right) + f(t_n, y_n)\right] + O\left(h^3\right) \\ \\ = y_n + \frac{h}{2}\left[f(t_{n+1}, y_n) + f(t_n, y_n)\right] + \frac{h}{2}y_{n+1}\frac{\partial f}{\partial y} - \frac{h}{2}y_n\frac{\partial f}{\partial y} + O\left(h^3\right) \end{matrix}$$

$$ y_{n+1} = y_n + \frac{ \frac{h}{2}\left[f(t_{n+1},y_n) + f(t_n, y_n)\right]}{1 - \frac{h}{2}\frac{\partial f}{\partial y}\biggr\vert_{t_{n+1}, y_n}} $$