Equação de Langevin

Artur Uhlik Fröhlich e Leonardo Dasso Migotto

O objetivo deste trabalho é resolver computacionalmente a equação de Langevin utilizando o método BAOAB[LEIMKUHLER.] Serão explorados os casos de partículas individuais livres ou sujeitas a um campo potencial, estudando os efeitos da variação do coeficiente de atrito no desvio quadrático médio e na transisão de fases.

Equação de Langevin

Esta equação diferencial estocástica descreve a evolução de um sistema quando sujeito a forças do tipo determinísticas e estocásticas simultâneamente. A sua aplicação mais popular é relativa ao movimento Browniano, o movimento de uma partícula macroscópica imersa em um fluído, sujeita à força de atrito excercida pelas partículas microscópicas do fluído. Neste caso, a equação pode ser escrita como:

${\frac {d{\vec {v}}}{dt}}=-\gamma {\vec {v}}/m+\mathrm {E} (t).$

Na equação acima, $\gamma$ é o coeficiente de atrito e $\mathrm {E} (t)$ é um ruído estocástico branco, que segue o Teorema Central do Limite com média 0 e desvio padrão relacionado à temperatura, a Constante de Boltzmann, $\gamma$ e a massa da partícula. A partir desta expressão, é possível descobrir a relação do coeficiente de difusão do fluído e os valores envolvidos na equação:

$D={\frac {2k_{B}T}{\gamma m}}$

Onde $D$ é o coeficiente de difusão do meio, $k_{B}$ é a constante de Boltzmann, $T$ é a temperatura e $m$ é a massa da partícula macroscópica. Outra relação presente no livro do Frenkel [FRENKEL], desenvolvida teoricamente, é a do coeficiente de difusão e o deslocamento quadrático médio de uma partícula no meio:

$\left\langle \left(r(t)-r_{0}(t)\right)^{2}\right\rangle =2dDt$

Método BAOAB

O método numérico escolhido para realizar a integração da equação é conhecido como BAOAB, desenvolvido por Leimkuhler e Mattews utilizado para resolver equações diferenciais estocásticas.

Ele é baseado na solução exata para o momentum,

$d{\vec {p}}=-\gamma {\vec {p}}dt+{\sqrt {2\gamma mk_{B}T}}d{\vec {W}},$

e faz o uso de um método de separação das equações entre as denominadas A, B e O, respectivamente representadas:

${\vec {r}}\left(t+{\frac {\Delta t}{2}}\right)={\vec {r}}(t)+{\frac {\Delta t}{2}}{\vec {p}}\left(t+{\frac {\Delta t}{2}}\right){\frac {1}{m}}$

${\vec {p}}\left(t+{\frac {\Delta t}{2}}\right)={\vec {p}}(t)+{\frac {\Delta t}{2}}{\vec {f}}(t)$

${\vec {p'}}\left(t+{\frac {\Delta t}{2}}\right)=exp(-\gamma \Delta t){\vec {p}}\left(t+{\frac {\Delta t}{2}}\right)+{\sqrt {1-exp(-2\gamma \Delta t)}}{\sqrt {mk_{B}T}}{\vec {G}}$

O ${\vec {G}}$ aqui representa um número aleatório Gaussiano que faz o papel da força estocástica.

A equação "A" realiza meio passo no tempo da distância, a "B" realiza um meio passo para o momentum e o "O" contabiliza a contribuição estocástica equação.

Essas equações podem formar vários algoritmos de integração mas o utilizado nesse trabalho será o BAOAB:

${\vec {p}}\left(t+{\frac {\Delta t}{2}}\right)={\vec {p}}(t)+{\frac {\Delta t}{2}}{\vec {f}}(t)\quad (1)$

${\vec {r}}\left(t+{\frac {\Delta t}{2}}\right)={\vec {r}}(t)+{\frac {\Delta t}{2}}{\vec {p}}\left(t+{\frac {\Delta t}{2}}\right){\frac {1}{m}}\quad (2)$

${\vec {p'}}\left(t+{\frac {\Delta t}{2}}\right)=exp(-\gamma \Delta t){\vec {p}}\left(t+{\frac {\Delta t}{2}}\right)+{\sqrt {1-exp(-2\gamma \Delta t)}}{\sqrt {mk_{B}T}}{\vec {G}}\quad (3)$

${\vec {r}}(t+\Delta t)={\vec {r}}\left(t+{\frac {\Delta t}{2}}\right)+{\frac {\Delta t}{2}}{\vec {p'}}\left(t+{\frac {\Delta t}{2}}\right){\frac {1}{m}}\quad (4)$

${\vec {p}}(t+\Delta t)={\vec {p'}}\left(t+{\frac {\Delta t}{2}}\right)+{\frac {\Delta t}{2}}{\vec {f}}(t)\quad (5)$

É importante lembrar que entre os dois últimos passos é necessário atualizar o termo ${\vec {f}}$, já que ele pode depender de termos já atualizados como ${\vec {p}}$ ou ${\vec {r}}$.

Implementação

Simulação de uma partícula livre sob o efeito da Equação de Langevin.