Método de Runge-Kutta 2ª e 4ª ordem: mudanças entre as edições

De Física Computacional
Ir para navegação Ir para pesquisar
Jhordan (discussão | contribs)
Sem resumo de edição
Jhordan (discussão | contribs)
Sem resumo de edição
Linha 118: Linha 118:
= Runge-Kutta 4ª ordem =
= Runge-Kutta 4ª ordem =


O método de Runge-Kutta de quarta ordem segue uma ideia similar e pode ser obtido utilizando a mesma técnica. Porém agora vamos ignorar termos de ordem <math display="inline">\Delta t^{5}</math> ou superior, então será necessário lidar com uma enorme quantidade de termos, o que torna a tarefa exaustiva e repetitiva. Logo não será feito esta demonstração aqui, mas o algoritmo de Runge-Kutta de quarta ordem pode ser dado por:
* <math display="inline"> k_{1}=f\left(y_{n},t_{n}\right)\Delta t </math>
* <math display="inline"> k_{2}=f\left(y_{n}+\frac{k_{1}}{2},t_{n}+\frac{\Delta t}{2}\right)\Delta t </math>
* <math display="inline"> k_{3}=f\left(y_{n}+\frac{k_{2}}{2},t_{n}+\frac{\Delta t}{2}\right)\Delta t </math>
* <math display="inline"> k_{4}=f\left(y_{n}+k_{3},t_{n}+\Delta t\right)\Delta t </math>
E por fim, temos então que o novo valor será dado por:
<math display="block">
y_{n+1}=y_{n}+\frac{1}{6}\left(k_{1}+2k_{2}+2k_{3}+k_{4}\right)
</math>


==Exemplo ==
==Exemplo ==
Vamos resolver o mesmo exemplo anterior, porém agora utilizando o Range-Kutta de quarta ordem.
<pre>
import matplotlib.pyplot as plt            #Biblioteca para plotar gráficos
import numpy as np                        #Biblitoeca de cálculos científicos
#Taxas de variação
def fv(x,w2):      #Velocidade
  return (-w2*x)
def fx(v):
  return (v)      #Posição


#Constantes
m=1  ; k= 1.; w2= k/m
#Parâmetros
dt  = 0.00001 ; tau = 2*np.pi; tf=4*tau ; Np= int(tf/dt)
#Valores iniciais
x=[1]; v=[0]; t=[0]; E=[k*x[0]**2/2+m*v[0]**2/2]


<pre>
#Método Range-Kutta de quarta ordem
for it  in range(Np):
  #Posição
  k1 = fx(v[it])*dt
  k2 = fx(v[it]+k1/2)*dt
  k3 = fx(v[it]+k2/2)*dt
  k4 = fx(v[it]+k3)*dt
  x.append(x[it]+(k1+2*k2+2*k3+k4)/6)
  #Velocidade
  k1 = fv(x[it]    ,w2)*dt
  k2 = fv(x[it]+k1/2,w2)*dt
  k3 = fv(x[it]+k2/2,w2)*dt
  k4 = fv(x[it]+k3  ,w2)*dt
  v.append(v[it]+(k1+2*k2+2*k3+k4)/6)
  #Energia
  E.append(k*x[it+1]**2/2+m*v[it+1]**2/2)
  #Tempo
  t.append(dt+it*dt)


plt.plot(t,x)
plt.plot(t,v)
plt.plot(t,E)
#plt.plot(x,v)
</pre>
</pre>


= Principais materiais utilizados =
= Principais materiais utilizados =


#[https://www.ufrgs.br/reamat/CalculoNumerico/livro-py/pdvi-metodos_de_runge-kutta_explicitos.html  Métodos de Runge-Kutta explícitos ] (REAMAT, UFRGS)
#[https://web.mit.edu/10.001/Web/Course_Notes/Differential_Equations_Notes/node5.html  Runge-Kutta Methods ] (Michael Zeltkevic, Instituto de Tecnologia de Massachusetts)
#[https://web.mit.edu/10.001/Web/Course_Notes/Differential_Equations_Notes/node5.html  Runge-Kutta Methods ] (Michael Zeltkevic, Instituto de Tecnologia de Massachusetts)
#[https://lpsa.swarthmore.edu/NumInt/NumIntSecond.html##section13  Second Order Runge-Kutta ] (Erik Cheever, Swarthmore)
#[https://lpsa.swarthmore.edu/NumInt/NumIntSecond.html##section13  Second Order Runge-Kutta ] (Erik Cheever, Swarthmore)

Edição das 19h48min de 20 de março de 2022

Runge-Kutta 2ª ordem

No método explícito de euler tínhamos:

yn+1=yn+f(tn,yn)Δt

Sendo dydt=f(t,y). Podemos reescrever como:

yn+1=yn+ak1

Onde a=1 e k1=f(tn,yn)Δt. Agor se supormos uma solução:

yn+1=yn+ak1+bk2(1) Com o termo adicional dependendo de uma posição genérica yn+cf(tn,yn)Δt em um tempo genérico t+dΔt, isto é k2=f(tn+dΔt,yn+cf(tn,yn)Δt)Δt. Usando o fato de que k1=f(tn,yn)Δt, podemos escrever então que:

  • k1=f(tn,yn)Δt
  • k2=f(tn+dΔt,yn+ck1)Δt

Agora lembrando a expansão em série de taylor que também vimos no método explícito e Euler:

y(t+Δt)=y(t)+y(t)Δt+y(t)Δt22+n=3y(n)(t)Δtnn!

Abrindo a segunda derivada, temos:

y(t)=d2dt2y(t)=ddtf(t,y)=ft+fydydt=ft+f(t,y)fy

Substituindo então, e escrevendo apenas 𝒪(Δ3)=n=3y(n)(t)Δtnn!, temos a seguinte expansão em série de Taylor:

y(t+Δt)=y(t)+y(t)Δt+(ft+f(t,y)fy)Δt22+𝒪(Δ3)(2)

Vamos expandir k2. Uma expansão de Taylor de primeira ordem para uma função de 2 variáveis em torno de (a,b) é dado por [1]:

f(x,y)f(a,b)+fx(a,b)(xa)+fy(a,b)(yb)

Onde fa denota a derivada da função f na variável a. Para o nosso caso, temos então para uma expansão em torno de (x,y):

f(x+Δx,y+Δy)f(x,y)+fx(x,y)Δx+fy(x,y)Δy

Expandindo então k2 em torno de (tn,yn) temos:

k2[f(tn,yn)+dΔttf(tn,yn)+ck1yf(tn,yn)]Δt

Aqui podemos notar que Δt multiplica a expansão da função, então quando desprezamos os termos de segunda ordem da expansão de f(x+Δx,y+Δy), deprezamos os termos de terceira ordem de k2. Substituindo então o k2aproximado e k1 na equação 1, temos:

yn+1=yn+af(tn,yn)Δt+b[f(tn,yn)+dΔttf(tn,yn)+cf(tn,yn)Δtyf(tn,yn)]Δt

Manipulando:

yn+1=yn+(a+b)f(tn,yn)Δt+[2bdtf(tn,yn)+2bcf(tn,yn)yf(tn,yn)]Δt22(3)

Comparando a aproximação 3 com a expansão 2 temos a seguinte relação:

a+b=1bd=12bc=12 Diferentes conjuntos de valore satisfazem este sistema. O método do ponto médio é obtido se ecolhermos: d=c=12, b=1 e a=1:

  • k2=f(tn+Δt2,yn+k12)Δt

Então:

yn+1=yn+f(tn+Δt2,yn+f(tn,yn)Δt2)Δt

O método de Heun é obtido se for escolhido a=b=12 e c=d=1:

  • k1=f(tn,yn)Δt
  • k2=f(tn+Δt,yn+k1Δt)Δt

yn+1=yn+(k1+k2)12

Uma observação, é que o erro global no algoritmo de Runge-Kutta de segunda ordem é 𝒪(Δ2) e o local é 𝒪(Δ3).

Exemplo

import matplotlib.pyplot as plt            #Biblioteca para plotar gráficos
import numpy as np                         #Biblitoeca de cálculos científicos

#Taxas de variação
def fv(x,w2):      #Velocidade
  return (-w2*x)
def fx(v):
  return (v)       #Posição

#Constantes
m=1  ; k= 1.; w2= k/m
#Parâmetros
dt  = 0.0001 ; tau = 2*np.pi; tf=4*tau ; Np= int(tf/dt)
#Valores iniciais
x=[1]; v=[0]; t=[0]; E=[k*x[0]**2/2+m*v[0]**2/2]

#Método Range-Kutta de segunda ordem, no método do ponto médio
for it  in range(Np):
  #Posição
  k1 = fx(v[it])*dt
  k2 = fx(v[it]+k1/2)*dt
  x.append(x[it]+k2)
  #Velocidade
  k1 = fv(x[it],w2)*dt
  k2 = fv(x[it]+k1/2,w2)*dt
  v.append(v[it]+k2)
  #Energia
  E.append(k*x[it+1]**2/2+m*v[it+1]**2/2)
  #Tempo
  t.append(dt+it*dt)

plt.plot(t,x)
plt.plot(t,v)
plt.plot(t,E)
#plt.plot(x,v)

Runge-Kutta 4ª ordem

O método de Runge-Kutta de quarta ordem segue uma ideia similar e pode ser obtido utilizando a mesma técnica. Porém agora vamos ignorar termos de ordem Δt5 ou superior, então será necessário lidar com uma enorme quantidade de termos, o que torna a tarefa exaustiva e repetitiva. Logo não será feito esta demonstração aqui, mas o algoritmo de Runge-Kutta de quarta ordem pode ser dado por:


  • k1=f(yn,tn)Δt
  • k2=f(yn+k12,tn+Δt2)Δt
  • k3=f(yn+k22,tn+Δt2)Δt
  • k4=f(yn+k3,tn+Δt)Δt

E por fim, temos então que o novo valor será dado por:

yn+1=yn+16(k1+2k2+2k3+k4)

Exemplo

Vamos resolver o mesmo exemplo anterior, porém agora utilizando o Range-Kutta de quarta ordem.

import matplotlib.pyplot as plt            #Biblioteca para plotar gráficos
import numpy as np                         #Biblitoeca de cálculos científicos

#Taxas de variação
def fv(x,w2):      #Velocidade
  return (-w2*x)
def fx(v):
  return (v)       #Posição

#Constantes
m=1  ; k= 1.; w2= k/m
#Parâmetros
dt  = 0.00001 ; tau = 2*np.pi; tf=4*tau ; Np= int(tf/dt)
#Valores iniciais
x=[1]; v=[0]; t=[0]; E=[k*x[0]**2/2+m*v[0]**2/2]

#Método Range-Kutta de quarta ordem
for it  in range(Np):
  #Posição
  k1 = fx(v[it])*dt
  k2 = fx(v[it]+k1/2)*dt
  k3 = fx(v[it]+k2/2)*dt
  k4 = fx(v[it]+k3)*dt
  x.append(x[it]+(k1+2*k2+2*k3+k4)/6)
  #Velocidade
  k1 = fv(x[it]     ,w2)*dt
  k2 = fv(x[it]+k1/2,w2)*dt
  k3 = fv(x[it]+k2/2,w2)*dt
  k4 = fv(x[it]+k3  ,w2)*dt
  v.append(v[it]+(k1+2*k2+2*k3+k4)/6)
  #Energia
  E.append(k*x[it+1]**2/2+m*v[it+1]**2/2)
  #Tempo
  t.append(dt+it*dt)

plt.plot(t,x)
plt.plot(t,v)
plt.plot(t,E)
#plt.plot(x,v)

Principais materiais utilizados

  1. Métodos de Runge-Kutta explícitos (REAMAT, UFRGS)
  2. Runge-Kutta Methods (Michael Zeltkevic, Instituto de Tecnologia de Massachusetts)
  3. Second Order Runge-Kutta (Erik Cheever, Swarthmore)

Citações