Quantização de sistemas hamiltonianos via método de diferenças finitas

Quantization of hamiltonian systems via finite differences method

ARTIGOS GERAIS

Quantização de sistemas hamiltonianos via método de diferenças finitas

Quantization of hamiltonian systems via finite differences method

G.A. Monerat¹; L.G. Ferreira Filho; E.V. Corrêa Silva; G. Oliveira-Neto; P.H.A.S. Nogueira; A.R.P. de Assumpção

Departamento de Matemática, Física e Computação, Faculdade de Tecnologia, Universidade do Estado do Rio de Janeiro, Resende, RJ, Brasil

RESUMO

Propomos a introdução do método de diferenças finitas como tópico a ser abordado na disciplina de mecânica quântica de um curso de graduação em física. O método é suficientemente simples para ser introduzido e exemplificado em seis horas de aula, permitindo a obten cão de resultados tanto qualitativamente corretos como de alta precisão quantitativa. Devido a sua grande aplicabilidade, seu entendimento por parte dos alunos de graduação em física, futuros pesquisadores, é essencial para a formação acadêmica destes. No presente trabalho, apresentamos o método em detalhes. Sua precisão é verificada calculando o espectro de energia para o oscilador harmônico e comparando-os com os resultados analí ticos conhecidos. Em seguida aplicamos o método a dois outros sistemas, o oscilador anarmônico quártico e o potencial linear. Para cada um destes sistemas, calculamos os dez ní veis de energia mais baixos, bem como seus respectivos auto-estados.

Palavras-chave: quantização de sistemas hamiltonianos, diferenças finitas.

ABSTRACT

We propose the introduction of the finite differences method as one topic to be inserted in the discipline of quantum mechanics in a physics undergraduate course. The method is simple enough to be introduced and exemplified in about six hours of class allowing both qualitatively correct and high-precision quantitative results. Due to the great applicability of the method, it is essential that the undergraduate students, future researchers, learn it. In the present work, we show the method in detail and verify its precision initially by calculating the energy spectrum of the harmonic oscillator and comparing it to its well-known analytical results. Then we apply it to two other systems, the anharmonic oscillator and the one with linear potential. For each of those systems We compute, for both systems, the ten lowest energy eigenvalues are calculated, as well as their corresponding eigenfunctions.

Keywords: quantization of hamiltonian systems, finite differences.

1. Introdução

Ao final do curso de graduação em física, o discente depara-se com a disciplina de mecânica quântica, na qual aprende que a dinâmica quântica é governada pela equação de Schrödinger [1]; ou seja, que a solução Ψ(r,t) desta equação, denominada função de onda, descreve o comportamento do sistema em questão, contendo toda informação sobre este.

Por outro lado, soluções analíticas exatas da equação de Schrödinger são conhecidas para poucos sistemas físicos como, por exemplo, a partícula livre, o oscilador harmônico e o átomo de hidrogênio. O discente mais atento rapidamente percebe que buscar soluções analí ticas para a equação de Schrödinger não é um trabalho trivial. Apresentam-se métodos aproximados, e dentre os mais comumente encontrados na maioria das ementas de mecânica quântica estão o método WKB [2] e as teorias de perturbação independente e dependente do tempo [3].

Além destes métodos, muitos autores têm trabalhado na busca de métodos mais precisos. Em 1971, Biswas e cols. [4] apresentaram um método para a obtenção do espectro de energia para osciladores anarmônicos com termos do tipo λx^2m. Falco e cols. [5] formularam um procedimento aproximado para calcular o espectro de energia para um oscilador anarmônico com termos cúbicos e quárticos no potencial. Em 1991, Fernandez [6] propôs um método para o cálculo do espectro de energia do estado fundamental para potenciais polinomiais. No mesmo ano, Chhajlany e Malnev [7] desenvolveram um método semi-analítico e não-perturbativo para o cálculo do espectro de energia, tanto para o estado fundamental, quanto para os estados excitados de osciladores anarmônicos. Em 2006, Monerat e cols. [8, 9] aplicaram o método proposto inicialmente por Chhajlany e Malnev para a quantização de modelos cosmológicos homogêneos e isotrópicos com radiação e constante cosmológica negativa.

O espectro de energia de um sistema quântico também pode ser obtido através de métodos numéricos de diferenças finitas [10]. Recentemente utilizamos o método de diferenças finitas no esquema de Crank-Nicolson[11] para obter a quantização de modelos cosmológicos homogêneos e isotrópicos com diferentes fontes de matéria [12–15]. Métodos espectrais [16] também têm se mostrado importantes ferramentas na determinação do espectro de energia de sistemas físicos com soluções analíticas desconhecidas.

No presente trabalho, faremos uso do método de diferenças finitas [10] para obter o espectro de energia e suas respectivas autofunções. Primeiramente o aplicamos ao oscilador harmônico unidimensional e comparamos os resultados obtidos com aqueles já bem conhecidos e determinados analiticamente, E_α = (α+1/2)ђω, tal que α ∈ {0,1,2, ... }. Em seguida, aplicamos o método a dois outros sistemas: uma partícula submetida a um potencial efetivo linear unidimensional e um oscilador anarmônico quártico.

Na seção 2, elaboramos algumas considerações sobre a equação de Schrödinger. Na seção 3 o método de diferenças finitas é apresentado em detalhes, e na seção 4 o aplicamos ao caso de um oscilador harmônico, calculando seus 10 níveis mais baixos de energia. Em seguida comparamos o resultado numérico para o espectro de energia com o resultado analítico, conhecido na literatura. Nas seções 5 e 6, aplicamos o método, respectivamente, aos casos de uma partícula submetida a um potencial efetivo linear unidimensional e a um oscilador anarmônico quártico. Na seção 7, apresentamos nossos comentários finais e conclusões.

2. A equação de Schrödinger

A mecânica quântica tem a função de onda Ψ(r,t), como principal objeto para descrever as propriedades de sistemas físicos. A forma da função de um dado sistema físico dependerá da função de Hamilton ou hamiltoniana do sistema. De um modo geral, a hamiltoniana de uma partí cula é da forma

em que p corresponde ao momentum canonicamente conjugado à posição da partí cula r, m é a massa da partícula e V(r) é a energia potencial desta sob a ação de um campo externo independente do tempo. Em sua versão unidimensional, temos

Elevando ao grau de operadores a posição , seu momentum canônico _x e a função de Hamilton , temos na representação de Schrödinger unidimensional a seguinte correspondência

em que fizemos ђ = 1. Assim, obtemos a equação de Schrödinger, que governa a dinâmica do sistema quântico

ou

Utilizando o método de separação de variáveis, fazemos Ψ(x,t) = φ(x) ξ(t), obtendo assim as equações

na qual introduz-se a constante de separação E. A solu cão da parte temporal é ξ(t) ∝ e^-iEt; a equa cão para a parte espacial φ(x) é da forma

denominada equação de Schrödinger independente do tempo, que equivale a uma equação de autovalores para o operador hamiltoniano . Os estados estacionários do sistema físico são representados pelas autofunções φ_α(x) do operador hamiltoniano, ou seja, funções que satisfazem

Os autovalores E_α correspondem aos ní veis de energia do sistema. Uma solução da equação de Schrödinger pode ser obtida através da superposição (combinação linear) destas soluções,

em que C_α, são constantes arbitrárias, tal que Ψ(x,t) seja normalizada.

Dependendo da forma do potencial V(x), soluções analíticas para a equação de Schrödinger não são conhecidas. Nesses casos, métodos numéricos são necessários.

3. O método de diferenças finitas

O método de diferenças finitas, essencialmente, é usado para transformar uma equação diferencial ordinária ou parcial, linear ou não, em um sistema de equações algébricas, onde as incógnitas são os valores da função em cada ponto de um conjunto discreto de pontos igualmente espaçados(denominado reticulado) obtido pela discretização da variável de integração x. O método apresenta certa ambigüidade, na medida em que, por construção, as derivadas podem ser aproximadas por diferenças finitas ascendente, descendente ou centrada, em relação ao ponto desejado. Consideremos a discretização mostrada na Fig. 1, na qual o intervalo de interesse [x₀, x_n+1] é particionado em n+1 sub-intervalos de mesma amplitude h = x_i+1– x_i. Os pontos x₀ e x_n+1 correspondem aos pontos extremos do intervalo.

Primeiramente expandimos uma função f(x) em série de Taylor em torno do ponto x_i do reticulado, até o termo em 2ª ordem, obtendo

em que f'(x) e f^"(x) indicam as derivadas ordinárias de primeira e segunda ordem, respectivamente, da função f(x) em relação à coordenada x. Considerando esta série somente até o termo de 1 ordem, obtemos

que é denominada diferenciação ascendente e é representada pelo sub-índice A. Por outro lado, também é possível escrever

o que nos conduz, após procedermos o mesmo truncamento, a

que é denominada diferenciação descendente, e é representada aqui pelo sub-índice D. Observe que o deslocamento da variável, neste caso, é x_i – x_i-1 = –h. Subtraindo a Eq. (10) da Eq. (12), obtemos a fórmula para a diferenciação centrada,

Somando a Eqs. (10) e (12) (nas quais os termos de ordem superiores aos de segunda ordem são desprezados), expressão para a derivada segunda centrada, que tem a forma

O mesmo resultado pode ser obtido se notarmos que

a substituição das Eqs. (11) e (13) na Eq. (16) conduz novamente à Eq. (15). Os erros cometidos nestas aproximações são da ordem de O(h²). Estes procedimentos podem ser facilmente estendidos para a determinação de expressões ascendentes, descendentes e centradas de derivadas de qualquer ordem [17].

Conforme mencionado na seção anterior, o comportamento de um sistema quântico é governado pela equação de Schrödinger. A busca de soluções para esta equação recai em um problema de autovalores e autovetores (8).

Aplicando a Eq. (15) na Eq. (6) para determinado ponto x_i do reticulado, obtemos,

Fazendo i ∈ {1,2,... n}, podemos aplicar a Eq. (17), desde o ponto x₁ até o ponto x_n. O resultado é um sistema algébrico de n equações com n+2 incógnitas - os valores da função φ_α nos pontos (x₀, x₁, ...,x_n+1). Entretanto, os valores de φ_α nestes pontos são determinados pelas condições de contorno, que aqui serão escolhidas como

Esta condição de contorno é compatível com a busca de autofunções localizadas de sistemas com potenciais finitos (como as autofunções correspondendo a ní veis discretos de energia de um poço de potencial, por exemplo) desde que o intervalo considerado seja suficientemente grande em rela cão à região classicamente permitida, onde as autofun cões são sensivelmente não-nulas.

Com esta informação adicional, o que obtemos finalmente é um sistema de equações algébricas com solução única. Este sistema pode ser convenientemente expresso na forma matricial, se escrevermos a função φ_α(x) como um vetor

Com isso, a Eq. (17) assume a forma

em que F(x) = [1/( mh²)]+V(x). Observamos que a matriz n×n é real e simétrica, e portanto diagonalizável; é garantida a existência de uma base de autoestados.

É claro que, numericamente, sempre utilizamos um valor finito para n, o que resulta em um número finito α_max de auto-energias e seus respectivos auto-estados. Isto porque, ao dividirmos o intervalo utilizando n pontos, gera-se uma matriz n×n, cuja equação característica é uma equação de grau n, cujas n raízes são os n auto-valores de valor mais baixo. No caso especifico de sistemas não-degenerados, temos α_max = n-1, com α variando de 0 até α_max. Quanto mais refinado for o esquema de discretização, ou seja, quanto maior for o número de pontos n do reticulado e menor o seu espaçamento h, conforme exibido na Fig. 1, mais preciso será o cálculo dos autovalores da energia. Obviamente, este refinamento está condicionado às limitações computacionais (e.g., erros de truncamento e arredondamento, quantidade de memória disponí vel).

Como uma ilustração simples, consideremos a discretização mostrada na Fig. 2, com apenas 5 pontos; o indice i varia de 1 a 3. Aplicando a Eq. (17) nos pontos 1, 2 e 3, obtemos o seguinte sistema de equações algébricas

Novamente utilizamos as condições de contorno φ_α(x₀) = φ_α(x₄) = 0, o que nos leva a um sistema linear de três equações. Este sistema pode ser escrito na forma matricial como segue,

Particularizemos a discussão para um sistema bem conhecido como, por exemplo, um oscilador harmônico, ou seja, para um potencial da forma

com k = m = 1 e integrando a equação no intervalo -10 a 10, de acordo com a Fig. 2, obtemos os três auto-valores E₀ = 0.03993600032767664, E₁ = 12.5400000000000 e E₂ = 12.5400639996723. É claro que esses resultados são extremamente ruins se comparados ao resultado analítico (veja a Tabela 1).

Este exemplo tem como finalidade exclusiva mostrar o funcionamento do método de uma forma mais didática, mesmo que em detrimento da precisão do cálculo numérico. Ele pode ser resolvido como um exercício simples de algebra linear, inclusive sem a necessidade de um computador, na medida em que o determinante da matriz conduz a uma equação característica de terceiro grau na variável E_α, solúvel em uma calculadora científica típica atual. As raizes desta equação são os três auto-valores apresentados acima.

Para se obter as auto-funções, basta substituir cada auto-valor nas Eqs. (21) e resolver o sistema linear correspondente. Por exemplo, para o auto-valor E₀, temos

O escalonamento deste sistema linear evidencia que ele tem, efetivamente, duas equações e três incógnitas e que portanto tem infinitas soluções ou, dito de outra maneira, uma solução com um parâmetro livre. Esta ambiguidade é removida quando impomos que cada auto-vetor tenha a norma euclidiana igual a 1, ou seja,

Dessa forma, a solução do sistema será φ₀(x₁) = φ₀(x₃) = 1.6 ×10^-3 e φ₀(x₂) = 1, como pode ser visto no primeiro gráfico da Fig. 3.

Repetindo este procedimento é possivel calcular os outros dois auto-vetores, correspondentes aos auto-valores E₁ e E₂, que também estão representados na Fig. 3.

Nas seções a seguir, o método exposto acima será aplicado a três sistemas físicos, utilizando agora um número muito maior de pontos no reticulado do que aquele do exemplo acima. Consideraremos novamente o (i) oscilador harmônico; em seguida, estudaremos o (ii) oscilador anarmônico quártico e, finalmente, (iii) uma partícula submetida a um potencial linear. O algoritmo para o método mencionado aqui foi implementado em linguagem FORTRAN e o código-fonte está disponível no Apêndice 1. Espera-se que o discente de graduação em física no último ano de graduação já possua conhecimentos básicos em alguma linguagem de programação; e que já tenha cursado disciplinas tais como Introdução ao Processamento de Dados, Cálculo Numérico e Métodos Computacionais.

4. O oscilador harmônico

Uma dúvida que pode surgir na mente do discente de graduação nesse momento é quanto à precisão do método. Para elucidar esta questão, vamos aplicar o método ao problema da quantização do oscilador harmônico, e com isso determinarmos numericamente o seu espectro de energia. O espectro de energia do oscilador harmônico é conhecido na literatura [1], o que permitirá compararmos o resultado numérico com o analítico e assim avaliarmos a precisão do método.

O operador hamiltoniano para o oscilador harmônico conduz à equação de Schrödinger estacionária (6) com um potencial na forma (20). De acordo com a Eq. (20) obtemos a equação matricial

em que F(x) = [1/( mh²)]+kx². Vamos considerar um reticulado com n = 1500 pontos e distância entre pontos adjacentes h = 0.013333.... Como no exemplo didático na seção anterior, usaremos aqui k = m = 1.

A partir da Tabela 1 observamos que os valores do espectro obtido através do método de diferenças finitas são muito próximos dos valores teóricos. O maior desvio encontrado é da ordem de δE ~ 10^-3. A precisão dos resultados é maior para os níveis de energia mais baixos, já que, conforme a energia aumenta, as auto-funções ficam cada vez menos localizadas, ou seja, elas ainda são não nulas nas proximidades dos pontos de "infinito" numérico, ou seja, dos pontos extremos do intervalo de trabalho. Nestes pontos, o potencial já não é mais um oscilador harmônico e sim uma barreira infinita, o que é uma das limitações do algoritmo computacional, já mencionadas anteriormente. Isso de forma alguma invalida o método, na medida em que possamos estender o intervalo de integração e refinar a discretização, até o limite da capacidade computacional. Vale enfatizar que a precisão ótima dos resultados resulta da mediania de fatores como a amplitude do intervalo de integração e o refinamento da malha, bem como de capacidade de memória capaz de lhes dar suporte.

Utilizando um método de interpolação [18] podemos estimar a dependência da energia E em relação ao parâmetro a, tal que α ∈ { 0,1,2,...} enumera os níveis de energia. Aplicando o método de interpolação, encontramos

Para o estado fundamental o desvio é da ordem de 10^-6; enquanto que para o primeiro estado excitado, por exemplo, o desvio é da ordem de 10^-5.

Os autoestados associados aos cinco autovalores mais baixos da energia mostradas na Tabela 1 são mostrados na Fig. 4.

Se compararmos o resultado dos três primeiros estados do oscilador harmônico exibidos na Fig. 4 com os mesmos autoestados calculados para o reticulado de cinco pontos (Fig. 2), observamos claramente uma grande melhoria na precisão destes.

Agora que a eficácia do método foi verificada, vamos considerar um sistema mais complicado, que não possui solução analítica conhecida.

5. O oscilador anarmônico quártico

Estudaremos a seguir o oscilador anarmônico quártico, caracterizado por um potencial na forma

em que L é um parâmetro do sistema. A aplicação do método de diferenças finitas a partir da Eq. (26) fornece os autoestados e as autovalores da energia com igual precisão. A malha para este caso terá 1500 pontos no intervalo [-10, 10], com h = 0.1 . A Fig. 5 exibe a variação do espectro de energia em relação à energia do respectivo estado fundamental com Λ = 0. Observamos também que à medida que Λ aumenta, os valores dos ní veis de energia também crescem.

Ao leitor interessado deixamos a construção das autofunções para os primeiros 3 níveis. Isso pode ser feito utilizando-se a rotina em FORTRAN disponível no apêndice deste trabalho.

6. Partícula sob a ação de um potencial linear

Agora aplicaremos o método a uma partícula submetida a um potencial linear (V(x) = x). Fazendo uso da expressão (20) encontramos o espectro de energia do sistema. Neste caso a malha construída utiliza 1500 pontos, no intervalo fechado [0, 20 ], com um h = 0.01. Na Tabela 2 exibimos os 10 primeiros níveis de energia do sistema.

Os cinco primeiros autoestados são mostrados na Fig. 6. Repare que as autofunções encontradas numericamente correspondem às conhecidas funções de Airy [19].

O resultado numérico obtido para o primeiro autoestado φ₀ pode ser interpolado [18] pela função φ₀(x) = xe^-xAi(x), em que Ai(x) é a função de Airy.

7. Conclusão e comentários finais

Observamos que o método de diferenças finitas, devido à sua simplicidade e ampla aplicabilidade, é ideal para a quantização de sistemas hamiltonianos cujas soluções analíticas não sejam conhecidas. Destacamos ainda a precisão do método, o que o torna uma ferramenta indispensável não somente para os discentes de graduação em física, mas também para pós-graduandos.

Neste artigo apresentamos o método mencionado acima em detalhes e posteriormente o aplicamos a três sistemas físicos, a saber: o oscilador harmônico simples, cuja solução analítica é conhecida; o oscilador anarmônico quártico; e uma partícula submetida a um potencial linear. Em todos os casos, o método aplicado forneceu o espectro de energia e os autoestados para os sistemas estudados. A dinâmica quântica desses sistemas será descrita por pacotes de onda construídos pela superposição dos respectivos autoestados para uma dada condição inicial. Deixamos isso como exercícios ao leitor interessado. No caso do oscilador harmônico simples, o espectro de energia obtido pelo método numérico concorda com o resultado analítido com uma precisão muito boa. O resultado para espectro de energia do segundo caso, oscilador quártico, mostrou apresentar um bom resultado; tendo em vista que para Λ = 0 reobtemos o estado fundamental do oscilador harmônico. No caso do sistema submetido a um potencial linear do tipo V(x) = x, os autoestados passam a ser descritos por funções de Airy, o que está de acordo com a literatura [1].

Agradecimentos

G.A. Monerat (Bolsista do CNPq, Brasil), E.V. Corrêa Silva (Bolsista do CNPq, Brasil), e G. de Oliveira-Neto (Bolsista do CNPq, Brasil) agradecem ao CNPq e a FAPERJ pelo apoio financeiro. Este trabalho foi parcialmente realizado no Laboratório de Computação Avançada (LCA) do Departamento de Matemática, Física e Computação da Faculdade de Tecnologia da UERJ.

Recebido em 28/5/2009; Aceito em 21/9/2009; Publicado em 26/3/2010

Apêndice 1

Código-fonte em Fortran para o calculo do espectro de energia e das autofunções

Reproduzimos a seguir o código-fonte em linguagem FORTRAN, com comentários, que implementa o algoritmo de diferenças finitas.

• [1] L.D. Landau and E.M. Lifshitz, Quantum Mechanics: Non-relativistic Theory (Pergamon, Oxford, 1977), 3rd. ed.
• [2] David J. Griffiths, Introduction Quantum Mechanics (Prentice Hall, Upper Saddle River, 1995).
• [3] W. Greiner, Quantum Mechanics - An Introduction (Springer, New York, 2001) 4th ed.
• [4] S.N. Biswas, K. Datta, R.P. Saxena, P.K. Srivastava and V.S. Varma. J. Math. Phys. 14, (1973).
• [5] L. De Falco, R. Mignani and R. Scipioni, Europhys. Lett. 36, 81 (1996).
• [6] F.M. Fernandez, Physical Review A 44, 3336 (1991).
• [7] S.C. Chhajlany, D.A. Letov and V.N. Malnev, J. Phys. A: Math. Gen. 24, 2731 (1991).
• [8] G.A. Monerat, E.V. Corrêa Silva, G. Oliveira-Neto, L.G. Ferreira Filho, and N.A. Lemos, Phys. Rev. D 73, 044022 (2006).
• [9] N.A. Lemos, G.A. Monerat, E.V. Corrêa Silva, G. de Oliveira Neto, L.G. Ferreira Filho. Physical Review D 75, 068504 (2007).
• [10] J.W. Thomas, Numerical Partial Differential Equations: Finite Difference Methods (Texts in Applied Mathematics) (Springer, New York, 1998).
• [11] J. Crank and P. Nicolson, Proc. Cambridge Philos. Soc. 43, 50 (1947).
• [12] G.A. Monerat, E.V. Corrêa Silva, G. de Oliveira Neto, L.G. Ferreira Filho, and N.A. Lemos, Brazilian Journal of Physics 35, 1106 (2005).
• [13] J. Acacio de Barros, E.V. Corrêa Silva, G.A. Monerat, G. Oliveira-Neto, L.G. Ferreira Filho and P. Romildo Jr., Phys. Rev. D 75, 104004 (2007).
• [14] G.A. Monerat, G. Oliveira-Neto, E.V. Corrêa Silva, L.G. Ferreira Filho, P. Romildo Jr., J.C. Fabris, R. Fracalossi, S.V.B. Gonçalves and F.G. Alvarenga, Phys. Rev. D 76, 024017 (2007).
• [15] L.G. Ferreira Filho, J. Acacio de Barros, E.V. Corrêa Silva, G.A. Monerat, G. de Oliveira-Neto and P.R. Pires JR, International Journal of Modern Physics E 16, 3014 (2007).
• [16] P. Pedram, M. Mirzaei, S. Jalalzadeh and S.S. Gousheh, Gen Relativ Gravit 40, 16631681 (2008).
• [17] S.C. Chapra and R.P. Canale, Numerical Methods for Engineers (McGraw-Hill, New York, 1988), 2nd ed.
• [18] D.F. Rogers and J.A. Adams, Mathematical for Computer Graphics (McGraw-Hill, New York, 1990).
• [19] M. Abramowitz and I.A. Stegun, Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables (Dover, New York, 1972).

Datas de Publicação

• Publicação nesta coleção
  02 Jul 2010
• Data do Fascículo
  Mar 2010

Histórico

• Recebido
  28 Maio 2009
• Aceito
  21 Set 2009