Métodos de Runge-Kutta explícitos

Cálculo Numérico - Versão Scilab

10.7 Métodos de Runge-Kutta explícitos

Nas seções anteriores, exploramos os métodos de Euler e Euler modificado para resolver problemas de valor inicial. Neste momento, deve estar claro ao leitor que o método de Euler melhorado produz soluções de melhor qualidade que o método de Euler para a maior parte dos problemas estudados. Isso se deve, conforme Seção 10.6, à ordem de precisão do método.

Os métodos de Runge-Kutta generalizam o esquema do método de Euler melhorado, inserindo mais estágios de cálculo e buscando ordens de precisão mais altas. Nesta seção, trataremos dos métodos de Runge-Kutta explícitos²

Para tal, considere o problema de valor inicial: $\begin{array}{rcl} u^{'} (t) & = & f (t, u (t)), & (10.135) \\ u (t_{0}) & = & a . & (10.136) \end{array}$

Integrando a EDO em $[t^{(n)}, t^{(n + 1)}]$ obtemos $\begin{array}{rcl} u^{(n + 1)} & = & u^{(n)} + \int_{t^{(n)}}^{t^{(n + 1)}} f (t, u (t)) d t & (10.137) \end{array}$

O método de Euler aproxima a integral no lado direito da expressão acima utilizando apenas o valor de $f (t, u)$ em $t = t^{(n)}$ , o método de Euler melhorado aproxima esta integral utilizando os valores de $f (t, u)$ $t = t^{(n)}$ e $t = t^{(n + 1)}$ . Os métodos de Runge-Kutta explícitos admitem estágios intermediários, utilizando outros valores de $f (t, u)$ nos pontos ${τ_{1}, τ_{2}, \dots, τ_{ν}}$ dentro do intervalo $[t_{n}, t^{(n + 1)}]$ . Veja esquema abaixo:

$u_{n}$			$u_{n + 1}$
$\|$			$\|$

$t^{(n)}$			$t^{(n + 1)}$
$τ_{1}$	$τ_{2}$	$\dots$	$τ_{ν}$
$t^{(n)}$	$t^{(n)} + c_{2} h$	$\dots$	$t^{(n)} + c_{ν} h$

Observe que $τ_{j} = t^{(n)} + c_{j} h$ com $0 = c_{1} \leq c_{2} \leq \dots \leq c_{ν} \leq 1$ , isto é, o primeira ponto sempre coincide com o extremo esquerdo do intervalo, mas o último ponto não precisa ser o extremo direito. Ademais, um mesmo ponto pode ser usado mais de uma vez com aproximações diferentes para $u (τ_{j})$ .

Desta forma, aproximamos a integral por um esquema de quadratura com $ν$ pontos: $\begin{array}{rcl} \int_{t^{(n)}}^{t^{(n + 1)}} f (t, u (t)) d t & \approx & u_{n} + h \sum_{j = 1}^{ν} b_{j} f (τ_{j}, u (τ_{j})) & (10.138) \end{array}$

onde $b_{j}$ são os pesos da quadratura numérica que aproxima a integral. Assim como na construção do método de Euler melhorado, não dispomos dos valores $u (τ_{j})$ antes de calculá-los e, por isso, precisamos estimá-los com base nos estágios anteriores:

\begin{matrix} \begin{aligned} ũ_{1} & = u^{(n)} \\ ũ_{2} & = u_{n} + h a_{21} k_{1} \\ ũ_{3} & = u_{n} + h [a_{31} k_{1} + a_{32} k_{2}] \\ ũ_{4} & = u_{n} + h [a_{41} k_{1} + a_{42} k_{2} + a_{43} k_{3}] \\ ⋮ \\ ũ_{ν} & = u_{n} + h [a_{ν 1} k_{1} + a_{ν 2} k_{2} + \dots + a_{ν ν} k_{ν}] \\ u^{(n + 1)} & = u_{n} + h [b_{1} k_{1} + b_{2} k_{2} + \dots + b_{ν} k_{ν}], \end{aligned} \end{matrix}

(10.139)

onde $k_{j} = f (τ_{j}, ũ_{j})$ , $A : = (a_{i j})$ é a matriz Runge-Kutta (triangular inferior com diagonal zero), $b_{j}$ são os pesos Runge-Kutta e $c_{j}$ são os nós Runge-Kutta. Estes coeficientes podem ser escritos de forma compacta em uma tabela conforme a seguir:

c

A

b^{T}

=

$c_{1}$	$0$	$0$	$0$
$c_{2}$	$a_{21}$	$0$	$0$
$c_{3}$	$a_{31}$	$a_{22}$	$0$

	$b_{1}$	$b_{2}$	$b_{3}$ .

=

$c_{1}$
$c_{2}$	$a_{21}$
$c_{3}$	$a_{31}$	$a_{22}$

	$b_{1}$	$b_{2}$	$b_{3}$ .

Na tabela mais à direita, omitimos os termos obrigatoriamente nulos.

Exemplo 10.7.1. O método de Euler modificado pode ser escrito conforme: $\begin{array}{rcl} ũ_{1} & = & u^{(n)} & (10.140) \\ ũ_{2} & = & u^{(n)} + h k_{1} & (10.141) \\ u^{(n + 1)} & = & u^{(n)} + h [\frac{1}{2} k_{1} + \frac{1}{2} k_{2}] & (10.142) \\ (10.143) \end{array}$

Identificando os coeficientes, obtemos $ν = 2$ , $c_{1} = 0$ , $c_{2} = 1$ , $a_{21} = 1$ , $b_{1} = \frac{1}{2}$ e $b_{2} = \frac{1}{2}$ . Escrevendo na forma tabular, temos:

$c$	$A$


	$b$

$0$
$1$	$1$

	$\frac{1}{2}$	$\frac{1}{2}$

Observação 10.7.1. Nos métodos chamados explícitos, os elementos da diagonal principal (e acima dela) da matriz $A$ devem ser nulos, pois a aproximação em cada estágio é calculada com base nos valores dos estágios anteriores. Nos métodos implícitos, essa restrição é removida, neste caso, o cálculo da aproximação em um estágio envolve a solução de uma equação algébrica.

Observação 10.7.2. Além da condição fixa $c_{1} = 0$ , para que um método seja de ordem pelo menos unitária, isto é, $p \geq 1$ .

10.7.1 Métodos de Runge-Kutta com dois estágios

Os métodos de Runge-Kutta com dois estágios $(ν = 2)$ são da seguinte forma: $\begin{array}{rcl} ũ_{1} & = & u^{(n)} & (10.144) \\ ũ_{2} & = & u^{(n)} + h a_{21} k_{1} & (10.145) \\ u^{(n + 1)} & = & u^{(n)} + h [b_{1} k_{1} + b_{2} k_{2}], & (10.146) \end{array}$

onde $k_{1} = f (t^{(n)}, u^{(n)})$ e $k_{2} = f (t^{(n)} + c_{2} h, ũ_{2})$

Assumindo suavidade suficiente em $f$ , usamos o polinônio de Taylor: $\begin{array}{rcl} k_{2} & = & f (t^{(n)} + c_{2} h, ũ_{2}) & (10.147) \\ = & f (t^{(n)} + c_{2} h, u^{(n)} + a_{21} h k_{1}) & (10.148) \\ = & f (t^{(n)}, u^{(n)}) + h [c_{2} \frac{\partial f}{\partial t} + a_{21} k_{1} \frac{\partial f}{\partial u}] + O (h^{2}) & (10.149) \\ = & k_{1} + h (c_{2} \frac{\partial f}{\partial t} + a_{21} k_{1} \frac{\partial f}{\partial u}) + O (h^{2}) & (10.150) \end{array}$

fazendo com que (10.146) se torne $\begin{array}{rcl} u^{(n + 1)} & = & u^{(n)} + h b_{1} k_{1} + h b_{2} k_{2} & (10.151) \\ = & u_{n} + h b_{1} k 1 + h b_{2} [k_{1} + h (c_{2} \frac{\partial f}{\partial t} + a_{21} k_{1} \frac{\partial f}{\partial u}) + O (h^{2})] & (10.152) \\ = & u_{n} + h (b_{1} + b_{2}) k_{1} + h^{2} b_{2} (c_{2} \frac{\partial f}{\partial t} + a_{21} k_{1} \frac{\partial f}{\partial u}) + O (h^{3}) & (10.153) \end{array}$

Usando a equação diferencial ordinária que desejamos resolver e derivando-a em $t$ , obtemos: $\begin{array}{rcl} u^{'} (t) & = & f (t, u (t)), & (10.154) \\ u^{″} (t) & = & \frac{\partial f}{\partial t} + \frac{\partial f}{\partial u} u^{'} (t) = \frac{\partial f}{\partial t} + \frac{\partial f}{\partial u} f (t, u (t)) . & (10.155) \end{array}$

Agora, expandimos em série de Taylor a solução exata $u (t)$ em $t = t^{(n)}$ , $\begin{array}{rcl} u (t^{(n + 1)}) & = & u (t^{(n)} + h) = u^{(n)} + h u^{'} (t) + \frac{h^{2}}{2} u^{″} (t) + O (h^{3}) & (10.156) \\ = & u^{(n)} + h f (t, u^{(n)}) + \frac{h^{2}}{2} \frac{d}{d t} f (t, u (t)) + O (h^{3}) & (10.157) \\ = & u^{(n)} + h f (t, u^{(n)}) + \frac{h^{2}}{2} (\frac{\partial f}{\partial t} + \frac{\partial f}{\partial u} u^{'} (t^{(n)})) + O (h^{3}) & (10.158) \\ = & u^{(n)} + h f (t, u^{(n)}) + \frac{h^{2}}{2} (\frac{\partial f}{\partial t} + \frac{\partial f}{\partial u} f (t, u^{(n)})) + O (h^{3}) & (10.159) \\ = & u^{(n)} + h f (t, u^{(n)}) + \frac{h^{2}}{2} (\frac{\partial f}{\partial t} + \frac{\partial f}{\partial u} k_{1}) + O (h^{3}) & (10.160) \end{array}$

Finalmente comparamos os termos em (10.153) e (10.160) de forma a haver concordância na expansão de Taylor até segunda ordem, isto é, restanto apenas o erro de ordem 3 e produzindo um método de ordem de precisão $p = 2$ :

b_{1} + b_{2} = 1, b_{2} c_{2} = \frac{1}{2} e a_{21} = c_{2}

(10.161)

Este sistema é formada por três equações e quatro incógnitas, pelo que admite infinitas soluções. Para construir toda a família de soluções, escolha um parâmetro $α \in (0, 1]$ e defina a partir de (10.161):

b_{1} = 1 - \frac{1}{2 α}, b_{2} = \frac{1}{2 α}, c_{2} = α, a_{21} = α

(10.162)

Portanto, obtemos o seguinte esquema genérico:

$c$	$A$


	$b$


$c_{2}$	$a_{21}$

	$b_{1}$	$b_{2}$


$α$	$α$

	$(1 - \frac{1}{2 α})$	$\frac{1}{2 α}$

0 < α \leq 1

Algumas escolhas comuns são $α = \frac{1}{2}$ , $α = \frac{2}{3}$ e $α = 1$ :

$0$
$\frac{1}{2}$	$\frac{1}{2}$

	$0$	$1$

$0$
$\frac{2}{3}$	$\frac{2}{3}$

	$\frac{1}{4}$	$\frac{3}{4}$

$0$
$1$	$1$

	$\frac{1}{2}$	$\frac{1}{2}$

Note que a tabela da direita fornece o método Euler modificado $(α = 1)$ . O esquema iterativo assume a seguinte forma: $\begin{array}{rcl} ũ_{1} & = & u^{(n)} & (10.163) \\ ũ_{2} & = & u^{(n)} + h α k_{1} & (10.164) \\ u^{(n + 1)} & = & u^{(n)} + h [(1 - \frac{1}{2 α}) k_{1} + \frac{1}{2 α} k_{2}], & (10.165) \end{array}$

onde $k_{1} = f (t^{(n)}, u^{(n)})$ e $k_{2} = f (t^{(n)} + h α, ũ_{2})$ . Ou, equivalentemente: $\begin{array}{rcl} k_{1} & = & f (t^{(n)}, u^{(n)}) & (10.166) \\ k_{2} & = & f (t^{(n)} + h α, u^{(n)} + h α k_{1}) & (10.167) \\ u^{(n + 1)} & = & u^{(n)} + h [(1 - \frac{1}{2 α}) k_{1} + \frac{1}{2 α} k_{2}], & (10.168) \end{array}$

10.7.2 Métodos de Runge-Kutta com três estágios

Os métodos de Runge-Kutta com 3 estágios podem ser descritos na forma tabular como:

$0$
$c_{2}$	$a_{21}$
$c_{3}$	$a_{31}$	$a_{32}$

	$b_{1}$	$b_{2}$	$b_{3}$

Seguindo um procedimento similar ao da Seção 10.7.1, podemos obter as condições equivalentes às condições (10.161) para um método com $ν = 3$ e ordem $p = 3$ , as quais são: $\begin{array}{rcl} b_{1} + b_{2} + b_{3} & = & 1, & (10.169) \\ b_{2} c_{2} + b_{3} c_{3} & = & \frac{1}{2}, & (10.170) \\ b_{2} c_{2}^{2} + b_{3} c_{3}^{2} & = & \frac{1}{3}, & (10.171) \\ b_{3} a_{32} c_{2} & = & \frac{1}{6}, & (10.172) \\ a_{21} & = & c_{2}, & (10.173) \\ a_{31} + a_{32} & = & c_{3} . & (10.174) \end{array}$

Assim, temos 6 condições para determinar 8 incógnitas, o que implica a existência de uma enorme família de métodos de Runge-Kutta com três estágios e ordem $p = 3$ . Se fixarmos os coeficientes $c_{2}$ e $c_{3}$ , podemos os outros de forma única: $\begin{array}{rcl} b_{1} & = & 1 - \frac{1}{2 c_{2}} - \frac{1}{2 c_{3}} + \frac{1}{3 c_{2} c_{3}} & (10.175) \\ b_{2} & = & \frac{3 c_{3} - 2}{6 c_{2} (c_{3} - c_{2})} & (10.176) \\ b_{3} & = & \frac{2 - 3 c_{2}}{6 c_{3} (c_{3} - c_{2})} & (10.177) \\ a_{21} & = & c_{2} & (10.178) \\ a_{31} & = & c_{3} - \frac{1}{6 b_{3} c_{2}} & (10.179) \\ a_{32} & = & \frac{1}{6 b_{3} c_{2}} & (10.180) \end{array}$

Alguns exemplos de métodos de Runge-Kutta de 3 estágios são o método clássico de Runge-Kutta $(c_{2} = \frac{1}{2} e c_{3} = 1)$ e o método de Nyström $(c_{2} = c_{3} = \frac{2}{3})$ :

$0$
$\frac{1}{2}$	$\frac{1}{2}$
$1$	$- 1$	$2$

	$\frac{1}{6}$	$\frac{4}{6}$	$\frac{1}{6}$

$0$
$\frac{2}{3}$	$\frac{2}{3}$
$\frac{2}{3}$	$0$	$\frac{2}{3}$

	$\frac{2}{8}$	$\frac{3}{8}$	$\frac{3}{8}$

10.7.3 Métodos de Runge-Kutta com quatro estágios

As técnicas utilizadas nas seções 10.7.1 e 10.7.2 podem ser usadas para obter métodos de quarta ordem ( $p = 4$ ) e quatro estágios ( $ν = 4$ ). As seguintes tabelas descrevem os dois esquemas mais conhecidos de Runge-Kutta quarta ordem com quatro estágios. O primeiro é denominado método de Runge-Kutta 3/8 e o segundo é chamado de método de Runge-Kutta clássico.

$0$
$\frac{1}{3}$	$\frac{1}{3}$
$\frac{2}{3}$	$- \frac{1}{3}$	$1$
$1$	$1$	$- 1$	$1$

	$\frac{1}{8}$	$\frac{3}{8}$	$\frac{3}{8}$	$\frac{1}{8}$

$0$
$\frac{1}{2}$	$\frac{1}{2}$
$\frac{1}{2}$	$0$	$\frac{1}{2}$
$1$	$0$	$0$	$1$

	$\frac{1}{6}$	$\frac{2}{6}$	$\frac{2}{6}$	$\frac{1}{6}$

O método de Runge-Kutta clássico é certamente o mais notório dos métodos de Runge-Kutta e seu esquema iterativo pode ser escrito como a seguir: $\begin{array}{rcl} k_{1} & = & h f (t^{(n)}, u^{(n)}) & (10.181) \\ k_{2} & = & h f (t^{(n)} + h ∕ 2, u^{(n)} + h k_{1} ∕ 2) & (10.182) \\ k_{3} & = & h f (t^{(n)} + h ∕ 2, u^{(n)} + h k_{2} ∕ 2) & (10.183) \\ k_{4} & = & h f (t^{(n)} + h, u^{(n)} + h k_{3}) & (10.184) \\ u^{(n + 1)} & = & u^{(n)} + \frac{k_{1} + 2 k_{2} + 2 k_{3} + k_{4}}{6} & (10.185) \end{array}$

A seguinte heurística, usando o método de Simpson para quadratura numérica, pode ajudar a compreender os estranhos coeficientes: $\begin{array}{rcl} u (t^{(n + 1)}) - u (t^{(n)}) & = & \int_{t^{(n)}}^{t^{(n + 1)}} f (t, u (s)) d s & (10.186) \\ \approx & \frac{h}{6} [f (t^{(n)}, u (t^{(n)})) + 4 f (t^{(n)} + h ∕ 2, u (t^{(n)} + h ∕ 2)) & (10.187) \\ + & f (t^{(n)} + h u (t^{(n)} + h))] & (10.188) \\ \approx & h \frac{k_{1} + 4 (\frac{k_{2} + k_{3}}{2}) + k_{4}}{6} & (10.189) \end{array}$

onde $k_{1}$ e $k_{4}$ representam os valores de $f (t, u)$ nos extremos; $k_{2}$ e $k_{3}$ são duas aproximações diferentes para a inclinação no meio do intervalo.

Exercícios resolvidos

ER 10.7.1. Construa o esquema iterativo o método clássico de Runge-Kutta três estágios cuja tabela é dada a seguir:

$0$
$\frac{1}{2}$	$\frac{1}{2}$
$1$	$- 1$	$2$

	$\frac{1}{6}$	$\frac{4}{6}$	$\frac{1}{6}$

Solução.

\begin{array}{rcl} k_{1} & = & f (t^{(n)}, u^{(n)}) & (10.190) \\ k_{2} & = & f (t^{(n)} + h ∕ 2, u^{(n)} + k_{1} ∕ 2) & (10.191) \\ k_{3} & = & f (t^{(n)} + h, u^{(n)} - k_{1} + 2 k_{2}) & (10.192) \\ u^{(n + 1)} & = & u^{(n)} + h \frac{k_{1} + 4 k_{2} + k_{4}}{6} & (10.193) \end{array}

$♢$

ER 10.7.2. Utilize o método clássico de Runge-Kutta três estágios para calcular o valor de $u (2)$ com passos $h = 1 0^{- 1}$ e $h = 1 0^{- 2}$ para o seguinte problema de valor inicial: $\begin{array}{rcl} u^{'} (t) & = & - u {(t)}^{2} + t, & (10.194) \\ u (0) & = & 0 . & (10.195) \end{array}$

Aplicando o processo iterativo obtido no Problema Resolvido 10.7.1, obtemos a seguinte rotina:

function [y]= f(t,u)
    y= t-u**2
endfunction

function [u] = RK3_classico(h,Tmax,u1)
    itmax = Tmax/h;
    u=zeros(itmax+1)
    u(1)=u1

    for i = 1:itmax
        t=(i-1)*h
        k1 = f(t, u(i))
        k2 = f(t+h/2, u(i) + h*k1/2)
        k3 = f(t+h,   u(i) + h*(2*k2-k1))

        u(i+1) = u(i) + h*(k1+4*k2+k3)/6
    end
endfunction

A qual pode ser invocada com:

->sol=RK3_classico(1e-2,2,0);sol(201)
ans  =

    1.1935760016451

-->sol=RK3_classico(1e-3,2,0);sol(2001)
ans  =

    1.1935759753635

Exercícios

E 10.7.1. Aplique o esquema de Runge-Kutta segunda ordem com dois estágios cujos coeficientes são dados na tabela a seguir

$0$
$\frac{2}{3}$	$\frac{2}{3}$

	$\frac{1}{4}$	$\frac{3}{4}$

para resolver o problema de valor inicial dado por:

\begin{array}{rcl} x^{'} (t) & = & sen (x (t)), & (10.196) \\ x (0) & = & 2 . & (10.197) \end{array}

para $t = 2$ com $h = 1 e - 1$ , $h = 1 e - 2$ e $h = 1 e - 3$ . Expresse sua resposta com oito dígitos significativos corretos.

Resposta.

2, 9677921

2, 9682284

2, 9682325

E 10.7.2. Resolva pelo método de Euler, Euler melhorado, Runge-Kutta clássico três estágios e Runge-Kutta clássico quatro estágios o problema de valor inicial tratados nos exercícios resolvidos 10.2.1 e 10.3.1 dado por: $\begin{array}{rcl} u^{'} (t) & = & - 0, 5 u (t) + 2 + t & (10.198) \\ u (0) & = & 8 & (10.199) \end{array}$

Usando os seguinte passos: $h = 1$ , $h = 1 0^{- 1}$ , $h = 1 0^{- 2}$ e $h = 1 0^{- 3}$ e compare a solução aproximada em $t = 1$ com as soluções obtidas com a solução exata dada por:

u (t) = 2 t + 8 e^{- t ∕ 2} \Rightarrow u (1) = 2 + 8 e^{- 1 ∕ 2} \approx 6, 85224527770107

(10.200)

Resposta.


Euler	6,0000000	6,7898955	6,8461635	6,8516386

$ε_{r e l}$	1,2e-01	9,1e-03	8,9e-04	8,9e-05

Euler mod,	7,0000000	6,8532949	6,8522554	6,8522454

$ε_{r e l}$	2,2e-02	1,5e-04	1,5e-06	1,5e-08

RK $_{3}$	6,8333333	6,8522321	6,8522453	6,8522453

$ε_{r e l}$	2,8e-03	1,9e-06	1,9e-09	1,8e-12

RK $_{4}$	6,8541667	6,8522454	6,8522453	6,8522453

$ε_{r e l}$	2,8e-04	1,9e-08	1,9e-12	1,3e-15

E 10.7.3. Aplique o método de Euler, o método de Euler melhorado, o método clássico de Runge-Kutta três estágios e o método clássico de Runge-Kutta quatro estágios para resolver o problema de valor inicial dado por $\begin{array}{rcl} u^{'} & = & u + t & (10.201) \\ u (0) & = & 1 & (10.202) \end{array}$

com passo $h = 1$ , $h = 1 0^{- 1}$ , $h = 1 0^{- 2}$ e $h = 1 0^{- 3}$ para obter aproximações para $u (1)$ . Compare com a solução exata dada do problema dada por $u (t) = 2 e^{t} - t - 1$ através do erro relativo e observe a ordem de precisão do método. Expresse a sua resposta com oito dígitos significativos para a solução e 2 dígitos significativos para o erro relativo.

Resposta.


Euler	2,0000000	3,1874849	3,4096277	3,4338479

$ε_{r e l}$	4,2e-01	7,2e-02	7,8e-03	7,9e-04

Euler mod	3,0000000	3,4281617	3,4364737	3,4365628

$ε_{r e l}$	1,3e-01	2,4e-03	2,6e-05	2,6e-07

RK $_{3}$	3,3333333	3,4363545	3,4365634	3,4365637

$ε_{r e l}$	3,0e-02	6,1e-05	6,5e-08	6,6e-11

RK $_{4}$	3,4166667	3,4365595	3,4365637	3,4365637

$ε_{r e l}$	5,8e-03	1,2e-06	1,3e-10	1,2e-14

Este texto é disponibilizado nos termos da licença Creative Commons Atribuição-CompartilhaIgual 3.0 Não Adaptada (CC-BY-SA 3.0). Página gerada em 19/8/2020 às 17:36:32.

Recursos

Projeto

IME - UFRGS
- Página do IME
- Página da UFRGS

UFRGS - IME - Recursos Educacionais Abertos de Matemática. Contato: reamat@ufrgs.br.

Informe erros ou edite você mesmo!