6.1.1 Equações Diferenciais Ordinárias — Métodos de Runge-Kuta — Método de Euler

1. Introdução

O método de Euler é um algoritmo destinado à solução de equações diferenciais na forma

$$ \dfrac{dy}{dx} = f(x,y) $$

onde $y = y(x)$ . A solução analítica da equação diferencial ordinária (EDO) é justamente encontrar a equação que expressa $y(x)$ . Por outro lado, quando utilizamos uma abordagem numérica, a solução retorna $y(x)$ avaliado em $x$.

Os métodos de Runge-Kuta são iterativos sobre a estimativa das variáveis e buscam aproximar a solução $y = y(x)$ através de sucessivos passos. Temos que

$$ y_i = y_{i-1} + \phi h $$

onde $\phi$ é a estimativa da inclinação (derivada), usada como estimadora dos novos valores ($y\_i$) a partir de valores prévios ($y\_{i-1}$) em uma distância $h$ (extrapolada). A figura abaixo mostra esta extrapolação.

Todos os métodos de passo único possuem esta abordagem, sendo que se diferenciam pela forma em que $\phi$ é estimado. Assim, a forma mais ingênua consiste em utilizar a própria equação diferencial (com o $y’$ isolado) como inclinação da reta. Isto é, a inclinação no início do intervalo é tomada como uma inclinação média de todo intervalo. Tal abordagem simplista é utilizada no método de Euler. Outras abordagem são conhecidas como métodos de Runge-Kuta de ordens superior.

 

2. Método de Euler

Conforme comentamos na seção anterior, o método de Euler utiliza a primeira derivada como estimativa direta da inclinação em $x_i$ . Isto é,

$$ \phi = \dfrac{dy}{dx} = f(x,y) $$

Assim, a solução será dada por

$$ y_i = y_{i-1} + f(x_{i-1}, y_{i-1}) h $$

 

2.1. Exemplo

Supondo que precisamos resolver a seguinte equação:

$$ \dfrac{dy}{dx} = x + y $$

Sabemos que esta equação está na forma

$$ \dfrac{dy}{dx} – ay = g(x) $$

que possui solução

$$ y = e^{-ax} \int_{x_0}^{x} e^{as} g(s) ds + c e^{-at} $$

onde $a=-1$ e $g(x)=x$ . Para $x_0 = 0$ , temos que a solução analítica da EDO em questão será

$$ y(x) = -1 -x + c e^x $$

Ou seja, a solução é sempre integral, por isso o método de Euler soluciona em um intervalo de integração. Assim, supondo que queremos integrar $f(x,y)$ no intervalo $x \in [0,1]$ com condições iniciais $x=0$, $y(x=0)=-1-0+ce^0=0$ e passo de tamanho 0.5, temos que

$$ y(0.5) = y(0) + f(0, 1) 0.5 $$

onde $y=1$ , então $f(0.1) = 0 + 1 = 1$ . Portanto, a aproximação no primeiro passo será

$$ y(0.5) = 1 + 1 * 0.5 = 1.5 $$

Novamente repetimos este processo, seguindo para $y(0.5 + h) = y(1)$, refinando a aproximação. Como esse processo é repetitivo e monótono, utilizamos ele no exemplo de implementação.

 

3. Implementação

 

def integral_euler(f, xi, xe, y0, n=1e6):
    """
    Numerical integration for solve ODE's using Euler's Method

    integral = integral_euler(f, xi, xe, n=1e6)

    INPUT:
      * f: derivative function f(x, y) = dy/dx
      * xi: beginning of integration interval
      * xe: end of integration interval
      * y0: initial estimative for y
      * n: number of points used

    return: \int_{xi}^{xe} f(x)

    Author: Pedro Garcia [sawp@sawp.com.br]
    see: http://www.sawp.com.br

    License: Creative Commons
             http://creativecommons.org/licenses/by-nc-nd/2.5/br/

    Sep 2012
    """
    def euler(x, y, h):
        return (x + h, y + f(x, y) * h)

    def integrator(x, y, h, n):
        for i in xrange(n):
            (x, y) = euler(x, y, h)
        return y

    (y, x) = (y0, xi)
    h = abs(xe - xi) / n
    return integrator(x, y, h, int(n))

Um exemplo de utilização dessa função pode ser obtido em http://www.sawp.com.br/code/ode/euler.py.

 

4. Copyright

Este documento é disponível sob a licença Creative Commons. As regras dos direitos de cópia deste conteúdo estão acessíveis em http://creativecommons.org/licenses/by-nc-nd/2.5/br/.

References

[1]Anthony Ralston and Philip Rabinowitz, A First Course in Numerical Analysis (2nd ed.), McGraw-Hill and Dover, (2001).
[2]N.B Franco, Cálculo Numérico, Pearson Prentice Hall (2006).