Termostato de Andersen

INTRODUÇÃO

A dinâmica molecular é uma técnica que naturalmente simula sistemas clássicos compostos por N partículas interagindo dentro de um volume V. Nesse contexto, as posições das partículas são atualizadas com base no potencial de interação escolhido. Sob a suposição de ergodicidade — ou seja, que as médias temporais equivalem às médias de ensemble —, as simulações resultam em amostragens do ensemble microcanônico (NVE). Nesse ensemble, o número de partículas N, o volume V, e a energia total E permanecem constantes (aproximadamente).

Ao colocar um sistema em contato com um reservatório térmico a uma temperatura T, mudamos do ensemble microcanônico (NVE), onde a energia é mantida constante, para o ensemble canônico (NVT), no qual a temperatura do sistema é constante. Nesse novo ensemble, a distribuição de probabilidade das velocidades das partículas segue a forma da distribuição de Maxwell-Boltzmann, uma distribuição gaussiana associada à temperatura T. Um dos métodos mais simples para realizar uma amostragem correta do ensemble canônico é o termostato de Andersen. Neste estudo, focaremos na análise desse termostato, explorando sua implementação, características e aplicação na simulação de sistemas termodinâmicos.

FUNDAMENTO TEÓRICO E DETALHES DA IMPLEMENTAÇÃO

Temperatura

A definição de temperatura em um sistema clássico em função das velocidades das partículas pode ser obtida pela equipartição de energia:

$\left\langle {\frac {1}{2}}mv_{\alpha }^{2}\right\rangle ={\frac {1}{2}}k_{B}T$

Temos então:

$T(t)={\frac {1}{k_{B}N_{f}}}\sum _{i=1}^{N_{f}}m_{i}v_{i}^{2}$

onde $N_{f}$ é o número de graus de liberdade do sistema. No ensemble microcanônico, como o momento é conservado, temos $N_{f}=3N-3$, mas no ensemble canônico, e portanto no termostato de Andersen, o momento não é conservado e utilizamos $N_{f}=3N$.

Condições de Contorno Periódicas

As condições de contorno periódicas são fundamentais para garantir que simulações computacionais de sistemas físicos representem com precisão o comportamento de sistemas grandes e infinitos. Elas ajudam a evitar efeitos de borda que não são representativos do comportamento real de partículas em um espaço muito grande, permitindo simulações mais realistas e precisas.

É utilizada a convenção da imagem mínima, que calcula a menor distância entre as partículas, sendo essa sempre menor ou igual a $L/2$.

Integração de Velocity-Verlet

Escolhido principalmente por não ser de ORDEM N² como o RK2, por exemplo. Precisamos calcular apenas 1 vez a força em cada passo do algoritmo. Assim, é o melhor "custo-benefício" dentre os possíveis algoritmos, tendo uma precisão razoável.

${\displaystyle }$

Condições Iniciais

Em todas as simulações, foram escolhidas as posições iniciais das partículas em uma grade cúbica igualmente espaçadas, para evitar que as partículas iniciem muito próximas, causando erros devido ao potencial entre as partículas, que no caso do Lennard-Jones cresce muito rápido quando as partículas se aproximam demais.

As velocidades foram escolhidas a partir de uma distribuição uniforme no intervalo $[-0.5,0.5]$, sendo reescaladas para a temperatura do reservatório, apesar de não fazer diferença para a convergência do algoritmo.

TERMOSTATO DE ANDERSEN

Para simular o contato do sistema com um reservatório térmico, a cada passo de tempo Δt é realizado um procedimento de Monte Carlo. Nesse processo, N partículas são selecionadas aleatoriamente, uma de cada vez, e suas velocidades são atualizadas. Essas novas velocidades são sorteadas a partir da distribuição gaussiana centrada em zero com $\sigma ={\sqrt {T}}$. Essa abordagem introduz um elemento de estocasticidade ao modelo, permitindo que a temperatura média do sistema oscile em torno da temperatura T do reservatório.

Na prática, a interação entre o sistema e o reservatório é controlada pela frequência de colisão $\nu$. Para cada partícula, é gerado um número aleatório $r$ no intervalo [0,1]. Caso $r<\nu *\Delta t$ a velocidade da partícula é atualizada conforme descrito.

Algoritmo

Podemos descrever o algoritmo em 4 passos:

1. Inicia-se com um conjunto de posições e velocidades.

2. Integra-se as equações do movimento para um passo $\Delta t$.

3. N partículas são selecionadas para colidir com o reservatório térmico.

4. Para cada partícula selecionada, definir nova velocidade a partir da distribuição de Maxwell-Boltzmann correspondente à temperatura T do reservatório.

RESULTADOS

Para verificar se o termostato está corretamente simulando o ensemble canônico, foram escolhidas 3 temperaturas diferentes. Para cada uma delas, foi simulado n=108 partículas, com densidade p=0.8442, nu=0.5, dt=0.001 e número de passos N=50.000. São mostradas abaixo a série temporal para cada uma das temperaturas. Para checar de forma mais consistente se o sistema no equilíbrio realmente corresponde ao ensemble canônico, foi comparada a curva teórica da distribuição de Maxwell-Boltzmann para a temperatura correspondente, com um histograma das velocidades em todas as dimensões obtidas nos últimos 40.000 passos (os primeiros 10.000 passos é o transiente para o equilíbrio). Além disso, foi calculada a média e o desvio-padrão para cada uma dessas medidas.

T=1

$\mu =1.011$ $\sigma =0.067$

Figura 1: Evolução da temperatura durante a simulação com temperatura do reservatório $T=1$.

Figura 2: Comparação da distribuição do valor absoluto das velocidades dos últimos $40.000$ passos com a distribuição de Maxwell-Boltzmann para $T=1$

T=2

$\mu =2.025$ $\sigma =0.166$

Figura 3: Evolução da temperatura durante a simulação com temperatura do reservatório $T=2$.

Figura 4: Comparação da distribuição do valor absoluto das velocidades dos últimos $40.000$ passos com a distribuição de Maxwell-Boltzmann para $T=2$.

T=3

$\mu =3.027$ $\sigma =0.247$

Figura 5: Evolução da temperatura durante a simulação com temperatura do reservatório $T=3$.

Figura 6: Comparação da distribuição do valor absoluto das velocidades dos últimos $40.000$ passos com a distribuição de Maxwell-Boltzmann para $T=3$.

CONCLUSÃO

Através deste estudo, foi possível implementar e analisar o funcionamento do termostato de Andersen como uma ferramenta eficaz para simular sistemas no ensemble canônico (NVT). A introdução de estocasticidade no modelo, via atualização aleatória das velocidades das partículas com base na distribuição de Maxwell-Boltzmann, mostrou-se eficiente para reproduzir a temperatura desejada do reservatório térmico e as propriedades estatísticas associadas ao equilíbrio termodinâmico.

Os resultados obtidos para diferentes temperaturas do reservatório confirmaram a eficácia do método, com as distribuições de velocidades convergindo para a forma teórica prevista pela estatística de Maxwell-Boltzmann. Além disso, os valores médios e os desvios padrão das temperaturas simuladas estiveram de acordo com os valores esperados, validando a robustez do termostato.

No entanto, foi observado um aumento das flutuações em temperaturas mais altas, como evidenciado pelos desvios padrão maiores. Isso é consistente com o comportamento esperado para sistemas térmicos e reflete o impacto da energia térmica sobre as velocidades das partículas.

Ao longo das simulações realizadas, foi necessário utilizar um valor de frequência de colisão $\nu =0.5$. Para valores menores de $\nu$, observou-se que o sistema não convergia adequadamente para a temperatura do reservatório térmico, apresentando desvios significativos em relação ao valor esperado. Esse comportamento pode ser explicado pelo $\Delta t$ escolhido. Se fosse menor, apesar de as colisões demorarem mais para acontecer, a dinâmica também seria menor, e talvez o valor de $\nu$ pudesse ser menor.

A frequência de colisão $\nu$ controla o nível de interação entre o sistema e o reservatório térmico. Quando $\nu$ é muito baixa, o número de colisões entre o sistema e o reservatório em um intervalo de tempo é insuficiente para equilibrar a temperatura. Nesse caso, o sistema permanece por períodos mais longos em estados que não refletem o equilíbrio térmico desejado, o que prejudica a convergência para o ensemble canônico.

Por outro lado, valores muito altos de $\nu$ podem introduzir um excesso de estocasticidade no modelo, fazendo com que a interação entre as partículas percam importância. Esse regime pode ser indesejado em simulações onde as interações interpartículas são críticas para o comportamento do sistema. Além disso, o termostato de Andersen introduz descontinuidades nas velocidades e, por consequência, no momento, sendo assim, não é possível medir propriedades de transporte, como difusão e viscosidade.

Em suma, o termostato de Andersen se mostrou uma ferramenta valiosa e relativamente simples para a simulação de sistemas interagindo com um reservatório térmico. A abordagem permite a obtenção de propriedades relevantes no ensemble canônico e abre caminho para estudos mais aprofundados, como medições de outras propriedades termodinâmicas, a comparação com outros termostatos ou a aplicação a sistemas mais complexos.