Análise de estabilidade de Von Neumann

Origem: Wikipédia, a enciclopédia livre.
Ir para: navegação, pesquisa

Na análise numérica, a Análise de estabilidade de Von Neumann (também conhecida como análise de estabilidade de Fourier) é um procedimento usado para verificar a estabilidade de métodos de diferenças finitas quando aplicados em equações diferenciais parciais.[1] A análise é baseada na decomposição de Fourier do Erro numérico e foi desenvolvida no Laboratório Nacional de Los Alamos depois de ter sido brevemente descrita em um artigo de 1947 pelos pesquisadores britânicos John Crank e Phyllis Nicolson.[2] Depois, foi dado um tratamento mais rigoroso ao método em um artigo [3] co-escrito por John von Neumann.

Estabilidade numérica[editar | editar código-fonte]

A estabilidade de métodos numéricos está intimamente associada ao erro numérico. Um método de diferenças finitas é estável se os erros produzidos em um passo de tempo do cálculo não provocam um aumento dos erros à medida que os cálculos avançam. Há 3 classes de métodos numéricos. Um método numérico condicionalmente estável depende de certos parâmetros para que os erros permaneçam limitados e não instabilizem. Se os erros diminuírem ou chegarem a desaparecer, o método numérico é dito como sendo estável. Se, de outra forma, os erros aumentarem com o tempo, a solução numérica irá divergir em relação à realidade (solução exata) e então o método numérico é dito como sendo instável. A estabilidade de métodos numéricos pode ser averiguada pela análise de estabilidade de von Neumann. Para problemas dependentes do tempo, a estabilidade garante que o método numérico produza uma solução limitada sempre que a solução da equação diferencial for limitada. Estabilidade em geral, pode ser dificilmente averiguada, especialmente se a equação em questão for não-linear.

A estabilidade de von Neumann é necessária e sucifiente (tal como utilizado no Teorema de Equivalência de Lax) apenas em certos casos: a EDP e o método de diferenças finitas devem ser lineares; a EDP precisa ter condições iniciais e de fronteira e ser no máximo de segunda ordem.[4] Devido à sua relativa simplicidade, a análise de von Neumann é geralmente usada no lugar de uma análise de estabilidade mais detalhada por ser um bom palpite em relação às restrições dos passos usados no método.

Ilustração do método[editar | editar código-fonte]

O método de von Neumann é baseado na decomposição dos erros em séries de Fourier. Para ilustrar o procedimento, considere Equação do calor unidimensional


  \frac{\partial u}{\partial t} = \alpha \frac{\partial^2 u}{\partial x^2}

definida no intervalo L, e discretizando pelo método FTCS

  
\quad (1) \qquad u_j^{n + 1} = u_j^{n} + r \left(u_{j + 1}^n - 2 u_j^n + u_{j - 1}^n \right)

onde

r = \frac{\alpha\, \Delta t}{\Delta x^2}

e a solução u_j^{n} da equação discretizada se aproxima da solução analítica u(x,t) da EDP na malha.

Definimos o Erro de arredondamento (também chamado de erro de truncamento) \epsilon_j^n como


  \epsilon_j^n = U_j^n - u_j^n

onde u_j^n é a solução da equação discretizada (1) que seria calculada sem os erros de truncamento, e U_j^n é a solução numérica obtida com precisão finita. Já que a solução exata u_j^n deve satisfazer a equação discretizada, o erro \epsilon_j^n também deve satisfazer a equação discretizada.[5] Portanto


  \quad (2) \qquad \epsilon_j^{n + 1} = \epsilon_j^n + r \left(\epsilon_{j + 1}^n - 2 \epsilon_j^n + \epsilon_{j - 1}^n \right)

é uma relação de recorrência para o erro. Equações (1) e (2) mostram que ambos os erros e a solução numérica tem o mesmo crescimento ou decaimento em relação ao tempo. Para equações diferenciais lineares com condições de contorno, a variação espacial do erro pode ser expandida em séries de Fourier, no intervalo L, como


  \quad (3) \qquad \epsilon(x) = \sum_{m=1}^{N/2} A_m e^{ik_m x}

onde o Número de onda k_m = \frac{\pi m}{L} com m = 1,2,\ldots,M, M = \frac{L}{\Delta x} e N sendo uma incógnita. A dependência no tempo do erro é incluída assumindo-se que a amplitude do erro A_m está em função do tempo. Já que o erro tende a crescer ou decair exponencialmente com o tempo, é sensato assumir que a amplitude varia exponencialmente com o tempo; daí


  \quad (4) \qquad \epsilon(x,t) = \sum_{m=1}^{N/2} e^{at} e^{ik_m x}

onde a é uma constante.

Já que a equação diferencial do erro é linear (o comportamento de cada tempo da série é o mesmo que o dá série), podemos então considerar que o crescimento do erro de um dado termo é:


  \quad (5) \qquad \epsilon_m(x,t) = e^{at} e^{ik_m x}

As característica da estabilidade podem ser estudadas usando apenas esta forma para o erro, sem grandes perdas em geral. Para descobrir como o erro varia nos espaços de tempo, substituimos (5) na equação (2), notando que

  • 
\begin{align}
 \epsilon_j^n & = e^{at} e^{ik_m x} \\
 \epsilon_j^{n+1} & = e^{a(t+\Delta t)} e^{ik_m x} \\
 \epsilon_{j+1}^n & = e^{at} e^{ik_m (x+\Delta x)} \\
 \epsilon_{j-1}^n & = e^{at} e^{ik_m (x-\Delta x)},
\end{align}

dando (depois de simplificar)


  \quad (6) \qquad e^{a\Delta t} = 1 + \frac{\alpha \Delta t}{\Delta x^2} \left(e^{ik_m \Delta x} + e^{-ik_m \Delta x} -  2\right).

Usando as identidades


  \qquad \cos(k_m \Delta x) = \frac{e^{ik_m \Delta x} + e^{-ik_m \Delta x}}{2} \qquad \text{and} \qquad \sin^2\frac{k_m \Delta x}{2} = \frac{1 - \cos(k_m \Delta x)}{2}

equação (6) pode ser reescrita como


  \quad (7) \qquad e^{a\Delta t} = 1 - \frac{4\alpha \Delta t}{\Delta x^2} \sin^2 (k_m \Delta x/2)

Definimos então o fator de amplitude


  G \equiv \frac{\epsilon_j^{n+1}}{\epsilon_j^n}

A condição necessária e suficiente para o erro se manter limitado é que, \vert G \vert \leq 1. Contudo,


  \quad (8) \qquad G = \frac{e^{a(t+\Delta t)} e^{ik_m x}}{e^{at} e^{ik_m x}} = e^{a\Delta t}

Daí, das equações (7) e (8), a condição para a estabilidade é dada por


  \quad (9) \qquad \left\vert 1 - \frac{4\alpha \Delta t}{\Delta x^2} \sin^2 (k_m \Delta x/2) \right\vert \leq 1

Para a condição acima manter todos os \sin^2 (k_m \Delta x/2), temos


  \quad (10) \qquad \frac{\alpha \Delta t}{\Delta x^2} \leq \frac{1}{2}

Equação (10) dá o requisito de estabilidade para o método FTCS quando aplicado a equação unidimensional do calor. Ela diz que para um determinado \Delta x,o valor permitido de \Delta t deve ser pequeno o bastante para satisfazer a equação (10).

Referências[editar | editar código-fonte]

  1. Analysis of Numerical Methods by E. Isaacson, H. B. Keller
  2. Crank, J.; Nicolson, P. (1947), "A Practical Method for Numerical Evaluation of Solutions of Partial Differential Equations of Heat Conduction Type", Proc. Camb. Phil. Soc. 43: 50–67, doi:10.1007/BF02127704 
  3. Charney, J. G.; Fjørtoft, R.; von Neumann, J. (1950), "Numerical Integration of the Barotropic Vorticity Equation", Tellus 2: 237–254, http://atoc.colorado.edu/~dcn/ATOC7500/members/Reading/Charney_ndbvm_1950.pdf 
  4. Smith, G. D. (1985), Numerical Solution of Partial Differential Equations: Finite Difference Methods, 3rd ed., pp. 67–68 
  5. Anderson, J. D., Jr.. Computational Fluid Dynamics: The Basics with Applications. [S.l.]: McGraw Hill, 1994.

Ver também[editar | editar código-fonte]