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.

numerical

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).