1. INTRODUÇÃO Este trabalho explica algumas das várias ferramentas do cálculo numérico para solucionar problemas matemáticos. Em análise numérica, os métodos de Runge–Kutta formam uma família importante de métodos iterativos implícitos e explícitos para a resolução numérica (aproximação) de soluções de equações diferenciais ordinárias. Em termos mais simples, os métodos numéricos correspondem a um conjunto de ferramentas ou métodos usados para se obter a solução de problemas matemáticos de forma aproximada, sendo aplicados a problemas que não apresentam uma solução exata. Os sucessores do método de Euler (Serie Taylor) foram sobretudo, Carl Runge e por Martin Wilhelm Kutta, estes métodos tornaram-se bastante populares devido às suas propriedades e fácil utilização.
página 1
2. HISTÓRICO DA ORIGEM DO MÉTODO Brook Taylor (Londres, 18 de agosto de 1685 — Londres, 30 de novembro de 1731) foi um matemático britânico. Publicou em 1719 o livro New Principles of Linear Perspective, uma versão melhorada do seu trabalho pioneiro intitulado Linear Perspective de 1715. Obra que foi revisada por John Colson em 1749 e reeditada em 1811. Taylor realizou a primeira investigação satisfatória sobre refração astronômico, é no entanto, mais conhecido pelo seu trabalho sobre as séries que hoje recebem seu nome, publicado em 1715 em Methodus incrementorum directa et inversa. Dito de outra maneira, uma série de taylor é uma expansão de uma série de funções ao redor de um ponto. No desenvolvimento em Taylor, tem-se: y(x0 + h) = y(x0) + y’ (x0).h + y’’(x0).h2/2!+ y’’’(x0).h3/3! + .... Leonhard Paul Euler (Basileia, 15 de abril de 1707 — São Petersburgo, 18 de setembro de 1783) foi um grande matemático e físico suíço de língua alemã que passou a maior parte de sua vida na Rússia e na Alemanha. Euler fez importantes descobertas em campos variados nos cálculos e grafos. Ele também fez muitas contribuições para a matemática moderna no campo da terminologia e notação, em especial para as análises matemáticas, como a noção de uma função matemática. Além disso ficou famoso por seus trabalhos em mecânica, óptica, e astronomia. Euler é considerado um dos mais proeminentes matemáticos do século XVIII. Uma declaração atribuída a Pierre-Simon Laplace manifestada sobre Euler na sua influência sobre a matemática: Em matemática e ciência computacional, o método de Euler, cujo nome relaciona-se com Leonhard Euler, é um procedimento numérico de primeira ordem para solucionar equações diferenciais ordinárias com um valor inicial dado. É o tipo mais básico de método explícito para integração numérica para equações diferenciais ordinárias. página 2
O Método de Euler segue o Método de Taylor, mas se limita à primeira derivada, ignorando da segunda derivada em diante. y(x0 + h) y(x0) + y (x0).h ’
Carl David Tolmé Runge (Bremen, 30 de Agosto de 1856 — Göttingen, 3 de Janeiro de 1927) foi um matemático alemão. Foi inicialmente professor em Hannover e em 1904, sob influência de Felix Klein, foi chamado para Göttingen para a nova cadeira de matemática aplicada (a primeira deste tipo na Alemanha). Já em Hannover contribuiu para a física da espectroscopia. Em Göttingen desenvolveu, juntamente com Martin Wilhelm Kutta, o método de Runge-Kutta para a resolução numérica de problemas de valores iniciais. Famosa é também sua observação de polinómios de interpolação e o seu comportamento quando se aumenta o grau do polinómio. Martin Wilhelm Kutta (Byczyna, 3 de novembro de 1867 — Fürstenfeldbruck, 25 de dezembro de 1944) foi um matemático alemão. De 1885 a 1890 estudou na Universidade de Wrocław, e depois, até
1894, na Universidade de Munique. De 1894 a 1897 foi assistente de Walther von Dyck na Universidade Técnica de Munique. Em 1898 passou meio ano na Universidade de Cambridge. De 1899 a 1909 foi novamente assistente de Walther von Dyck. De 1909 a 1910 foi professor na Universidade Friedrich Schiller de Jena, e de 1910 a 1912 foi professor na RWTH Aachen. Em 1912 foi professor na Universidade de Stuttgart, onde permaneceu até aposentar-se em 1935. Em 1901, baseado em um artigo de Carl Runge, desenvolveu o Método de Runge-Kutta, utilizado para resolver equações diferenciais ordinárias. Runge – Kutta é um dos métodos mais populares para a integração da equação diferencial de primeira ordem. Esses métodos de aproximação de uma função se usam da expansão por série de Taylor. Desta forma, o método página 3
de Runge - Kutta de primeira ordem se utiliza da expansão de Taylor de primeira ordem, o método de Runge - Kutta de segunda ordem se utiliza da expansão de Taylor de segunda ordem, e, assim por diante. Lembrando que o método de Euler é equivalente ao método de Runge - Kutta de primeira ordem. 3. MÉTODO DE RUNGE-KUTTA DE 1ª ORDEM – MÉTODO DE EULER O método de Euler é um método de série de Taylor de 1ª ordem, como
Então,
é o método de Euler que satisfaz as três propriedades mencionadas no item anterior que o caracterizam como um método de Runge-Kutta de ordem p =1. O método de Euler utilizado para resolver EDO com condições iniciais, é o método numérico mais simples. Ele consiste em aproximar a solução y ( x ), no sentido de uma linearização, por meio de suas tangentes.
Podemos, porém, melhorar esta aproximação se subdividirmos o intervalo [x 0 ; z ] em subintervalos de amplitude constante, genericamente página 4
chamada h, e como sabemos calcular a direção da função incógnita y ( x ) em cada ponto, substituiremos tal função por um segmento de reta, em cada um destes subintervalos. Estes segmentos terão a direção que ela (função) tem no início de cada dos subintervalos. Método de Euler considerando dois subintervalos
página 5
4. MÉTODO DE RUNGE-KUTTA CARACTERÍSTICAS
A idéia dos métodos que estudaremos é aproveitar as qualidades dos métodos das séries de Taylor eliminando o seu maior defeito que é o cálculo de derivadas de f(x,y). Os métodos de Runge-Kutta de ordem p caracterizam-se pelas propriedades: o
o
São de passo um (para calcular yi usamos apenas y i-1 ); Não exigem o cálculo de qualquer derivada de f (x, y); no entanto, pagam, por isso, o preço de calcular f (x, y) em vários pontos;
Após expandir f (x, y) por Taylor para função de duas variáveis em torno de (xn, Yn) e agrupar os termos semelhantes, sua expressão coincide com a do método de série de Taylor de mesma ordem. O método de Euler é um método de série de Taylor de 1ª ordem:
Então e assim o método de Euler satisfaz as propriedades acima que o caracteriza como um método de RungeKutta de ordem p=11 . 5. PROBLEMA DE VALOR INICIAL EM EQUAÇÕES ORDINÁRIAS Uma equação diferencial na qual a variável dependente é função de apenas uma variável é dita equação diferencial ordinária. Um problema de valor inicial (PVI) é um problema de evolução, no qual a informação inicial (conhecida) é propagada para o interior do domínio unidimensional pela equação diferencial. Matematicamente, o mais simples dos problemas de valor inicial pode ser apresentado na forma: y f ( x, y ) (1.1) y ( a) ,
página 6
onde
f : IR 2 IR é
função incógnita e
uma função contínua. A função
y y ( x) ( x a) é
a
é o seu valor inicial no ponto a .
Um ponto importante a ser citado é a questão da existência e unicidade de soluções de problema de valor inicial para uma equação diferencial, visto que só faz sentido buscar a sua solução numérica, ou solução aproximada, se de antemão estiver garantida a existência e unicidade de sua solução, pois ao determinar uma solução numérica de um problema que a princípio não possua solução ou que possua, mas não seja única, o processo numérico poderá não convergir, ou convergir para algo que não seja confiável. No caso de um problema de valor inicial de equações ordinárias esta questão está bem resolvida. Suponhamos que a função
f ( x, y )
satisfaça as seguintes condições:
é contínua, onde E ( x, y), x [a, b], y IR n ; Existe uma constante K tal que para todo x [a, b] e quais quer dois valores y e y em IR n f : E IR n
f ( x, y ) f ( x, y ) K y y
Teorema 1.1: Seja f ( x, y ) satisfazendo as condições anteriores e seja um vetor dado. Então, existe exatamente uma função y( x) com as seguintes propriedades:
i. ii. iii.
y y ( x) é contínua e diferenciável para x y( x) f ( x, y ( x)) , x [ a, b] ; y (a) .
em [a, b] ;
Vamos tratar agora de um método conhecido por Método de Diferenças Finitas. O método de diferenças finitas é a discretização do domínio e a substituição das derivadas presentes na equação diferencial, por aproximações destas envolvendo somente valores numéricos da função. Na prática substituise as derivadas pela razão incremental que converge para o valor da derivada quando o incremento tende a zero. Dizemos então que o problema foi discretizado. Vamos utilizar como ferramenta matemática no cálculo de aproximações para as derivadas a série de Taylor que relaciona valores da função e suas derivadas num ponto x com valores dessa mesma função numa vizinhança de página 7
x , ou seja, em ( x h ) com h 0 . Se podemos escrever :
y ( x) , tem derivadas de ordem n 1, em x ,
y ( x h) y ( x) h y ( x)
h ( n 1) (n 1)!
h2 2!
y ( x ) .. .
y ( n 1) ( ),
hn n!
y ( n ) ( x)
x x h.
(1.2).
Tomando n 2 em (1.2), e reescrevendo (1.2) para h e h , respectivamente, obtemos :
y ( x h) y ( x) hy( x)
h2
y( x)
2!
h3 3!
y(1 )
e y ( x h) y ( x) hy ( x)
h2 2!
y ( x)
h3 3!
y ( 2 ).
Subtraindo a última expressão da penúltima obtemos a fórmula centrada ou fórmula de diferenças centradas: y ( x)
y ( x h) y ( x h) 2h
h3 3!
y (1 ) y ( 2 ) .
O primeiro método para solução aproximada de PVI´s que trataremos é o Método de Euler Explícito que utiliza como base a fórmula progressiva para a aproximação da derivada primeira. O primeiro passo é dividir o intervalo [a, b] em N subintervalos iguais, cada um de comprimento h , definindo uma malha Rh xi a (i 1) h, i 1,..., N 1 . Com h
(b a)
N
.
Sejam yi y( xi ) aproximações para y ( xi ) , i 2 , ..., N 1 e y1 . Em cada um dos pontos x2 , x3 , ... , x N 1 , aproximamos a equação diferencial (1.1) por : f ( xi , yi )
yi 1 yi
,
i 1, ..., N
yi 1 yi hf ( xi , yi ),
i 1, ... , N
h
ou,
Perceba que o valor de mesmo quando a função
f
yi
1
é calculado explicitamente em função de
yi ,
é não linear. página 8
Quando aproximamos a derivada primeira em (1.1) por diferenças regressivas obtemos: f ( xi , yi )
yi yi 1 h
,
i 2, ..., N 1
ou
yi 1 yi hf ( xi 1 , yi 1 ),
i 1, ... , N
Esse Método é chamado Método de Euler Implícito. Note que se uma função não linear teremos que o vetor
yi
1
f
for
será definido implicitamente e,
portanto para a sua obtenção deve-se lançar mão de método eficiente para a solução de sistemas de equações algébricas não lineares A comparação entre um método explícito e outro implícito é que em geral o método implícito é mais estável, mas exige um passo h menor para ser obter uma boa aproximação da solução. Vamos mostrar os Métodos de Euler Explícito e Implícito para resolver numericamente o PVI da forma (1.1) : Método de Euler Explícito: Vamos ilustrar os resultados obtidos pelos Métodos de Euler Explícito e Implícito para resolver numericamente o seguinte PVI: y x y . y(0) 1
(1.3)
A solução exata do PVI é y ( x) x 1 2e x , que por sua vez foi utilizada na comparação da eficiência dos métodos para x entre 0 e 1 e com h = 0,1. Veja comparação na tabela: Solução exata: x
0,0
y
1,0
0,1
0,2
0,3
0,4
0,5
0,6
0,7
0,8
0,9
1,0
0,9097 0,8375 0,7816 0,7406 0,7131 0,6976 0,6932 0,6987 0,7131 0,7358
Solução Numérica: Método de Euler Explícito x
0,0
y
1,0
0,1
0,2
0,3
0,4
0,5
0,6
0,7
0,8
0,9
1,0
0,9000 0,8200 0,7580 0,7122 0,6810 0,6629 0,6566 0,6609 0,6748 0,6974
página 9
Solução Numérica: Método de Euler Implícito x
0,0
y
1,0
0,1
0,2
0,3
0,4
0,5
0,6
0,7
0,8
0,9
1,0
0,9182 0,8529 0,8026 0,7660 0,7418 0,7289 0,7263 0,7330 0,7482 0,7711
Graficamente obtivemos :
Também usamos o método das diferenças centradas para obter a solução numérica do PVI (1.3).
página 10
6. MÉTODO DE RUNGE-KUTTA DE 2ª ORDEM O objetivo deste método é aproximar a solução de uma equação diferencial de 1º ordem y f x, y , sendo y x a função solução da equação e x a variável desta função solução. Neste método considera-se uma solução inicial conhecida do tipo y0 , a partir desta solução inicial calculam-se as k iteradas que nos vão devolver o valor aproximado da função solução y x nos pontos xn xn1 h y x0
com
n 1,..., k
e h o espaçamento considerado entre as variáveis, obtendo-se
desta forma a construção da função solução y x por pontos. Ao conjunto da equação diferencial e da solução inicial: y f x, y y x0 y 0 Método de Runge-Kutta de ordem 2
y n 1 y n h 1 f x n , y n 2 f x n 3 h, y n 3 hf x n , y n 3 4 4 3 x x h n 1 n
Os métodos de Runge-Kutta de 2ª or-dem devem ter fórmulas que devem ser semelhantes às fórmulas do Método de Taylor até termos de segunda ordem em h. Consideremos o método de Euler aprimorado ou fórmula de Heun: yn 1 yn
f xn , yn f xn 1 , yn 1 h 2
(2)
é semelhantes às fórmulas do Método de Taylor até termos de segunda ordem em h. Da fórmula de Taylor de y(x) em x=xn+1 h2 y ( xn1 ) y ( xn ) y ( xn ) h y( xn )
yn 1 yn y n h yn
2!
h2 2!
página 11
Logo, o Método de Euler Aprimorado é um método de Taylor de 2ª ordem e portanto, devido às 3 propriedades verificadas, também é um Método de Runge-Kutta de 2ª ordem.
7. PARA QUE SERVE O MÉTODO Trajetórias balísticas, trajetória dos satélites artificiais,
Estudo de redes elétricas, Resolver problema deflexão de uma viga, Estabilidade de aviões, Teoria das vibrações, Reações químicas Exemplo 01 - Trajetórias balísticas
Procedimento numérico para calcular o movimento da partícula consiste em interar as equações a partir das condições iniciais. Dados x0 e v0 em t0, obtemos x1 e v1 em t1, e daí x2 e v2 em t2, e assim por diante. Se o movimento é em duas ou três dimensões, o método continua o mesmo devemos apenas escrever as equações na forma vetorial.
página 12
Exemplo - 02 - O problema da deflexão de uma viga.
O problema da deflexão de uma viga foi estudado por meio do método do disparo linear, um método para soluções de problemas de equações diferenciais com valores de limite, ou seja, equações diferenciais com condições impostas em diferentes pontos. Teorema: Seja o problema linear com valor de limite y' ' p( x) y' q( x) y r ( x), a x b y(a) , y(b)
(P)
satisfazendo (i ) p ( x ) , q ( x ) e (ii) q( x) >0
r ( x)
são contínuas em
[a, b] ,
em [a, b] .
Então o problema tem uma única solução. O problema do disparo linear para equações lineares consiste na substituição do problema linear por dois problemas de valor inicial. Dado o problema com valor de limite ( P)
y' ' p( x) y'q( x) y r ( x) y(a) , y(b)
Considere os dois problemas de valor inicial
Se
( P1 )
y' ' ( x) p( x) y' q( x) y r ( x), a x b y ( a ) , y ' (a ) 0
( P2 )
y' ' ( x) p( x) y'q( x) y, a x b y(a) 0, y' (a) 1
y1 ( x )
é a única solução de
( P1 )
e y2 ( x) é a única solução de ( P2 ) então página 13
(1)
y ( x ) y1 ( x )
é solução de
y1 ( x ) y2 ( x ) y2 (b)
( P) .
Foram utilizados métodos de partida para aproximar a solução de
y1 ( x )
e y2 ( x) , assim podemos aproximar a solução do problema com valor de limite y ( x) .
O método de partida utilizado foi o método de Runge-Kutta para resolver um problema de valor inicial y' f (t , y) , y (t0 ) y0 , que envolve uma média ponderada de valores de f (t , y) no intervalo
tn t tn 1 , (n 0,1,2,...) .
É dado
por kn1 2kn 2 2kn 3 kn 4 6
yn 1 yn h
onde h é o tamanho do passo e k n1 f tn , yn ,
k n 2 f tn
h h , yn kn1 , 2 2
k n 3 f tn h , yn h k n 2 , 2 2
kn 4 f tn h, yn hkn 3 .
O método do Disparo para equações lineares se baseia na substituição do problema linear com valor de limite por dois problemas com valor inicial
( P1 )
e ( P2 ) . Temos muitos métodos com os quais se podem aproximar as soluções y1 ( x )
e y2 ( x) e, uma vez que contamos com essas aproximações, a solução
do problema com valor de limite é aproximada por meio da equação
y ( x) .
página 14
8. EXEMPLO APLICAÇÃO DO MÉTODO RUNGE-KUTTA 1ª E 2ª ORDEM MODELAMENTO MATEMÁTICO PARA PLANEJAMENTO REGIONAL. A capacidade de avaliar quantitativamente a evolução de um ecossistema, principalmente no tocante a alterações potencialmente adversas para a natureza e para o ser humano, permite que a sociedade possa intervir de forma precoce no sentido de evitar transtornos futuros. Como exemplo, apresenta-se aqui o seguinte modelo de poluição de um lago. O problema consiste em fazer previsões para futuros investimentos industriais numa região composta por uma lagoa ou represa de água salgada, onde a principal fonte financeira é o cultivo de frutos do mar. Este problema foi formulado e resolvido por Boynton , Hawkins e Gray e consiste em planejar que tipos de investimentos devem ser alocados para um crescimento adequado das cidades ao redor da lagoa sem prejudicar o meio biológico existente. O modelo é bastante complexo e envolve 16 variáveis que medem desde aspectos biológicos até aspectos sobre o comprometimento do turismo regional. Para efeito de simplificação, apresentaremos uma situação hipotética com 13 variáveis, excluindo as variáveis que representam o turismo regional. Os valores também são hipotéticos, onde a finalidade aqui é apenas mostrar o comportamento e inter-relações destas variáveis e não seu resultado numérico. O modelo é bastante interessante e nele são consideradas as variáveis econômicas, biológicas e populacional. A idéia é apresentar, sobre determinados parâmetros conhecidos, o comportamento de cada variável frente às variações globais durante um determinado período de tempo. O modelo é da forma:
página 15
onde: x1:área urbanizada; x2:área possível de ser urbanizada; x3:capital local; x4:estrutura da cidade; x5:imagem da cidade; x6:residentes; x7:capital industrial; x8:estrutura industrial; x9:frutos do mar; x10:matéria orgânica despejada; x11:toxinas; x12:nutrientes; x13:coliforme fecal. Os parâmetros I e k são parâmetros de proporcionalidade entre as variáveis e dependem das séries históricas e estatísticas sobre o comportamento das variáveis durante um longo período de tempo Simulação página 16
Foram adotadas como condições iniciais para este exemplo valores hipotéticos e admensionais, servindo apenas para uma demonstração didática do modelo. Assim:
Olhando para as variáveis biológicas na Figura 1 percebe-se que a estratégia utilizada não é recomendável, uma vez que, apesar da matéria orgânica despejada no lago e os nutrientes diminuírem radicalmente a zero em 4 anos, a quantidade de toxinas e coliforme fecal aumentaram abruptamente em conseqüência do grande aumento de residentes (Figura 3) atraídos pelos investimentos industriais (Figura 3), não estando a região adequadamente preparada do ponto de vista do saneamento básico (Figura 2) onde o valor final da parte estrutural da cidade foi menor do que no início. Outro responsável pela poluição deste lago é apontado no indicador de estrutura industrial, que de zero teve um excelente aumento indo a mais de 3000 (Figura 3). Porém, como este desenvolvimento não foi auto-sustentável, ou seja, não houve uma adequação com o meio ambiente, percebe-se que as industrias instaladas tiveram prejuízo ao cabo de 4 anos (Figura 3), por estas estarem inteiramente ligadas aos produtos do lago, que foram a zero (Figura 2). Isto contrasta com o capital local que teve aumento ao fim de 4 anos, mas um aumento que sairá caro para a região, uma vez que sua fonte srcinal de renda voltada à biosfera do lago foi exterminada pelo mau planejamento urbano ao receber industrias e novos residentes.
página 17
A Importância na Escolha do Método Numérico A escolha de um método numérico para a simulação de sistemas dinâmicos é de fundamental importância para uma conclusão segura e coesa sobre o evento. Mostraremos primeiramente como funciona o método de Euler. Suponhamos que tenhamos que resolver a equação diferencial.
página 18
Neste caso, teremos como fazer uma verificação sobre a precisão ou não do método escolhido uma vez que a solução analítica é facilmente obtida e será:
Como já foi vista, a fórmula para o método de Euler é:
e escolhendo a variação no tempo
Poderemos passo a passo obter, escolhendo como h = 0.1 de t = 0 até t = 1,
página 19
9. CONCLUSÃO O melhor método numérico em termos de precisão, para o cálculo de soluções aproximadas de PVI’s foi o método de Runge-Kutta; No Método de Euler Explícito conseguimos estabelecer a convergência do método, fornecendo um limite superior para o erro. No Método das Diferenças Finitas para PVC o método implementado usando diferenças centradas no interior da malha e diferenças progressivas e regressivas no contorno proporcionou resultados numéricos com ótima concordância com as soluções analíticas.
10. REFERÊNCIAS BIBLIOGRAFICAS 1. BOYCE, W. e DIPRIMA, R. C.; Equações Diferenciais Elementares e Problemas de Valores de Contorno. 7ºed. Rio de Janeiro: LTC 2002. 2. RICHARD L. BURDEN e J. DOUGLAS FAIRES; Análise Numérica; São Paulo: Pioneira Thomson Learning, 2003. 3. CONTE, S.D. Elementos de Análise Matemática, 3ª Edição. Editora Globo, 1977. p. 244-317. 4. www.wikipedia.org 5. PMR 2420 – Mecânica Computacional 33 CAPÍTULO III MÉTODOS DE RUNGE-KUTTA - Autor não informado 6. Renato S. Silva, Regina C. Almeida, Métodos Numéricos Equações Diferenciais Ordinárias, p. 10-29
página 20
