Física Computacional - Contribuições do usuário [pt-br]

Algoritmo de Wang-Landau

2021-11-28T17:44:03Z

Williampantaleao:

==Introdução==
Ao contrário dos métodos convencionais de Monte Carlo, que geram diretamente uma distribuição canônica <math>g(E) e^{-E/k_B T}</math> a uma dada temperatura <math>T</math>, a abordagem de Wang-Landau estima a densidade de estados <math>g(E)</math> diretamente por meio de um passeio aleatório, que produz um histograma plano no espaço de energia <math>H(E)</math>.

==Amostragem de Wang-Landau==
No início da simulação, <math>g(E)</math> é desconhecido e fazemos uma estimativa inicial para ele. A abordagem mais simples é definir <math>g(E) = 1</math> para todas as energias possíveis <math>E</math>. A configuração de spin inicial para toda a rede pode ser escolhida arbitrariamente. Então, uma caminhada aleatória no espaço de energia é iniciada pela formação de estados de teste, cada um dos quais é produzido escolhendo aleatoriamente um spin e alterando aleatoriamente seu estado.

Cada vez que uma energia <math>E</math> é visitada, a entrada correspondente em <math>H(E)</math> é incrementada em 1. A estimativa de <math>g(E)</math> é então modificada por um fator multiplicativo <math>f</math>, e o valor atualizado realiza um passeio aleatório adicional no espaço de <math>E</math>.

Se <math>E_1</math> e <math>E_2</math> são energias antes e depois de um valor de spin ser alterado, a probabilidade de transição da energia <math>E_1</math> para <math>E_2</math> é

<math>p(E_1 \rightarrow E_2) = min \left( \frac{g(E_1)}{g(E_2)}, 1\right).</math>

A razão das probabilidades de transição de <math>E_1</math> para <math>E_2</math> e de <math>E_2</math> a <math>E_1</math> podem ser calculados como

<math>\frac{p(E_1 \rightarrow E_2)}{p(E_2 \rightarrow E_1)} = \frac{g(E_1)}{g(E_2)}.</math>

Logo, o algoritmo de passeio aleatório satisfaz o equilíbrio detalhado:

<math>\frac{1}{g(E_1)}p(E_1 \rightarrow E_2) = \frac{1}{g(E_2)}p(E_2 \rightarrow E_1),</math>

onde <math>1 / g (E1)</math> é a probabilidade na energia <math>E_1</math> e <math>p (E_1\rightarrow E_2)</math> é a probabilidade de transição de <math>E_1</math> para <math>E_2</math>.

====Algoritmo====
Resumindo, o passo a passo do algoritmo pode ser escrito como:

1. Seto <math>g(E) = 1</math> e um fator de modificação <math>f=e</math>;

2. Aleatoriamente, flipo um spin com probabilidade <math> p(E_1 \to E_2) = min(g(E_1)/g(E_2), 1) </math>;

3. Modifico a densidade de estados <math> g(E) \to g(E) \times f </math> e atualizo o histograma <math>H(E)</math>;

4. Continuo até o histograma estar reto, então diminuo o valor de <math>f</math> e reseto o valor de <math>H(E)</math>;

5. Repito 2-4 até <math> \ln f \approx 1 </math>.

==Aplicação ao Modelo de Ising 2D==

====Modelo de Ising====
Uma rede 2D que consiste de uma variável discreta em cada sítio que pode ser usada para representar o momento de dipolo magnético de um átomo.

Cada sítio pode ter o valor de spin <math>+1</math> ou <math>-1</math>.

O hamiltoniano pode ser calculado por

<math> \mathcal{H} = -J \sum_{\langle ij \rangle}\sigma_i \sigma_j </math>

A soma ocorre sobre todos sítios vizinhos e a interação é ferromagnética para <math>J>0</math> e antiferromagnética para <math>J<0</math>.

====Implementação====
====Resultados====

[[Arquivo:fig1.png]]
[[Arquivo:fig2.png]]
[[Arquivo:fig3.png]]

Arquivo:Fig3.png

2021-11-28T17:10:13Z

Williampantaleao:

Arquivo:Fig2.png

2021-11-28T17:09:50Z

Williampantaleao:

Arquivo:Fig1.png

2021-11-28T17:09:22Z

Williampantaleao:

Algoritmo de Wang-Landau

2021-11-28T17:06:21Z

Williampantaleao: /* Resultados */

Algoritmo de Wang-Landau

2021-11-28T17:01:44Z

Williampantaleao: Criou página com ' ==Introdução== Ao contrário dos métodos convencionais de Monte Carlo, que geram diretamente uma distribuição canônica <math>g(E) e^{-E/k_B T}</math> a uma dada tempera...'

Trabalhos 2021-1

2021-11-28T16:30:31Z

Williampantaleao:

===[[Algoritmo de Wang-Landau]] ===

===[[Motility-Induced Phase Separation(MIPS)]]===

Equação de Burgers

2021-10-04T13:18:22Z

Williampantaleao:

(EM EDIÇÃO)

'''Grupo: Eduardo Pedroso, Luis Gustavo Lang Gaiato e William Machado Pantaleão'''

Um dos maiores desafios no campo dos [https://en.wikipedia.org/wiki/Complex_system sistemas complexos] é a compreensão do fenômeno de turbulência. Simulações computacionais contribuíram bastante para o entendimento dessa área, no entanto, ainda não existe nenhuma teoria que explique com sucesso esse comportamento e permita prever outros importantes fenômenos como misturas, convecção e combustão turbulentas, com base nas equações fundamentais da dinâmica de fluidos. Isso se deve ao fato de que a equação para os fluidos mais simples (incompressíveis) já deve levar em consideração propriedades não lineares <ref name=White> F. M. White, Viscous Fluid Flow, 3rd ed. New York, NY: McGraw-Hill Professional, 2005.
</ref>. Da [https://en.wikipedia.org/wiki/Navier%E2%80%93Stokes_equations equação de Navier–Stokes]:

<math> \frac{\partial}{\partial t}v(x,t) + \nabla \cdot v(x,t) = -\nabla p(x,t) + \nu\Delta v(x,t) </math>

<math> \nabla \cdot v(x,t) = 0 </math>

Devido à incompressibilidade, a pressão <math>p</math> é definida pela [https://en.wikipedia.org/wiki/Poisson%27s_equation equação de Poisson]:

<math> \Delta p(x,t) = -\nabla \cdot v(x,t) \cdot \nabla v(x,t)) </math>

Em 1939 o cientista alemão Johannes Martinus Burgers simplificou a equação de Navier-Stokes, removendo o termo de pressão. A equação ficou conhecida como equação de Burgers. Em uma dimensão, onde <math>u</math> corresponde ao campo de velocidades e <math>\nu</math> ao coeficiente de difusão:

<math> \frac{\partial}{\partial t}u(x,t) + u(x,t)\frac{\partial}{\partial x}u(x,t) = \nu\frac{\partial^{2}}{\partial x^{2}}u(x,t) </math>

A equação de Burgers é não linear e, portanto, espera-se um comportamento similar ao da turbulência. No entanto, foi demonstrado posteriormente que a equação de Burgers homogênea, não possuí a propriedade mais importante atribuída ao fenômeno de turbulência: o comportamento caótico em relação à pequenas mudanças nas condições iniciais. Utilizando a ''transformação de Hopf-Cole'', que transforma a equação de Burgers em uma equação parabólica linear é possível observar essa característica <ref name=Hopf>Hopf, E. (1950). The partial differential equation ut + uux = μxx. Communications on Pure and Applied Mathematics, 3(3), 201–230. doi:10.1002/cpa.3160030302</ref>.

Do ponto de vista numérico, isso é bastante importante, pois permite a comparação das soluções da equação não linear numérica com o resultado analítico. Dessa forma, pode-se investigar a qualidade do método numérico utilizado.

== A transformação de Hopf-Cole ==
A transformação de Hopf-cole mapeia a solução da equação de Burgers na equação do calor <ref name=Hopf></ref> <ref name=Evans>Evans, Lawrence C. (2010). Partial differential equations. [S.l.]: Providence, R.I. : American Mathematical Society. pp. 175–176</ref> <ref name='Meyland'>
Meylan, M., 2020. Nonlinear PDE Meylan. Lecture 12. [online] Youtube.com. Disponível em: <https://youtu.be/CsnUKrLjtyQ> [Acessado em 30 Setembro de 2021] </ref>:

<math> \frac{\partial}{\partial t}\psi(x,t)=\nu \Delta \psi(x,t) </math>

Escrevendo a equação, para uma certa condição inicial <math>F(x)</math>:

<math> \frac{\partial}{\partial t}u(x,t) + u(x,t)\frac{\partial}{\partial x}u(x,t) = \nu\frac{\partial^{2}}{\partial x^{2}}u(x,t) </math>

<math> u(x,0)=F(x) </math>

Assim, reescrevendo a equação:

<math> \frac{\partial}{\partial t}u(x,t) + \frac{\partial}{\partial x}\left(\frac{(u(x,t))^{2}}{2}-\nu \frac{\partial}{\partial x}u(x,t)\right) = 0 </math>

Busca-se uma solução <math> \psi(x,t) </math> que satisfaça:

<math> \left\{\begin{array}{l}
\cfrac{\partial \psi}{\partial x}=u\\
\cfrac{\partial \psi}{\partial t}=\nu \cfrac{\partial}{\partial x}u - \cfrac{u^{2}}{2}
\end{array}\right. </math>

Como:

<math> \frac{\partial}{\partial x}\frac{\partial \psi}{\partial t} = \frac{\partial}{\partial t}\frac{\partial \psi}{\partial x} </math>

Tem-se, que:

<math> \frac{\partial\psi}{\partial t} = \nu\frac{\partial^{2}\psi}{\partial x^{2}} -\frac{1}{2}\left(\frac{\partial\psi}{\partial x}\right)^{2} </math>

Aplicando a transformação de Hopf-Cole:

<math> \psi=-2\nu\ln{(\phi)} </math>

Assim, os seguintes resultados são obtidos:

<math> \frac{\partial\psi}{\partial x} = -\frac{2\nu}{\phi}\left(\frac{\partial \phi}{\partial x}\right) </math>

<math> \frac{\partial^{2}\psi}{\partial x^{2}} = \frac{2\nu}{\phi^{2}}\left(\frac{\partial \phi}{\partial x}\right)^{2} - \frac{2\nu}{\phi}\left(\frac{\partial^{2} \phi}{\partial x^{2}}\right) </math>

<math> \frac{\partial\psi}{\partial t} = -\frac{2\nu}{\phi}\left(\frac{\partial \phi}{\partial t}\right) </math>

Dessa forma:

<math> \frac{\partial\psi}{\partial t} = \nu \frac{\partial^{2} \psi}{\partial x^{2}} - \frac{1}{2}\left(\frac{\partial \psi}{\partial x}\right)^{2} \rightarrow -\frac{2\nu}{\phi}\left(\frac{\partial \phi}{\partial t}\right) = \frac{2\nu^{2}}{\phi^{2}}\left(\frac{\partial \phi}{\partial x}\right)^{2}-\frac{2\nu^{2}}{\phi}\left(\frac{\partial^{2} \phi}{\partial x^{2}}\right)-\frac{1}{2}\left[\frac{2\nu}{\phi}\left(\frac{\partial \phi}{\partial x}\right)\right]^{2} </math>

<math> \implies \frac{\partial \phi}{\partial t} = \nu\frac{\partial^{2} \phi}{\partial x^{2}} </math>

A equação se torna a própria [https://en.wikipedia.org/wiki/Diffusion_equation equação de difusão]. É necessário, no entanto, transformar também as condições de contorno; assim:

<math> F(x)=u(x,0)=-\frac{2\nu}{\phi(x,0)}\left(\frac{\partial \phi(x,0)}{\partial x}\right) </math>

Podendo ser reescrito como:

<math> \rightarrow \frac{d}{dx}\ln{(\phi)}=-\frac{1}{2\nu}F(x)</math>

Cuja solução:

<math> \phi(x,0)=\Phi(x)=\exp{\left(-\frac{1}{2\nu}\int_{0}^{x}F(s)ds\right)}</math>

Dessa forma, é preciso resolver:

<math> \left\{\begin{array}{l}
\cfrac{\partial \phi}{\partial t}=\nu\frac{\partial^{2} \phi}{\partial x^{2}}\\
\phi(x,0)=\Phi(x)
\end{array}\right. </math>

Utilizando a [https://en.wikipedia.org/wiki/Fourier_transform transformada de Fourier]:

<math> \left\{\begin{array}{l}
\cfrac{\partial \hat{\phi}}{\partial t}=-k^{2}\nu\hat{\phi}\\
\hat{\phi}(k,0)=\hat{\Phi}(k)
\end{array}\right. </math>

Cuja solução:

<math> \hat{\phi}(k,t)=\hat{\Phi}(k)e^{-k\nu t} </math>

Aplicando o [https://en.wikipedia.org/wiki/Convolution_theorem teorema da convolução]:

<math> \phi(x,t)=\Phi(x)\mathcal{F}^{-1}\left[e^{-k^{2}\nu t}\right] = \frac{1}{2\sqrt{\pi\nu t}}\int_{-\infty}^{\infty}\Phi(y)e^{\frac{(x-y)^{2}}{4\nu t}}dy </math>

<math> \implies \phi(x,t)=\frac{1}{2\sqrt{\pi\nu t}}\int_{-\infty}^{\infty}e^{-\frac{f}{2\nu}}dy </math>

Onde:

<math> f(x,y,t)=\frac{1}{2\nu}\int_{0}^{y}F(s)ds + \frac{(x-y)^{2}}{2t} </math>

Em <math>u</math>:

<math> \rightarrow u(x,t)=-\frac{2\nu}{\phi(x,t)}\left(\frac{\partial \phi(x,t)}{\partial x}\right) </math>

<math> \implies u(x,t)=\frac{\int_{-\infty}^{\infty}\left(\frac{x-y}{t}\right)e^{-\frac{f}{2\nu}}dy}{\int_{-\infty}^{\infty}e^{-\frac{f}{2\nu}}dy} </math>

Onde:

<math> f(x,y,t)=\frac{1}{2\nu}\int_{0}^{y}F(s)ds + \frac{(x-y)^{2}}{2t} </math>

== Modelo de Deposição - Crescimento de Interfaces ==
[[Arquivo:kpzgrow.gif|275px|thumb|right|'''Figura 1''': Modelo de deposição balística (pertencente à família dos modelos KPZ). O modelo é como um tetris onde os blocos caem e se fixam no primeiro contato <ref name='Halpin'> Halpin-Healy, T., 2015. KPZ growth model - ballistic deposition (BD). [online] Youtube.com. Disponível em: <https://youtu.be/pdeswgu9rS8> [Acessado em 30 Setembro de 2021] </ref>.]]
Um exemplo de aplicação é no crescimento de interfaces por deposição, já que a equação de Burgers é equivalente a equação conhecida como, [https://en.wikipedia.org/wiki/Kardar%E2%80%93Parisi%E2%80%93Zhang_equation ''equação de Kardar-Parisi-Zhang''] (equação ''KPZ''), um modelo de crescimento de uma superfície sólida por deposição de vapor (ou erosão de material de uma superfície sólida), que mostra a evolução da altura da camada com o tempo <ref name=Reis>REIS, F. D. A. A. Depinning transitions in interface growth models. Brazilian journal of physics, v. 33, n. 3, p. 501–513, 2003</ref>.

<math>\frac{\partial}{\partial t} h(x, t)-\frac{1}{2}\left(\nabla h(x, t)\right)^{2}=\nu \frac{\partial^{2}}{\partial x^{2}} h(x, t)+F(x, t)</math>

A equação é obtida a partir da equação de [https://en.wikipedia.org/wiki/Advection advecção simples] para uma superfície <math>z=h(x,t)</math> se movimentando com velocidade <math>U(x,t)</math>:

<math>\frac{\partial}{\partial t} h(x, t) + U \cdot \nabla h(x,t)=0</math>

A velocidade é assumida como sendo proporcional ao gradiente de <math>h(x,t)</math> (superfície evolui na direção do gradiente). A difusão da superfície é descrita pelo termo de difusão.

A equação KPZ é obtida a partir da equação Burgers, no passo imediato à aplicação da transformação de Hopf-Colem.

== Solução de Onda Viajante ==
Pode-se ainda, encontrar uma solução para a equação de Burgers na forma de uma onda viajante. Assim:

<math>u(x,t)=w(\Delta x - st)\equiv w(y)</math>

Dessa forma:

<math> \frac{\partial}{\partial t} u(x,t) = -sw^{\prime} </math>

<math> \frac{\partial}{\partial x} u(x,t) = w^{\prime}
</math>

<math> \frac{\partial^{2}}{\partial x^{2}} u(x,t) = w^{\prime\prime} </math>

Substituindo na equação de Burgers:

<math> -s w^{\prime}+w w^{\prime}=\nu u^{\prime \prime}
</math>

<math> -s w^{\prime}+\left(\frac{w^{2}}{2}\right)^{\prime}=\nu
w^{\prime\prime} \rightarrow -s w+\frac{w^{2}}{2}=\nu w^{\prime}+C </math>

Impondo condições em <math>\pm\infty</math>, tais que:

<math> \left\{\begin{array}{l}
w(-\infty)=u_{E}\\
w(\infty)=u_{D}\\
w^{\prime}(\pm\infty)=0
\end{array}\right. </math>

Onde <math>u_{E}>u_{D}</math>. Dessa forma, tem-se que:

<math> -s u_{E}+\frac{u_{E}^{2}}{2}=C=-s u_{D}+\frac{u_{D}^{2}}{2} \implies s=\frac{u_{E}+u_{D}}{2}</math>

<math> C=-u_{E}u_{D}/2 </math>

Desse modo, a velocidade de choque é a mesma para o caso de viscosidade nula <math>\left(\nu=0\right)</math>. Continuando:

<math> \nu w^{\prime}=\frac{w^{2}}{2}-\frac{u_{E}+u_{D}}{2} w+\frac{u_{E} u_{D}}{2} </math>

<math> \frac{d y}{2 \nu}=\frac{d w}{w^{2}-\left(u_{E}+u_{D}\right) u+u_{E} u_{D}} </math>

<math> \implies \frac{d y}{2 \nu}=\frac{d w}{\left(w-\cfrac{u_{E}+u_{D}}{2}\right)^{2}-\cfrac{\left(u_{E}-u_{D}\right)^{2}}{4}} </math>

Integrando ambas as partes, utilizando que:

<math> \int \frac{d w}{(w-a)^{2}-b^{2}}=\frac{1}{2 b} \ln \left|\frac{w-a-b}{w-a+b}\right| </math>

<math> \implies \frac{y}{2 \nu}+C=\frac{1}{u_{E}-u_{D}} \ln \left|\frac{w-\cfrac{u_{E}+u_{D}}{2}-\cfrac{u_{E}-u_{D}}{2}}{w-\cfrac{u_{E}+u_{D}}{2}+\cfrac{u_{E}-u_{D}}{2}}\right| </math>

<math> =\frac{1}{u_{E}-u_{D}} \ln \left|\frac{w-u_{E}}{w-u_{D}}\right|=\frac{1}{u_{E}-u_{D}} \ln{\left(\frac{u_{E}-w}{w-u_{D}}\right)} </math>

Utilizando o fato de que <math>u_{E}>w>u_{D}</math>:

<math> \frac{y}{2 \nu}+C = \frac{1}{u_{E}-u_{D}} \ln{\left(\frac{u_{E}-w}{w-u_{D}}\right)} </math>

Portanto:

<math> \frac{u_{E}-w}{w-u_{D}}=\exp{\left[{\frac{y\left(u_{E}-u_{D}\right)}{2 \nu}+C}\right]} </math>

<math> u_{E}-w=w e^{\alpha}-u_{D} e^{\alpha} \rightarrow w\left(e^{\alpha}+1\right)=u_{E}+u_{D}e^{\alpha} </math>

<math> \implies w=\frac{u_{E}+u_{D} e^{\alpha}}{e^{\alpha}+1}=u_{D}+\frac{u_{E}-u_{D}}{2} \frac{2}{e^{\alpha}+1} </math>

Onde:

<math> \alpha=\frac{y\left(u_{L}-u_{R}\right)}{2 \nu}+C </math>

Multiplicando e dividindo por <math>\exp{(-\alpha/2)}</math>:

<math> \frac{2 e^{-\alpha / 2}}{e^{\alpha / 2}+e^{-\alpha / 2}}=1-\frac{e^{\alpha / 2}-e^{-\alpha / 2}}{e^{\alpha / 2}+e^{-\alpha / 2}}=1-\tanh{\left(\frac{\alpha}{2}\right)} </math>

Dessa forma:

<math> w(y)=\frac{u_{D}+u_{E}}{2}-\frac{u_{E}-u_{R}}{2} \tanh \left(\frac{y\left(u_{E}-u_{D}\right)}{4 \nu}+C\right) </math>

Portanto:

<math> u(x, t)=\frac{u_{D}+u_{E}}{2}-\frac{u_{E}-u_{D}}{2} \tanh \left(\frac{\left(\Delta x-s t\right)\left(u_{E}-u_{D}\right)}{4 \nu}\right) </math>

=== Gráfico ===

O perfil da função <math>w(y)</math> pode ser plotado para diferentes viscosidades. Sabe-se, que conforme a viscosidade se aproxima de zero (fluido invíscido) a função se aproxima de uma função degrau de argumento <math>x-st</math> (ver [[#Apêndice|Apêndice]]).

[[Arquivo:Shockwaveee.gif|1024px|thumb|center|'''Figura 2''': Os perfis da solução da equação de Burgers víscida para, uD=0, uE=1, x0=0 e ν variando de 0.4 a 0.02 em 50 passos. Nota-se, que a função se aproxima de uma função degrau quando ν->0.]]

O gráfico pode ser construido utilizando o código em [https://www.python.org/ ''Python''] abaixo.

<source lang="python">

# -*- coding: utf-8 -*-
import matplotlib.pyplot as plt
import matplotlib.cm as cm
from matplotlib.animation import FuncAnimation
from matplotlib.lines import Line2D
import numpy as np
from numpy import ndarray
from typing import Any

u, v, x0 = 0, 1, 0 # define-se os parametros iniciais u=u_D, v=u_E

def w(y:ndarray, nu:float): # funcao w(y)
arg = (y*(v-u))/(4*nu)
return((u+v)/2 - ((v-u)/2) *np.tanh(arg))

fig = plt.figure(figsize=(12.8, 3.7), dpi=80) # resolucao da gráfico
ax = plt.axes(xlim=(-2, 2), ylim=(0, 1), # características dos eixos
xlabel=r'$y$', ylabel=r'$w(y)$')
line, = ax.plot([], [], lw=3) # cria-se um obejto 'Line' de 2D
plt.title('Perfis da Equação '+r'$w(y)$'+ # titulo do grafico
' para Diferentes Viscosidades')

def init(): # funcao que inicia a animacao caso o frame=0 seja vazio
line.set_data([], [])
return line,

def animate(frame, *fargs: Any) -> tuple[Line2D]: # função que plota os gráficos
c = 0 # inicializa a variavel
y = np.linspace(-2,2,100) # cria um array de valores para y
nu_list = np.linspace(0.4,0.02,50) # cria um array de valores de viscosidade
w_arr = w(y, nu_list[frame]) # cria o array de w(y, nu)
colour = iter(cm.rainbow(np.linspace(0,2,110))) # iteravel utilizado para as cores

for k in range(frame + 1): # avanca o objeto iteravel em (frame + 1) vezes
c = next(colour)
plt.plot(y, w_arr, color=c) # faz com que as linhas sejam mantidas entre frames
ax.text(1, 0.9, r'$\nu = $'+f'{nu_list[frame]:.3f}', fontsize=15, # adiciona nu
bbox=dict(edgecolor='white', facecolor='white', alpha=1))
line.set_color(c) # define a cor da linha com base no iteravel atualizado
line.set_data(y, w_arr) # define o novo conjunto de dados
return(line,) # retorna o objeto iteravel [Linha2D]

anim = FuncAnimation(fig, animate, init_func=init, # funcao que produz a animacao
frames=50, interval=100, blit=True) # 50 frames com intervalo de 0.1s
anim.save('shockwaveee.gif', writer='imagemagick') # gera um .GIF (e renderizacao)

</source>

== Velocidade de Choque ==
Tomando a viscosidade como nula, a equação de Burgers pode ser reescrita como<ref name=Cameron>Cameron, M. Notes on Burgers's Equation. University of Maryland, 2021</ref>:

<math> \frac{\partial u(x,t)}{\partial t} + \frac{\partial}{\partial x}\frac{u^{2}}{2} = 0 </math>

Cuja solução é dada na forma e uma onda viajante:

<math> u(x,t)=V(\Delta x - st) </math>

Onde <math>V(y)</math> é uma função degrau:

<math> \left\{\begin{array}{l}
u_{E}\text{ para } y<0\\
u_{D}\text{ para } y>0
\end{array}\right. </math>

Onde <math>u_{L}>u_{R}</math>. Essa onda é chamada de [https://en.wikipedia.org/wiki/Shock_wave onda de choque], e possui velocidade de propagação <math>s</math>. Por conservação, [https://en.wikipedia.org/wiki/Rankine%E2%80%93Hugoniot_conditions#The_jump_condition The condição de salto]:

<math> \frac{d}{d t} \int_{-M}^{M} u(x, t) d x=\int_{-M}^{M}-u u_{x} d x=-\left.\frac{u^{2}}{2}\right|_{-M} ^{M}=\frac{u_{L}^{2}}{2}-\frac{u_{R}^{2}}{2} </math>

Por outro lado:

<math> \int_{-M}^{M} u(x, t) d x=(M+s t) u_{L}+(M-s t) u_{R} </math>

Portanto:

<math> \frac{d}{d t} \int_{-M}^{M} u(x, t) d x=s\left(u_{L}-u_{R}\right) </math>

Dessa forma:

<math> s=\left(\frac{u_{L}^{2}}{2}-\frac{u_{R}^{2}}{2}\right) /\left(u_{L}-u_{R}\right)=\frac{u_{L}+u_{R}}{2} </math>

O caso generalizado da velocidade de propagação de uma onda de choque é conhecido como [https://en.wikipedia.org/wiki/Rankine%E2%80%93Hugoniot_conditions condição de Rankine-Hugoniot], e pode ser obtido, desse resultado, pela aplicação das leis de conservação hiperbólicas:

<math> \frac{\partial u}{\partial t}+\frac{\partial}{\partial x} f(u)=0 </math>

Assim:

<math> s=\frac{f\left(u_{L}\right)-f\left(u_{R}\right)}{u_{L}-u_{R}} = \frac{\text{salto em }f(u)}{\text{salto em }u} </math>

A condição de Rankine-Hugoniot, onde o salto da função é definido como choque.

== Métodos Numéricos ==
Discretizando a equação de Burgers para um fluido víscido, a partir da abordagem FTCS ([https://en.wikipedia.org/wiki/FTCS_scheme ''Foward Time Central Space'']), com a derivada à direita:

<math> \cfrac{\partial}{\partial{t}}u(x,t) = \cfrac{u(x,t+\Delta t) - u(x,t)}{\Delta t} +\mathcal{O}(\Delta t^{2}) </math>

<math> \cfrac{\partial^{2}}{\partial{x^{2}}}u(x,t) = \cfrac{u(x+\Delta x, t) - 2u(x,t) + u(x-\Delta x,t)}{(\Delta t)^{2}} +\mathcal{O}(\Delta x^{3}) </math>

Se:

<math> \left\{\begin{array}{l}
u(x, t) := U_{i}^{n}\\
u(x, t+\Delta t) := U_{i}^{n+1}\\
u(x+\Delta x, t) := U_{i+1}^{n}
\end{array}\right. </math>

<math> \implies \cfrac{\partial{u}}{\partial{t}}+u\cfrac{\partial{u}}{\partial{x}}=\nu\cfrac{\partial^{2}{U}}{\partial{x^{2}}} \longrightarrow \cfrac{U_{i}^{n+1}-U_{i}^{n}}{\Delta t} + U_{i}^{n}\cfrac{U_{i}^{n}-U_{i-1}^{n}}{\Delta x} = \nu\cfrac{U_{i+1}^{n}-2U_{i}^{n}+U_{i-1}^{n}}{(\Delta x)^{2}} </math>

<math> \longrightarrow U_{i}^{n+1} = U_{i}^{n} - U_{i}^{n}\cfrac{\Delta t}{\Delta x}(U_{i}^{n}-U_{i-1}^{n})+\nu\cfrac{\Delta t}{(\Delta x)^{2}}(U_{i+1}^{n}-2U_{i}^{n}+U_{i-1}^{n}) </math>

=== Conservação ===
Um dos problemas no cálculo de soluções descontinuas para equações como a de Burgers pode ser exemplificado da seguinte forma; supondo a equação de Burgers na forma ''quasi''-linear com viscosidade nula:<ref name=Cameron></ref>

<math> \frac{\partial u}{\partial t} + u\frac{\partial u}{\partial x} = 0 </math>

Com condições iniciais:

<math> u(x, 0)=u_{0}(x)= \begin{cases}1, & x<0 \\ 0, & x \geq 0\end{cases} </math>

Pela discretização feita anteriormente:

<math> U_{i}^{n+1} = U_{i}^{n} - U_{i}^{n}\cfrac{\Delta t}{\Delta x}(U_{i}^{n}-U_{i-1}^{n})+\nu\cfrac{\Delta t}{(\Delta x)^{2}}(U_{i+1}^{n}-2U_{i}^{n}+U_{i-1}^{n}) </math>

Pelas condições iniciais, <math>U_{j}^{0}=1</math> para <math>j<0</math> e <math>U_{j}^{0}=0</math> para <math>j\geq 0</math>:

<math> U_{j}^{1}= \begin{cases}1-\frac{\Delta t}{\Delta x} 1(1-1)=1, & j<0 \\ 0-\frac{\Delta t}{\Delta x} 0\left(0-U_{j-1}^{0}\right)=0, & j \geq 0\end{cases} </math>

Como resultado, <math>U_{j}^{1}=U_{j}^{0}</math>. Dessa forma, <math>U_{j}^{n}=U_{j}^{0}</math>; isso significa que o método propaga a descontinuidade (onda de choque) com a velocidade incorreta, <math>s=0</math>; ou seja, não ocorre propagação após o choque.

==== Conservação no Caso Invíscido ====
A partir da discretização para o caso víscido, se <math>\nu=0</math>, então:

<math>U_{i}^{n+1} = U_{i}^{n} - U_{i}^{n}\cfrac{\Delta t}{\Delta x}(U_{i}^{n}-U_{i-1}^{n})</math>

Assim como para o caso anterior, esse método não é conservativo e é adequado apenas para soluções que não possuem descontinuidades. Considerando a [https://en.wikipedia.org/wiki/Conservation_law lei de conservação]:

<math>\frac{\partial u}{\partial t} + f(u)\frac{\partial u}{\partial x} = 0</math>

Discretizando:

<math>U_{i}^{n+1}=U_{i}^{n}-\frac{\Delta t}{\Delta x}\left[f\left(U_{i}^{n}\right)-f\left(U_{i-1}^{n}\right)\right]</math>

Onde a função <math>f</math> corresponde à uma função de <math>p+q+1</math> argumentos, chamada de fluxo numérico <ref name=Cameron></ref><ref name=Mikel>Landajuela, M. Burgers Equation. Basque Center for Applied Mathematics, 2011</ref> <ref name=LeVeque>LeVeque, Randall J. (1992). Numerical Methods for Conservation Laws (PDF). Boston: Birkhäuser. p. 125. ISBN 0-8176-2723-5.</ref>. Assim, para <math>p=1</math>, <math>q=0</math> e <math>F(U,V)=f(u)</math>:

<math>U_{i}^{n+1}=U_{i}^{n}-\frac{\Delta t}{\Delta x}\left[\frac{1}{2}\left(U_{i}^{n}\right)^{2}-\frac{1}{2}\left(U_{i-1}^{n}\right)^{2}\right]</math>

=== Lax-Friedrichs ===
O método de [https://en.wikipedia.org/wiki/Lax%E2%80%93Friedrichs_method Lax-Friedrichs] para sistemas não lineares possui a forma:

<math>U_{i}^{n+1}=\frac{1}{2}\left(U_{i-1}^{n}+U_{i+1}^{n}\right)-\frac{k}{2 h}\left[f\left(U_{i+1}^{n}\right)-f\left(U_{i-1}^{n}\right)\right]</math>

Para o caso da equação de Burgers invíscida:

<math>U_{i}^{n+1}=\frac{1}{2}\left(U_{i-1}^{n}+U_{i+1}^{n}\right)-\frac{\Delta t}{2 \Delta x}\left[\frac{1}{2}\left(U_{i+1}^{n}\right)^{2}-\frac{1}{2}\left(U_{i-1}^{n}\right)^{2}\right]</math>

=== Lax-Wendroff ===
O método de [https://en.wikipedia.org/wiki/Lax%E2%80%93Friedrichs_method Lax-Wendroff] é um método de segunda ordem e possui a forma <ref name=LeVeque></ref><ref name=Mikel></ref>:

<math>U_{i}^{n+1}=U_{i+1}^{n}-\frac{k}{2 h}\left(f\left(U_{i+1}^{n}\right)-f\left(U_{i-1}^{n}\right)\right)+\frac{k^{2}}{2 h^{2}}\left[A_{i+\frac{1}{2}}\left(f\left(U_{i+1}^{n}\right)-f\left(U_{i}^{n}\right)\right)-A_{i-\frac{1}{2}}\left(f\left(U_{i}^{n}\right)-f\left(U_{i-1}^{n}\right)\right)\right]</math>

Onde <math>A_{j\pm\frac{1}{2}}</math> corresponde à matriz jacobiana <math>A(u)=f^{\prime}(u)</math>, avaliada em <math>1/2(U_{i}^{n}+U_{j\pm1}^{n})</math>. Para a equação de Burgers (invíscida) <math>f^{\prime}(u)=u</math>, assim:

<math>U_{i}^{n+1}= U_{i+1}^{n}-\frac{k}{2 h}\left(\frac{1}{2}\left(U_{i+1}^{n}\right)^{2}-\frac{1}{2}\left(U_{i-1}^{n}\right)^{2}\right)+</math>

<math>\frac{k^{2}}{2 h^{2}}\left[\left(\frac{1}{2}\left(U_{i}^{n}+U_{i+1}^{n}\right)\right)\left(\frac{1}{2}\left(U_{i+1}^{n}\right)^{2}-\frac{1}{2}\left(U_{i}^{n}\right)^{2}\right)-\left(\frac{1}{2}\left(U_{i}^{n}+U_{i-1}^{n}\right)\right)\left(\frac{1}{2}\left(U_{i}^{n}\right)^{2}-\frac{1}{2}\left(U_{i-1}^{n}\right)^{2}\right)\right]</math>

=== Método parabólico ===
Considerando a equação de Burgers víscida com <math>f(u)=u^{2}/2</math>, na forma:

<math>\frac{\partial u}{\partial t} + f(u)\frac{\partial u}{\partial x} = \nu\frac{\partial^{2} u}{\partial x^{2}}</math>

Integrando de <math>x_{j-1/2}</math> até <math>x_{j+1/2}</math>:

<math>\frac{\partial u}{\partial t} + f(u)\frac{\partial u}{\partial x} = \nu\frac{\partial^{2} u}{\partial x^{2}}</math>

== Implementação ==
Foram implementados quatro métodos: para o caso víscido, FTCS não conservativo e o método parabólico; para o caso invíscido, FTCS conservativo e Lax-Friedrichs. Definindo a condição de contorno periódica <math>U_{i=0}^{n}=u_{i_{max}}^{n}</math>, e utilizando as condições iniciais:

Caso 1:

<math>u(x, t)= \begin{cases}1, & x \leq 1 \\ 0, & x>1\end{cases}</math>

Caso 2:

<math>\left\{
\begin{array}\\
u_{i}^{n} = \cfrac{-2\nu}{\phi(x,\nu)}\left(\cfrac{\partial\phi(x,\nu)}{\partial x}\right)+4 \\
\phi=e^{\left(\cfrac{-(x-4 t)^{2}}{4 \nu(t+1)}\right)}+e^{\left[\cfrac{-(x-4 t-2 \pi)^{2}}{4 \nu(t+1)}\right]}
\end{array}
\right.</math>

Derivando <math>\phi</math> e substituindo na equação:

<math>\implies u_{i}^{n=0} = -\cfrac{x}{2\nu}e^{\left(\cfrac{-x^{2}}{4\nu}\right)}-\cfrac{(x-2\pi)}{2\nu}e^{\left[\cfrac{-(x-2\pi)^{2}}{4\nu}\right]}</math>

Cuja solução analítica é

<math>\left.u_{i}^{n}=-\frac{(x-4 t)}{2 \nu(t+1)} e^{-\left(\frac{(x-4 t)^{2}}{4 \nu(t+1)}\right)}-\frac{(x-4 t-2 \pi)}{2 \nu(t+1)} e^{\left[\frac{-(x-4 t-2 \pi)^{2}}{4 \nu(t+1)}\right.}\right]</math>

=== Importando os Módulos Utilizados ===

<source lang="python">

# -*- coding: utf-8 -*-
from math import pi, exp
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm

Array2D = np.ndarray # typing...
</source>

== Solução

== Objetivos Futuros ==
Lorem ipsum dolor sit amet, consectetur adipiscing elit. Vestibulum metus.

== Apêndice ==
O limite bilateral da função não existe no ponto na qual é centrada (assim como a função [https://en.wikipedia.org/wiki/Heaviside_step_function degrau de Heavisde]). Dessa forma, pode-se calcular o os dois limites unilaterais para entender o comportamento da função, assim:

<math>\lim _{\nu \rightarrow 0^{-}}\left(\frac{u_{D}+u_{E}}{2}-\frac{1}{2}(u_{E}-u_{D}) \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)\right) \text{ para } u_{D}<u_{E}</math>

<math>\lim _{\nu \rightarrow 0^{-}}\left(\frac{u_{E}+u_{D}}{2}-\frac{1}{2}(u_{E}-u_{D}) \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)\right)</math>

<math>\lim _{\nu \rightarrow 0^{-}}\left(\frac{u_{E}+u_{D}}{2}-\frac{1}{2} \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)(u_{E}-u_{D})\right)</math>

<math>\lim _{\nu \rightarrow 0^{-}}-\frac{1}{2} \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)(u_{E}-u_{D})+\left(\frac{u_{E}+u_{D}}{2}\right)</math>

Aplicando a regra do produto:

<math>\frac{\lim _{\nu \rightarrow 0^{-}} \tanh \left(\cfrac{y(u_{E}-u_{D})}{4 \nu}\right)(-u_{E}+u_{D})}{2}+\frac{u_{E}+u_{D}}{2}</math>

<math>\lim _{\nu \rightarrow 0^{-}} \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)=-1 \text{ para } u<u_{E} </math>

<math>\rightarrow \frac{u_{E}+u_{D}}{2}+\frac{(-u_{E}+u_{D})(-1)}{2}=\frac{u_{E}-u_{D}}{2}+\frac{u_{E}+u_{D}}{2}</math>

<math>\frac{u_{E}-u_{D}}{2}+\frac{u_{E}+u_{D}}{2}=u_{E} \text{ para } u_{D}<u_{E} </math>

<math>\implies \lim _{\nu \rightarrow 0^{-}}\left(\frac{u_{D}+u_{E}}{2}-\frac{1}{2}(u_{E}-u_{D}) \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)\right)=u_{E} \text{ para } u_{D}<u_{E}</math>

Para o segundo caso:

<math>\lim _{\nu \rightarrow 0^{+}}\left(\frac{u_{E}+u_{D}}{2}-\frac{1}{2}(u_{E}-u_{D}) \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)\right)</math>

<math>\lim _{\nu \rightarrow 0^{+}}\left(\frac{u_{E}+u_{D}}{2}-\frac{1}{2} \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)(u_{E}-u_{D})\right)</math>

<math>\lim _{\nu \rightarrow 0^{+}}-\frac{1}{2} \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)(u_{E}-u_{D})+\left(\frac{u_{E}+u_{D}}{2}\right)</math>

Aplicando a regra do produto:

<math>\frac{\lim _{\nu \rightarrow 0^{+}} \tanh \left(\cfrac{y(u_{E}-u_{D})}{4 \nu}\right)(-u_{E}+u_{D})}{2}+\frac{u_{E}+u_{D}}{2}</math>

<math>\lim _{\nu \rightarrow 0^{+}} \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)=1 \text{ para } u_{D}<u_{E} </math>

<math>\rightarrow \frac{u_{E}+u_{D}}{2}+\frac{(-u_{E}+u_{D})(1)}{2}=\frac{u_{E}+u_{D}}{2}+\frac{u_{D}-u_{E}}{2}</math>

<math>\frac{u_{E}+u_{D}}{2}+\frac{u_{D}-u_{E}}{2}=u_{D} \text{ para } u_{D}<u_{E} </math>

<math>\implies \lim _{\nu \rightarrow 0^{+}}\left(\frac{u_{E}+u_{D}}{2}-\frac{1}{2}(u_{E}-u_{D}) \tanh \left(\frac{y(u_{E}-u_{D})}{4 \nu}\right)\right)=u_{D} \text{ para } u_{D}<u_{E}</math>

Ao analisar gráfico, verifica-se esse comportamento (ver [[#Gráfico|figura 2]])..

== Referências ==

<references/>