Modelo de Potts 2D: mudanças entre as edições

De Física Computacional
Ir para navegação Ir para pesquisar
Sem resumo de edição
Sem resumo de edição
Linha 52: Linha 52:
O algoritmo de banho térmico, assim como o metropolis, é um algoritmo em que mudamos um spin por vez. O algoritmo segue as seguintes etapas: primeiro escolhemos um spin na rede (<math>i</math>), e independente do seu valor atual, escolhemos um novo valor para <math>s_i</math>. Esse novo valor será aceito, ou não, de acordo com os pesos de Boltzmann. Temos que no algoritmo de Banho Térmico nós atribuímos um valor <math>n</math>, entre 1 e <math>q</math>, ao spin com uma probabilidade
O algoritmo de banho térmico, assim como o metropolis, é um algoritmo em que mudamos um spin por vez. O algoritmo segue as seguintes etapas: primeiro escolhemos um spin na rede (<math>i</math>), e independente do seu valor atual, escolhemos um novo valor para <math>s_i</math>. Esse novo valor será aceito, ou não, de acordo com os pesos de Boltzmann. Temos que no algoritmo de Banho Térmico nós atribuímos um valor <math>n</math>, entre 1 e <math>q</math>, ao spin com uma probabilidade


<math>A(\mu \rightarrow \nu) = \frac{e^{-\beta E_{\nu}}}{\sum_{m=1}^q e^{-\beta E_m}}  </math>,
<math>A(\mu \rightarrow \nu) = \frac{e^{-\beta E_{\nu}}}{\sum_{m=1}^q e^{-\beta E_m}}  ,</math>


onde <math>E_n</math> é a energia do sistema quando <math>s_i = n</math> e o somatório é dado em todas energias possíveis. Pelo fato desse algoritmo permitir que o spin assuma qualquer valor ele satisfaz a condição de ergodicidade. Temos que <math>A(\mu \rightarrow \nu) = a_{\nu}</math> e <math>A(\nu \rightarrow \mu) = a_{\mu}</math>, assim, pela descrição do algoritmo temos
onde <math>E_n</math> é a energia do sistema quando <math>s_i = n</math> e o somatório é dado em todas energias possíveis. Pelo fato desse algoritmo permitir que o spin assuma qualquer valor ele satisfaz a condição de ergodicidade. Temos que <math>A(\mu \rightarrow \nu) = a_{\nu}</math> e <math>A(\nu \rightarrow \mu) = a_{\mu}</math>, assim, pela descrição do algoritmo temos


<math>\frac{A(\mu \rightarrow \nu)}{A(\nu \rightarrow \mu)} = \frac{a_{\nu}}{a_{\mu}} = \frac{e^{-\beta E_{\nu}}}{\sum_{m=1}^q e^{-\beta E_m}}  \times \frac{\sum_{m=1}^q e^{-\beta E_m}}{e^{-\beta E_{\mu}}} = e^{-\beta \Delta E}</math>
<math>\frac{A(\mu \rightarrow \nu)}{A(\nu \rightarrow \mu)} = \frac{a_{\nu}}{a_{\mu}} = \frac{e^{-\beta E_{\nu}}}{\sum_{m=1}^q e^{-\beta E_m}}  \times \frac{\sum_{m=1}^q e^{-\beta E_m}}{e^{-\beta E_{\mu}}} = e^{-\beta \Delta E}.</math>


Ou seja, o algoritmo de banho térmico respeita a condição de balanço detalhado.  
Ou seja, o algoritmo de banho térmico respeita a condição de balanço detalhado.  
Linha 64: Linha 64:
=== Aplicação ===
=== Aplicação ===


Neste trabalho, o modelo de Potts foi estudado em uma rede quadrada 2D com vizinhança de von Neumann para primeiros vizinhos e condições de contorno periódicas. A quantidade de spins no modelo é <math>N = L\times L</math> com interações ferromagnéticas e <math>J = 1</math>, favorecendo vizinhanças de spins que compartilham o mesmo valor de <math>q</math> para minimizar a energia do sistema. Em cada passo temporal, são realizadas <math>L^2</math> iterações do algoritmo de aceitação do Metropolis ou Banho Térmico. Definimos Passo de Monte Carlo (MCS) cada um desses passos temporais. Em cada passo de Monte Carlo, todos os sítios tem a probabilidade de  tentar realizar uma troca.
Neste trabalho, o modelo de Potts foi estudado em uma rede quadrada 2D com vizinhança de von Neumann para primeiros vizinhos e condições de contorno periódicas. A quantidade de spins no modelo é <math>N = L\times L</math> com interações ferromagnéticas e <math>J = 1</math>, favorecendo vizinhanças de spins que compartilham o mesmo valor de <math>q</math> para minimizar a energia do sistema. Em cada passo temporal, são realizadas <math>L^2</math> iterações do algoritmo de aceitação do Metropolis ou Banho Térmico. Definimos Passo de Monte Carlo (MCS) como cada um desses passos temporais. Em cada passo de Monte Carlo, todos os sítios tem a probabilidade de  tentar realizar uma troca.


== Resultados ==
== Resultados ==
Linha 76: Linha 76:


<source lang=C>
<source lang=C>
sitio[L*L]
int sitio[L*L];
</source>
</source>


Linha 85: Linha 85:
</source>
</source>


Uma praticidade que ganhamos com isso é que diminuímos as componentes que precisamos utilizar para diferenciar cada sítio na rede, sendo necessário apenas um número aleatório identificar o sítio. Outra sugestão é criar uma matriz de vizinhos, pois nas simulações desses sistemas geralmente precisamos saber a configuração deles. Nos vizinhos criaremos um "array" bidimensional, uma coordenada indicará o sitio que queremos e a outra indicará qual dos quatro vizinhos estamos lidando:
Uma praticidade que ganhamos com isso é que diminuímos as componentes que precisamos utilizar para diferenciar cada sítio na rede, sendo necessário apenas um número aleatório para identificar o sítio. Outra sugestão é criar uma matriz de vizinhos, pois nas simulações desses sistemas geralmente precisamos saber a configuração deles. Nos vizinhos criaremos um "array" bidimensional, uma coordenada indicará o sitio que queremos e a outra indicará qual dos quatro vizinhos estamos lidando:


<source lang=C>
<source lang=C>
Linha 91: Linha 91:
</source>
</source>


Então por lidarmos com os nossos sítios em um "array" unidimensional e condições periódicas utilizamos de fórmulas que guardarão a posição de cada vizinho nessa matriz.
Então por lidarmos com os nossos sítios em um "array" unidimensional e condições periódicas utilizamos as seguintes relações que determinam a posição de cada vizinho nessa matriz.


<source lang=C>
<source lang=C>
Linha 113: Linha 113:
=== Redução de "if"s ===
=== Redução de "if"s ===


No Modelo de Potts nós possuímos uma delta de Kronecker. Iremos subtrair <math>J</math> da nossa energia se o valor do sítio for igual ao de um vizinho. O jeito mais lógico de se pensar em simular essa parte do sistema é utilizando de condicionais "if/else", da forma
No Modelo de Potts nós possuímos uma delta de Kronecker. Iremos subtrair <math>J</math> da nossa energia se o valor do sítio for igual ao de um vizinho. A maneira mais direta de se pensar em simular essa parte do sistema é utilizando condicionais "if/else", da forma


<source lang=C>
<source lang=C>
Linha 168: Linha 168:
</source>
</source>


Assim obtemos um melhor desempenho computacional, já que não temos uma grande quantidade de "if"" no código. Também utilizamos de uma forma mais clara a delta de kronecker do nosso hamiltoniano.
Assim obtemos um melhor desempenho computacional, já que não temos uma grande quantidade de "if"" no código. Também utilizamos de forma mais clara a delta de kronecker do nosso hamiltoniano.


=== Função de Visualização ===
=== Função de Visualização ===
Linha 195: Linha 195:
</source>
</source>


Onde a função recebe o sitio utilizando um "pipe" redirecionamos a saida diretamente para a o gnuplot. Podemos também não querer plotar sempre, assim definimos uma diretiva que responde a compilação do programa
onde a função recebe o sitio. Utilizando um "pipe" redirecionamos a saida diretamente para a o gnuplot. Podemos também não querer plotar sempre, assim definimos uma diretiva que responde a compilação do programa


<source lang=C>
<source lang=C>
Linha 207: Linha 207:
</source>
</source>


Assim, se na linha de compilação tiver a opção <source lang=sh> -DVIEW </source> o programa ira, a cada 10mcs, plotar a configuração do sistema. Essa opção as vezes se torna importante para visualizar se nosso sistema esta evoluindo da forma que estamos esperando.
Assim, se na linha de compilação tiver a opção <source lang=sh> -DVIEW </source> o programa irá, a cada <math>10 mcs</math>, plotar a configuração do sistema. Essa opção as vezes se torna importante para visualizar se nosso sistema esta evoluindo da forma que estamos esperando.




Linha 218: Linha 218:
</source>
</source>


Onde ''prog.c'' é o programa que deseja utilizar. É necessário que o programa e a biblioteca com funções para MC estejam na mesma pasta no momento da compilação. E execute da seguinte maneira
Onde ''prog.c'' é o programa que deseja utilizar. É necessário que o programa e a biblioteca com funções para simulações de Monte Carlo estejam no mesmo diretório no momento da compilação. Execute da seguinte maneira


<source lang=sh>
<source lang=sh>
Linha 236: Linha 236:
</source>
</source>


Entretanto com a visualização no gnuplot o programa pode demorar mais para executar, então é recomendado diminuir os tempos de transiente (''TRAN'') e medidas (''TMAX'').
Entretanto com a visualização no gnuplot o programa demora mais para executar, então é recomendado diminuir os tempos de transiente (''TRAN'') e medidas (''TMAX'') caso se queira observar a simulação inteira.


O script possui alguns exemplos de execução através de diferentes estados e temperaturas do sistema. Copie o conteudo em um arquivo e então o deixar executável:
O script possui alguns exemplos de execução através de diferentes estados e temperaturas do sistema. Copie o conteudo em um arquivo e então o deixar executável:

Edição das 14h30min de 17 de maio de 2021

Modelo de Potts

O Modelo

O "modelo de Potts de Q-estados" trata de um sistema de rede com N spins interagentes , onde um spin pode assumir valores discretos . Cada spin do sistema está limitado a interagir com outros spins em sua vizinhança e a energia da interação entre dois spins e é dada pelo potencial

onde é a função delta de Kronecker e é a constante de interação entre os spins. Dessa maneira, a interação entre dois spins vizinhos contabiliza um valor de energia ao sistema apenas se . A hamiltoniana do sistema é dada pela soma entre todas as interações entre spins vizinhos:

Este modelo é tido como uma generalização natural do Modelo de Ising e para o caso ambos modelos são equivalentes a menos de uma constante:

Nesse caso, a interação entre dois spins e assume a mesma dinâmica do modelo de Ising a contribuição para a energia do sistema será

Motivações

Método de Monte Carlo

O método de Monte Carlo é aplicado ao modelo de Potts com o objetivo de gerar estados de equilíbrio para medir os observáveis do sistema. Os spins são iniciados com valores aleatórios de Q na rede e o método de Monte Carlo escolhe arbitrariamente um spin e gera um novo valor de para o spin. A partir disso, através de algum algoritmo específico, se escolhe como os estados serão gerados e quais serão aceitos ou não para o sistema transicionar, respeitando as condições de balanço detalhado e ergodicidade das cadeias markovianas. Para este trabalho, foram estudados os algoritmos de Metropolis-Hasting e o algoritmo de banho térmico.

Algorítmo de Metropolis-Hasting

O primeiro algoritmo utilizado para gerar as configurações do sistema foi o algoritmo de Metropolis-Hasting. O algoritmo escolhe repetidamente um novo estado para o sistema e aceitando ou rejeitando ele de acordo com uma probabilidade de aceitação de transitar de um estado antigo para o novo estado . O algoritmo que iremos descrever utiliza a dinâmica de inversão única de spins.

Temos que a condição de balanceamento detalhado é dada por [1]:

onde é a diferença de energia entre o novo e o antigo estado.

Vamos supor que tenhamos os estados e e que temos a relação de energias: . Então, a maior das duas chances de aceitação é , portanto iremos igualar essa probabilidade a 1. Para que seja respeitada, iremos definir o valor de como . Temos, assim, o algoritmo de Metropolis-Hasting:

Dessa forma, sempre que tivermos um estado cuja energia seja menor do que a do estado atual, iremos aceitar a transição, mas se a energia for maior, teremos uma pequena probabilidade de trocarmos de estado.

Algoritmo de Banho Térmico

O algoritmo de Metropolis-Hasting para inversão única de spins é eficaz para o modelo de Potts em baixos valores de ou temperaturas acima da temperatura crítica, entretanto para valores altos de ou baixas temperaturas o algoritmo falha em convergir o sistema rapidamente para o estado estacionário.

Considerando um caso onde e um spin que possui 4 vizinhos, se todos os vizinhos do spin possuem valores diferentes uns do outro e do próprio spin, poderá levar em média passos de Monte Carlo para sortear um novo valor de que tem a transição aceita, e dessa forma o algoritmo irá demorar mais tempo para alcançar a configuração de equilíbrio do sistema. A dificuldade de aceitar transições é maior ainda para baixas temperaturas, onde a probabilidade de transicionar para um novo estado tem um peso maior para qualquer diferente dos spins da vizinhança, dessa maneira poderá demorar 96 passos para gerar um spin que seja igual a algum spin da vizinhança e realizar a transição. Para contornar este problema podemos utilizar o algoritmo de banho térmico.

O algoritmo de banho térmico, assim como o metropolis, é um algoritmo em que mudamos um spin por vez. O algoritmo segue as seguintes etapas: primeiro escolhemos um spin na rede (), e independente do seu valor atual, escolhemos um novo valor para . Esse novo valor será aceito, ou não, de acordo com os pesos de Boltzmann. Temos que no algoritmo de Banho Térmico nós atribuímos um valor , entre 1 e , ao spin com uma probabilidade

onde é a energia do sistema quando e o somatório é dado em todas energias possíveis. Pelo fato desse algoritmo permitir que o spin assuma qualquer valor ele satisfaz a condição de ergodicidade. Temos que e , assim, pela descrição do algoritmo temos

Ou seja, o algoritmo de banho térmico respeita a condição de balanço detalhado.

Veremos que para valores pequenos de (como , que é o modelo de Ising), o algoritmo de Metropolis é mais eficiente. Porém, para valores altos de ou altas temperaturas o algoritmo de banho térmico é mais eficiente.

Aplicação

Neste trabalho, o modelo de Potts foi estudado em uma rede quadrada 2D com vizinhança de von Neumann para primeiros vizinhos e condições de contorno periódicas. A quantidade de spins no modelo é com interações ferromagnéticas e , favorecendo vizinhanças de spins que compartilham o mesmo valor de para minimizar a energia do sistema. Em cada passo temporal, são realizadas iterações do algoritmo de aceitação do Metropolis ou Banho Térmico. Definimos Passo de Monte Carlo (MCS) como cada um desses passos temporais. Em cada passo de Monte Carlo, todos os sítios tem a probabilidade de tentar realizar uma troca.

Resultados

Otimizações Computacionais

Essa seção é focada em algumas otimizações computacionais que podemos fazer ao simular sistemas de rede. Os "snippets" de códigos estão em C, mas a ideia pode ser adaptada para qualquer linguagem de programação.

Redes

Uma boa sugestão ao criar a topologia do sistema é utilizar um "array" unidimensional ao invés de bidimensional. Como estamos falando de uma rede quadrada de tamanho podemos fazer

int sitio[L*L];

no lugar de

int sitio[L][L];

Uma praticidade que ganhamos com isso é que diminuímos as componentes que precisamos utilizar para diferenciar cada sítio na rede, sendo necessário apenas um número aleatório para identificar o sítio. Outra sugestão é criar uma matriz de vizinhos, pois nas simulações desses sistemas geralmente precisamos saber a configuração deles. Nos vizinhos criaremos um "array" bidimensional, uma coordenada indicará o sitio que queremos e a outra indicará qual dos quatro vizinhos estamos lidando:

int viz[L*L][4];

Então por lidarmos com os nossos sítios em um "array" unidimensional e condições periódicas utilizamos as seguintes relações que determinam a posição de cada vizinho nessa matriz.

for(i=0; i<L*L; i++)
{
        viz[i][0] = (i-L+(L*L))%(L*L);
        viz[i][1] = (i+1)%L + (i/L)*L;
        viz[i][2] = (i+L)%(L*L);
        viz[i][3] = (i-1+L)%L + (i/L)*L;
}

Assim, se necessitarmos comparar com apenas um vizinho utilizamos o numero correspondente a ele, e se precisarmos comparar com todos podemos fazer um "loop". Lembrando que esse "array" guarda a posição dos vizinhos, se quisermos, por exemplo, associar um valor ao vizinho 3 do sítio fazemos da forma

sitio[viz[i][3]] = q;


Redução de "if"s

No Modelo de Potts nós possuímos uma delta de Kronecker. Iremos subtrair da nossa energia se o valor do sítio for igual ao de um vizinho. A maneira mais direta de se pensar em simular essa parte do sistema é utilizando condicionais "if/else", da forma

for(int j=0; j<_Q; j++)
{
        for(int k=0; k<4; k++)
        {
                E1 = 0;
                if(j==spin[neigh[site][k]])
                {
                        E1-=J;
                }
        }
                 E2 += exp(-E1/(KB*_TEMP));
}

Entretanto, Para valores grandes de esses condicionais "if"s resultarão em um alto custo computacional. Como temos uma delta de Kronecker no sistema podemos utilizar dela dentro do programa. Primeiro definimos o array bidimensional como a delta de kronecker da forma

for(i=0; i<_Q; i++)
{
        for(j=0; j<_Q; j++)
        {
                if(i==j)
                {
                        kronecker[i][j] = 1;
                }
                else
                {
                        kronecker[i][j] = 0;
                }
        }
}

Note que utilizamos "if/else", porém essa parte do código executa apenas uma vez. E então usamos esse "array" para otimizar as comparações:

for(j=0; j<_Q; j++)
{
        E1 = 0;

        for(k=0; k<4; k++)
        {
                E1-=J*kronecker[j][spin[neigh[site][k]]];
        }

        E2 += exp(-E1/(KB*_TEMP));
}

Assim obtemos um melhor desempenho computacional, já que não temos uma grande quantidade de "if"" no código. Também utilizamos de forma mais clara a delta de kronecker do nosso hamiltoniano.

Função de Visualização

Podemos estar interessados em ver a evolução do nosso sistema. Podemos plotar diretamente do programa sem gerar os dados e plotar depois. Definimos essa função da forma

void visualize(int *spin)
{
        int i;
        printf("unset xtics\nunset ytics\n");
        printf("set size square\n");
        printf("set xrange [0:%d]\nset yrange [0:%d]\n",L-1,L-1);
        printf("pl '-' matrix w image\n");
        for(i=L2-1; i>=0; i--)
        {
                printf("%d ", spin[i]);
                if(i%L == 0)
                {
                        printf("\n");
                }
        }

        printf("e\n\n");
}

onde a função recebe o sitio. Utilizando um "pipe" redirecionamos a saida diretamente para a o gnuplot. Podemos também não querer plotar sempre, assim definimos uma diretiva que responde a compilação do programa

#ifdef VIEW     
        if(mcs % 10 == 0)
        {
                printf("set title 'T = %d MCS'\n", mcs);
                visualize(spin);
        }
#endif

Assim, se na linha de compilação tiver a opção

 -DVIEW

o programa irá, a cada , plotar a configuração do sistema. Essa opção as vezes se torna importante para visualizar se nosso sistema esta evoluindo da forma que estamos esperando.


Programas Utilizados

Programas na linguagem C. Para utilizar os programas, abra o terminal e compile da forma

$ gcc prog.c -lm

Onde prog.c é o programa que deseja utilizar. É necessário que o programa e a biblioteca com funções para simulações de Monte Carlo estejam no mesmo diretório no momento da compilação. Execute da seguinte maneira

$ ./a.out Q TEMP

onde o segundo termo é o estado Q do sistema e o terceiro é a temperatura do sistema. Os programas possuem uma diretiva de compilação para visualização do sistema ao decorrer da execução. Para utilizar é necessário ter o gnuplot [2] instalado e compilar da forma

$ gcc -DVIEW prog.c -lm

e então executar da maneira

$ ./a.out Q TEMP | gnuplot

Entretanto com a visualização no gnuplot o programa demora mais para executar, então é recomendado diminuir os tempos de transiente (TRAN) e medidas (TMAX) caso se queira observar a simulação inteira.

O script possui alguns exemplos de execução através de diferentes estados e temperaturas do sistema. Copie o conteudo em um arquivo e então o deixar executável:

$ chmod +x script.sh

e então com o programa compilado pode executar da forma

$ ./script.sh


Potts Metropolis

Potts Banho Térmico

Biblioteca com funções em C para simulações de Monte Carlo

Script em bash para gerar os dados

Referências

  1. M. E. J. Newman, G. T. Barkema, "Monte Carlo Methods in Statistical Physics". Oxford University Press Inc., New York, 1999.
  2. https://fiscomp.if.ufrgs.br/index.php/Gnuplot