Partícula em um poço de potencial infinito e o método variacional de Monte Carlo

Particle in an infinity potential well and the variational Monte Carlo method

Partícula em um poço de potencial infinito e o método variacional de Monte Carlo

Particle in an infinity potential well and the variational Monte Carlo method

Rodrigo Canabrava; Silvio A. Vitiello

Instituto de Física, Universidade Estadual de Campinas. 13083 Campinas - SP, Brasil

^{Endereço para correspondência}

RESUMO

Neste trabalho, apresentamos o método variacional de Monte Carlo, aplicado a um sistema de uma partícula confinada em um poço de potencial infinito. A amostragem da densidade de probabilidade é feita utilizando-se o algoritmo de Metropolis. O objetivo do trabalho é apresentar de forma didática o método e sua implementação em computadores.

ABSTRACT

We apply the variational Monte Carlo method to study a particle in an infinite potential well. The probability density is sampled using the Metropolis algorithm. The aim of this work is to present in a pedagogical way both the method and its computer implementation.

I Introdução

A simulação de sistemas físicos em computadores vem se tornando essencial para o estudo da física e da ciência dos materiais, pois os modelos teóricos dos sistemas da natureza, cada vez mais realistas, são também cada vez mais complexos.

Nessas simulações são avaliadas grandezas de interesse do sistema, sem que as equações que as descrevem sejam explicitamente resolvidas. A idéia desses métodos é obter respostas numéricas sujeitas apenas a incertezas de natureza estatística em situações em que não é possível determinar resultados de forma analítica.

Neste artigo, apresentamos o método variacional de Monte Carlo e o aplicamos num problema físico extremamente simples de uma partícula confinada em um poço de potencial infinito.

O método de Monte Carlo possibilita o cálculo de diversas grandezas de interesse do sistema, entre elas, a sua energia total. Para tanto, diversas configurações (posições da partícula no espaço), são utilizadas. Cada configuração possui uma probabilidade de ocorrer. Assim, a energia do sistema é dada pela média das energias associadas a cada configuração, e ponderada pela densidade de probabilidade P de esta configuração ocorrer.

Essa ponderação é feita escolhendo-se as configurações por essa densidade de probabilidade. Dentre os algoritmos que permitem fazer essa amostragem, está o algoritmo de Metropolis[1].

O algoritmo de Metropolis permite, a cada passo, obter uma nova configuração a partir da anterior, utilizando a densidade de probabilidade P. O algoritmo permite essencialmente que configurações mais prováveis tenham mais chance de ocorrer na média das grandezas de interesse.

Na próxima seção, descrevemos o modelo da partícula no poço de potencial. Na seção 3 discutimos o método variacional. Nas duas seções seguintes, apresentamos o método variacional, e o algoritmo de Metropolis, respectivamente. Em seguida, discutimos detalhes da implementação, e apresentamos os resultados obtidos.

II Partícula confinada num poço de potencial infinito

Para uma partícula em estado estacionário, isto é, em que não há variação na sua energia, a equação de Schrödinger é dada por

em que H é a Hamiltoniana, em um problema unidimensional, escrita como

Para simplificarmos, utilizamos unidades onde = 2m = 1. Para uma partícula em um poço de potencial infinito, podemos ter:

Assim, as condições de contorno satisfeitas por y(x) devem ser:

A Eq. (1) é uma equação diferencial de segundo grau homogênea e permite soluções exatas da forma

Considerando as condições de contorno, obtemos as soluções

com n inteiro e maior que zero. Utilizando este resultado na Eq. (1), obtemos os auto-valores da energia E_n:

Portanto, para o estado fundamental, temos E₀ = p²

III Método variacional

Em situações onde não é possível determinarmos a solução exata da equação de Schrödinger, o método variacional é de grande utilidade. Isso porque, como sabemos[2], ele fornece um limite superior para o valor da energia. Ao construirmos uma função de onda tentativa y_T, é importante que ela incorpore da melhor maneira possível os conhecimentos que temos do sistema.

A energia variacional de um sistema pode ser escrita como:

Se na construção de y_T utilizamos parâmetros, E_T é minimizada com respeito aos mesmos.

Para ilustrar a aplicação do método variacional ao problema da partícula no poço de potencial, vamos adotar como função tentativa, para o estado fundamental:

em que N é a constante de normalização. Na Fig. 1 apresentamos uma comparação desta função tentativa com a função de onda exata do problema.

Portanto, para esse problema, se substituírmos y_T dada pela Eq. (9) na Eq. (8), temos que E_T = 10. Este valor é bem próximo do valor exato da Eq. (7).

IV Método variacional de Monte Carlo

A integração que acabamos de realizar analiticamente sem dúvida pode ser feita numericamente. Caso estivéssemos tratando um problema de muitas partículas e que envolvessem funções não integráveis analiticamente, provavelmente a única alternativa seria empregarmos o método de Monte Carlo para realizarmos as integrações necessárias. Apesar de nosso problema não justificar a utilização deste método, iremos empregá-lo aqui para ilustrar sua aplicação.

Na aplicação do método variacional de Monte Carlo, é útil reescrevermos a Eq. (8) como

onde P(R) é uma densidade de probabilidade dada por:

e

é chamada de energia local. Obviamente, para as auto-funções da Hamiltoniana (1), E_L será uma constante.

A probabilidade de encontrarmos uma partícula numa certa posição, P(x) = |y₀(x)|², é dada, considerando-se a função variacional da Eq. (9), por

Note que o método de Metropolis não envolve o conhecimento explícito do valor da constante N. A energia local, neste caso, é dada por

A energia variacional, como escrita na Eq. (10), mostra a utilidade dos métodos de Monte Carlo na sua determinação. É necessário amostrarmos a densidade de probabilidade P(x) para obter a seqüência de configurações x₁,x₂,...,x_M e compor a média:

com as M configurações x_i.

No processo de amostragem construimos uma seqüência x_i distribuida de acordo com P(x), isto é, teremos mais pontos onde P assume seus maiores valores. Como já mencionamos, os resultados estão sujeitos a incertezas estatísticas, o que torna essencial a estimativa da variância:

e conseqüentemente do erro padrão s/.

A amostragem de P(x), no caso particular deste trabalho, pode ser feita sem dificuldades utilizando diversos algoritmos. O algoritmo de Metropolis é um dos mais empregados na literatura por ser de aplicação geral. Ele será utilizado aqui e é descrito a seguir.

V Algoritmo de Metropolis

O algoritmo de Metropolis [1] é muito poderoso porque permite amostrarmos qualquer distribuição de probabilidade P. Uma nova configuração aleatória para o sistema é proposta a partir da anterior. As probabilidades da posição atual e da nova posição proposta são utilizadas para determinar se a nova configuração é aceita. São permitidos movimentos tanto para configurações de maior como de menor probabilidade que a atual.

O algoritmo é o seguinte:

1) Escolha uma configuração arbitrária x_i onde i = 0.

2) Faça x_novo = x_i + ( – )D,

onde é um número aleatório[3] no intervalo [0,1], e D é um parâmetro do qual os resultados são independentes.

3) Se

> ¢ então x_i+1 = x_novo

senão x_i+1 = x_i;

¢ é um novo número aleatório.

4) Faça i = i + 1.

5) Se i < M vá para 2, (M é o número total de amostras).

A nova configuração obtida a cada passo pode, eventualmente, ser igual à anterior. Novas configurações são calculadas e propostas no passo 2. D é um parâmetro que determina a amplitude do deslocamento da partícula. Seu valor é ajustado para que cerca de 50% das novas configurações sejam aceitas.

No passo 3, temos os testes para aceitar ou não essa configuração. Se o sistema vai para uma nova configuração mais provável, ela é aceita. Mesmo que a probabilidade da nova configuração seja menor, ela ainda pode ser aceita, com uma probabilidade dada pelo número aleatório ¢.

Observa-se que a configuração inicial é arbitrária. Levando em conta esse fato, descartamos os primeiros passos para que o sistema se equilibre.

VI Implementação

Para realizar a simulação, escrevemos um programa que, a partir de uma configuração inicial arbitrária x₀, com 0 < x₀ < 1, amostrava a distribuição de probabilidade dada pela Eq. (13) utilizando o algoritmo de Metropolis. Com as amostras obtidas, formamos as médias da energia E_T, de acordo com a Eq. (15).

O algoritmo de Metropolis introduz correlações entre as configurações amostradas. Isso é conseqÉência direta do passo 2 do algoritmo. Uma maneira óbvia de resolver esse problema é utilizar no cálculo da energia apenas uma amostragem a cada K posições obtidas. Por exemplo, com K = 200, formaríamos a média da energia apenas com a seqÉência X₀, X₂₀₀, X₄₀₀, ....

Com isso teríamos eliminado o problema de correlações. Uma maneira mais eficiente de tratar este problema envolve o que chamamos de média por blocos. Utilizamos um valor de K pequeno, obtendo uma seqüência X₀, X_K, X_2K, X_3K, .... Calculamos as médias parciais das energias a cada grupo de N_E dessas posições. Após obtermos as energias de N_b blocos, fazemos a média final, com o desvio padrão referente a essa média com N_b valores.

Assim, por exemplo, para K = 5 e N_E = 40, obtemos X₀, X₅, X₁₀, .... E a cada 40 dessas posições, agrupamos numa média parcial E_b. A média final E_T e o seu erro padrão são calculados assumindo que cada energia de bloco (E_b) seja estatisticamente independente das demais.

Esse procedimento resolve o problema de correlação entre as configurações amostradas, com maior eficiência no tempo de uso da máquina.

VII Resultados obtidos

O primeiro passo é escolher o valor do parâmetro D, que controla a amplitude do movimento da partícula. Para determinarmos o melhor valor para o parâmetro, consideramos a fração de movimentos aceitos como medida, isto é, a razão entre o número de vezes em que uma nova configuração é aceita pelo número total de tentativas. Os resultados podem ser vistos na Tabela 1.

Observamos que quanto menor o valor de D, maior é a fração de aceitação. Quando D é pequeno, a partícula se move muito pouco a cada passo. Para valores maiores de D, a razão das probabilidades (passo 3 do algoritmo) pode ser pequena e, além disto, o algoritmo tenta mover a partícula para fora do intervalo com muita freqÉência, resultando em rejeitarmos a configuração proposta. Uma regra empírica sugere que a fração de aceitação ideal deve ser aproximadamente 0,5. Note que tanto para o maior como para o menor valor de D apresentados na tabela, os cálculos perdem muito de sua eficiência; nestas situações o número de amostras utilizado deixa de ser suficiente. Fato este, que se reflete no valor irrealista do erro padrão.

Em seguida, realizamos a implementação de média por blocos. Mostramos na Tabela 2 alguns resultados. Esses testes foram realizados sempre com 2000 passos, que eram divididos igualmente num certo número de blocos. Cada um dos testes foi realizado usando-se a mesma sequência de configurações dos outros, de forma que obtivemos o mesmo valor de energia média em todos os casos. Assim, podemos ver mais diretamente a influência do número de blocos no valor do erro.

Conforme esperávamos, os valores de desvio padrão apresentam um máximo em função do número de blocos (N_b). Como o número de passos total é sempre o mesmo, quanto maior for o número de blocos, menor será o número de passos por bloco (N_E). Por conseqüência, menor será a variância.

De fato, estaremos tratando como avaliações independentes da energia resultados que, na realidade, não são independentes. Como já mencionamos, esta é uma decorrência do algoritmo de Metropolis. No outro extremo, a medida que diminuirmos excessivamente o número de blocos, a variância volta a diminuir.

A escolha do número de blocos mais adequado nem sempre é fácil[4]. Neste problema, 10 a 20 blocos dão uma boa idéia do erro do valor calculado.

VIII Conclusão

Neste artigo, apresentamos o método variacional de Monte Carlo. Esse método é muito útil em simulações de sistemas físicos em computadores, pois a equação de Schrödinger desses sistemas muitas vezes não permitem soluções analíticas.

De forma a exemplificar o uso desse método, utilizamos um modelo muito simples, que tem solução analítica fácil, de forma que podemos comparar o resultado obtido na simulação com a solução da equação de Schrödinger.

Além da discussão teórica, também analisamos aspectos de implementação do método variacional de Monte Carlo. Exemplificamos o algoritmo de Metropolis como um dos possíveis algoritmos que fazem a amostragem de configurações neste problema. Analisamos o problema de correlação entre as configurações, e mostramos uma maneira de se resolver esse problema.

Os resultados obtidos com a nossa função variacional são muito bons, quando comparados ao resultado exato da Eq. (7). O método variacional de Monte Carlo é de fácil implementação e fornece a energia variacional em excelente acordo com o resultado analítico. Este valor está sujeito apenas a incertezas de natureza estatística, as quais podem ser sistematicamente reduzidas pelo aumento do número de amostras.

Agradecimentos

Agradecemos ao Centro Nacional de Processamento de Alto Desempenho em São Paulo (CENAPAD-SP) na Unicamp pelo uso de suas facilidades de processamento em paralelo.

Agradecemos também à Fundação de Amparo à Pesquisa no Estado de São Paulo (FAPESP) por financiar a pesquisa com bolsa de Iniciação Científica.

Silvio A. Vitiello

E-mail: [e-mail]

Recebido em 9 de julho, 2001

Manuscrito revisado recebido em 16 de dezembro, 2002

Aceito em 13 de janeiro, 2003

• Referęncias [1] N. Metropolis et al., J. Chem. Phys. 21, 1087 (1953).
• [2] D. Griffiths, Introduction to Quantum Mechanics, Prentice Hall, New Jersey, 1995.
• [3] G. Marsaglia e A. Zaman, Comp. Phys. 8, 117 (1994).
• [4] M. P. Allen e D. J. Tildesley, Computer Simulation of Liquids, Oxford, 1989.

Datas de Publicação

• Publicação nesta coleção
  20 Maio 2003
• Data do Fascículo
  2003

Histórico

• Revisado
  16 Dez 2002
• Recebido
  09 Jul 2001
• Aceito
  13 Jan 2003