Acessibilidade / Reportar erro

Sobre o método semi-analítico de Chhajlany e Malnev para soluções aproximadas não-perturbativas da Equação de Schrödinger com potencial polinomial par

On the semi-analytical Chhajlany and Malnev method for aproximate non-perturbative solutions of the Schrödinger equation with even polynomial potential

Resumos

Propomos, como parte do estudo introdutório de Mecânica Quântica na graduação, o uso do método de Chhajlany e Malnev (MCM) para obtenção de soluções aproximadas semi-analíticas da equação de Schrödinger unidimensional (autovalores e autofunções de energia), como uma alternativa a métodos perturbativos (que trazem questões de convergência de séries infinitas), a métodos como o WKB (de caráter semi-clássico, e não totalmente quântico), e a métodos como o de diferenças finitas (puramente numéricos), no estudo de potenciais polinomiais de potências pares. Potenciais polinomiais surgem, por exemplo, como potenciais efetivos no estudo de oscilações em torno de mínimos locais de um dado potencial. Neste trabalho, desenvolvemos em detalhe o MCM para o osciladores quânticos harmônico e anarmônico quártico, utilizando um código Fortran que implementa o método, nestes casos. Pré-requisitos para a compreensão do MCM são a familiaridade com a Equação de Schrödinger, bem como ferramentas básicas do cálculo íntegro-diferencial e álgebra linear.

Palavras-chave:
quantização de sistemas hamiltonianos; método de Chhajlany e Malnev; oscilador harmônico quântico; oscilador quártico quântico


We propose, as part of the introductory study of Quantum Mechanics at undergraduate level, to employ the Chhajlany and Malnev method (MCM) for obtaining approximate semi-analytical solutions of the unidimensional Schrödinger equation (energy eigenvalues and eigenfunctions), as an alternative to perturbative methods (which bring issues of convergence of infinite series), to methods like WKB (of semi-classical character, and not completely quantum), and to methods such as the finite difference method (purely numerical), for the study of polynomial potentials with even powers. Polynomial potentials appear, for instance, as effective potentials in the study of oscillations about local minima of a given potential. In the present work, we develop in detail the MCM for the harmonic and quartic anharmonic quantum oscillators, using a Fortran code that implements that method, for those cases. Pre-requisites for understanding the MCM are the familiarity with Schrödinger equation, as well as with basic tools of integral-differential calculus and linear algebra.

Keywords:
quantization of hamiltonian systems; method of Chhajlany and Malnev; quantum harmonic oscillator; quantum quartic oscillator


1. Introdução

A equação de Schrödinger é ponto fundamental em um curso de introdução à Mecânica Quântica na graduação, tanto para a compreensão conceitual desta disciplina quanto para a capacitação do aluno no uso efetivo da teoria quântica. Para obterem-se soluções fisicamente razoáveis desta equação de uma partícula, são explorados e exemplificados métodos de solução através de um elenco de potenciais (e.g., potencial constante, a barreira de potencial, o poço de potencial, o oscilador harmônico e o átomo de hidrogênio) concomitantemente à introdução de diversos conceitos para a solução sob estes potenciais (e.g., separação de variáveis, onda plana, princípio de superposição, limites assintóticos, continuidade da função de onda e suas consequências, coeficientes de transmissão e reflexão, a quantização propriamente dita de grandezas físicas em decorrência das condições de contorno, e o operador momentum angular) [1[1] W. Greiner, Quantum Mechanics: An Introduction (Springer-Verlag, Berlin, 2001), 4th ed.4[4] David J. Griffiths, Introduction to Quantum Mechanics (Prentice-Hall, New Jersey, 1995).]. Para os potenciais mais simples introduzem-se métodos analíticos de solução, visando a obtenção de resultados exatos. Para os casos em que isto não é possível estudam-se métodos aproximados tradicionais como a teoria da perturbação [2[2] L.D. Landau and E.M. Lifshitz, Quantum Mechanics: Non-Relativistic Theory (Pergamon, Oxford, 1977), 3rd ed.] e o método WKB [4[4] David J. Griffiths, Introduction to Quantum Mechanics (Prentice-Hall, New Jersey, 1995).]. (Geralmente a aplicação de métodos puramente numéricos como o de diferenças finitas [5[5] J.W. Thomas, Numerical Partial Differential Equations: Finite Difference Methods (Springer, New York, 1995).] não fazem parte do escopo desta disciplina introdutória.)

No presente trabalho sugerimos o aprendizado, no contexto de um curso de graduação ou atividade de iniciação científica, do método semi-analítico de Chhajlany e Malnev (MCM) [6[6] S.C. Chhajlany, D.A. Letov and V.N. Malnev, J. Phys. A 24, 2731 (1991).], que se aplica a potenciais V(x) polinomiais com potências pares de x, para sistemas com domínio em x(-,) (reta real) ou x[0,) (semi-reta, apropriado a potenciais radiais). Potenciais polinomiais surgem, por exemplo, no estudo de oscilações em torno de um ponto de mínimo do potencial, aproximado via série de Taylor, o que revisaremos na seção 2. Na reta real, o termo assintoticamente dominante (i.e., de potência mais alta) deve possui potência par e coeficiente positivo, para que limx±V(x)=+ e o espectro de energia do sistema seja limitado inferiormente (em outras palavras, para que exista um estado fundamental). A restrição do método a potenciais com potências exclusivamente pares limita-o a potenciais simétricos em torno do ponto de mínimo (assimetrias seriam introduzidas pelas potências ímpares). O método é descrito em detalhes na seção 3., para o potencial quártico da forma V(x)=ax2+bx4 na semi-reta. Alguns dos autores do presente trabalho abordaram este potencial via MCM, em um estudo de Cosmologia [7[7] 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).].

Uma vantagem da utilização deste método reside no emprego de autofunções de energia manifestamente normalizáveis — um polinômio multiplicado por uma exponencial de um outro polinômio par que tenda a - nos limites assintóticos, fazendo com o que a autofunção tenda assintoticamente a zero — evitando as divergências do método perturbativo [6[6] S.C. Chhajlany, D.A. Letov and V.N. Malnev, J. Phys. A 24, 2731 (1991).,8[8] C.M. Bender and T.T. Wu, Phys. Rev. 184, 1231 (1969).] e sutilezas no tratamento de uma série infinita [9[9] J.E. Drummond, J. Phys. A: Math. Gen. 14, 1651 (1981).] que não caberiam em um curso introdutório. A forma das autofunções aproximadas deste método facilita o cálculo de valores esperados de grandezas físicas do sistema. Ele também permite um tratamento totalmente quântico, e não semi-clássico como o WKB.

Para um eficiente acompanhamento do método pelo aluno, recomenda-se a familiaridade com a equação de Scrhödinger, e com as ferramentas básicas do cálculo diferencial e álgebra linear (em especial o do cálculo de autovalores e autovetores de matrizes de dimensão finita).

O presente trabalho é organizado da seguinte maneira. Na seção 2. faremos uma rápida revisão de oscilações harmônicas e anarmônicas clássicas. Na seção 3., apresentamos o MCM em detalhes. Na seção 4. exemplificamos este método com dois osciladores quânticos (harmônico e quártico) na semi-reta, comparando-os com resultados analíticos exatos (no caso harmônico) e resultados obtidos por dois outros métodos, o de diferenças finitas [5[5] J.W. Thomas, Numerical Partial Differential Equations: Finite Difference Methods (Springer, New York, 1995).] e o espectral de Galerkin [10[10] E.V. Corrêa Silva a, G.A. Monerat, G. de Oliveira Neto and L.G. Ferreira Filho, Computer Physics Communications 185, 380 (2014).] (no caso quártico). A seção 5. apresenta nossas conclusões. O apêndice A Apêndice A O Programa CHHAJLANY Implementamos o MCM, descrito na seção 3., em um programa em linguagem Fortan, aqui denominado CHHAJLANY.FOR, cujo código fonte encontra-se disponível no apêndice deste artigo. Ele permite também a construção de pacotes de onda e do cálculo do valor esperado da posição. Em relação à avaliação da precisão dos resultados, foram criadas duas ferramentas, tendo por base a verificação da conservação da norma da função de onda e a inspeção de como o potencial U(x) (10) se aproxima do potencial original V(x) (9), comportamento este determinado pela variável auxiliar c. Em relação a c, observamos que seu valor diminui à medida que aumentamos o grau de truncamento Ω do polinômio v(x) (17), garantindo assim que U(x)→V(x). Para calcularmos o valor de c, usamos a rotina zplrc(coef, zero, ndeg) da biblioteca IMSL embarcada no compilador FORTRAN Profissional. O parâmetro coef é um vetor com dimensão ndeg+1, com ndeg=3, ou seja, coef contém os quatro coeficientes da Eq.(27). O parâmetro zero é um vetor de dimensão ndeg contendo as raízes da referida equação. Na sequência o programa seleciona a raiz real e positiva resultante e atribui este valor à variável auxiliar c. Esse resultado é mostrado na Fig. 5. Os auto-valores e auto-vetores da matriz dada pela Eq.(21) são obtidos através da subrotina evcrg(a, eval, evec), onde eval é um vetor de dimensão (Ω+1)∕2 contendo os auto-valores e evec é uma matriz de dimensão (Ω+1)∕2×(Ω+1)∕2 contendo na n-esima coluna os coeficientes D correspondentes à expansão Eq.(17) do auto-estado n. O parâmetro SSE é a soma dos quadrados dos erros. Além disso, o programa também faz a evolução de um pacote de onda em particular. Não comentaremos aqui esta sua funcionalidade. descreve o programa Fortran CHHAJLANY.FOR que utilizamos em nossos exemplo, cujo código fonte se encontra no apêndice B Apêndice B Código fonte CHHAJLANY.FOR .

2. Oscilações Harmônicas e Anarmônicas

Uma situação física importante discutida nos cursos de Mecânica Clássica e Quântica são oscilações de uma partícula em torno de um ponto de equilíbrio clássico x0, o que corresponde a um mínimo local do potencial V(x). Considere a situação exemplificada na Fig.1. Classicamente, uma partícula com uma certa energia mecânica E, em uma vizinhança suficientemente pequena em torno de um ponto de mínimo do potencial ficará confinada à esta região, entre os dois pontos de retorno A e B (nos quais a energia cinética é nula, logo E=V(x)). A força resultante sobre esta partícula é F=-dV(x)dx; portanto ela é positiva em A (a derivada do potencial aqui é negativa, então a força neste ponto também é positiva) e negativa em B (a derivada aqui é positiva, portanto a força é negativa), fazendo com que a particula sempre retorne em direção ao ponto de mínimo. Neste ponto temos dVdxx=x0=0; portanto, nas vizinhanças deste, a expansão em série de Taylor do potencial é da forma

(1)V(x)V0+V2(x-x0)2+V3(x-x0)3+O(x-x0)4,
em que V0=V(x0), Vk=1k!dkVdxkx=x0, para k{1,2,3,}, e O(x-x0)4 indica termos em potências iguais ou superiores a 4 da distância x-x0 ao ponto de equilíbrio. O termo na derivada de primeira ordem é nulo, (V1=0). A força resultante sobre esta partícula é
(2)F(x)-2V2(x-x0)-3V3(x-x0)2+O(x-x0)3.

Figure 1
Região de oscilação clássica em torno de um mínimo do potencial V(x). Nos pontos A e B a energia mecânica E da partícula é igual à energia potencial V(x), e consequentemente a sua energia cinética é nula. Além disso, a força resultante sobre a partícula é F=-dVdx; em A esta derivada é negativa, logo esta força é orientada positivamente, enquanto que em B a derivada é positiva, portanto a força é orientada negativamente. Deste modo, ao alcançar qualquer um destes pontos, denominados pontos de retorno, a partícula retornará em direção ao ponto de equilíbrio, permanecendo perpetuamente na região entre A e B. As chamadas oscilações harmônicas ocorrem quando a energia mecânica E for suficientemente próxima do mínimo do potencial, de tal modo que na região de oscilação [xA,xB] seja pequena o suficiente para que aproximemos o potencial como uma função quadrática de (x-x0), equivalente à função potencial do oscilador harmônico. Caso contrário, teremos oscilações ditas anarmônicas.

Para a descrição do comportamento dinâmico da partícula, o termo constante V0 da expansão (1) não é relevante, pois sua derivada é nula e assim não contribui para a força sobre a partícula. No que segue, faremos V0=0, o que equivale a redefinirmos o potencial como V(x)V(x)-V0 e portanto sem perda de generalidade.

Se (x-x0) for suficientemente pequeno para que possamos consideramos a expansão do potencial somente até 2a ordem,

(3)V(x)V2(x-x0)2+O(x-x0)3,
teremos um potencial do tipo do oscilador harmônico, gerando uma força restauradora de Hooke da forma
(4)F(x)-2V2(x-x0)+O(x-x0)2.
proporcional ao deslocamento em relação ao ponto de equilíbrio 1 1 Podemos definir uma constante elástica efetiva k tal que k=2V2. . Esta situação aproximada é denominada regime de oscilações harmônicas. No caso em que (x-x0) não seja tão pequeno e que tenhamos que considerar termos de ordens superiores ao quadrático na expansão (1), estaremos no regime de oscilações anarmônicas.

A equação de Schrödinger [1,3,4,11-13]

(5)-22m2x2+V(x)Ψ(x,t)=itΨ(x,t),
governa a evolução temporal da função de onda complexa Ψ(x,t), representativa do estado de uma partícula quântica unidimensional de massa m sob a ação de um potencial V(x) contínuo por partes (sistema conservativo), em um certo intervalo de tempo (ta,tb) sem medições do sistema. Juntamente com a equação de Schrödinger independente do tempo
(6)-22md2dx2+V(x)ψ(x)=Eψ(x),
forma o ponto-chave de um curso introdutório à Mecânica Quântica 2 2 Analogamente, nas disciplinas de Mecânica clássica, a 2a Lei de Newton, que governa a dinâmica de uma partícula unidimensional, é a base para a modelagem de situações mais complexas, envolvendo maiores dimensões e conjuntos de partículas. , para os objetivos de compreensão e aplicação da teoria quântica que se almejam nesta disciplina.

Conhecer uma variedade de métodos para abordar esta equação é de grande valor para capacitar o aluno a prever, dado um certo potencial V(x), o comportamento quântico de uma partícula. O caso de oscilações harmônicas pode ser resolvido exatamente e suas soluções são bem conhecidas [1,3,4,11-13], sendo este sistema comumente abordado nos cursos de graduação. Já o caso anarmônico não possui solução exata, e é uma oportunidade para o aprendizado de métodos aproximados como a teoria da perturbação e o método WKB, também comumente abordados. Neste artigo, propomos o estudo de um certo tipo de potencial anarmônico através do MCM. que discutimos na seção 3.

3. O Método de Chhajlany e Malnev

3.1. O tipo de potencial em questão

O método de Chhajlany e Malnev (MCM) [6[6] S.C. Chhajlany, D.A. Letov and V.N. Malnev, J. Phys. A 24, 2731 (1991).] aplica-se a potenciais polinomiais pares,

(7)V(x)=n=1Mpnx2n,
em que x, M, pn para n{1,2,,M}, e pM>0. (O termo de ordem mais alta deve ter o coeficiente positivo, para que o espectro de energia do sistema seja limitado inferiormente.) Define-se um potencial auxiliar
(8)U(x)=V(x)+c2x2(M+1)
em que c é uma variável auxiliar real e positiva. O potencial auxiliar, também um polinômio par, resulta portanto do acréscimo do termo não-negativo c2x2(M+1) ao potencial original, aumentando-se o grau do polinômio em duas unidades.

3.2. O potencial considerado neste artigo

Consideraremos neste artigo um oscilador anarmônico com potencial quártico

(9)V(x)=ax2+bx4,
em que a e b são parâmetros dados do sistema, tal que b>0. O potencial auxiliar é da forma
(10)U(x)=ax2+bx4+c2x6,
de tal modo que a equação de Schrödinger independente do tempo assume a forma
(11)d2dx2η(x)+ϵ-ax2-bx4-c2x6η(x)=0.

3.3. Ansatz de solução

Em geral, desejam-se soluções normalizáveis. Para tanto, buscam-se soluções aproximadas contínuas da forma

(12)η(x)=v(x)e-c4x4-γ2x2,
em que v(x) é uma função polinomial (sem paridade definida, a princípio), e γ>0 uma nova variável auxiliar, em princípio independente de c. Os fatores 12 e 14 no expoente são convenientes para a simplicidade das expressões que seguirão. Observe que este ansatz garante que
(13)limx±η(x)=0,
e portanto estas soluções são normalizáveis em (-,). Substituindo (12) em (11) temos:
(14)d2v(x)dx2-2(cx3+γx)dv(x)dx+ϵ-γ+γ2-a-3cx2+(2cγ-b)x4v(x)=0.

As raízes da parte polinomial da solução (12) corresponderão aos nodos aproximados das autofunção de energia.

Uma escolha simplificadora para γ é

(15)γ=b2c,
que reduz (14) a
(16)d2v(x)dx2-2(cx3+γx)dv(x)dx+ϵ-γ+γ2-a-3cx2v(x)=0.

3.4. Soluções polinomiais para v(x) e relações de recorrência entre coeficientes

Escrevamos v(x) como uma série infinita (a princípio) de potências inteiras,

(17)v(x)=n=0Dnxn,
em que os coeficientes Dn são constantes a serem determinadas. Substituindo (17) em (16), obtemos
(18)n=2n(n-1)Dnxn-2+n=1-2cnDnxn+2-2γnDnxn+n=0(ϵ-γ)Dnxn+n=0γ2-a-3cDnxn+2=0.

O termo de menor potência comum a todos os somatórios em (18) é o termo cúbico. Extraindo os termos de ordem zero, lineares e quadráticos de cada somatório, obtemos

(19)2D2+6D3x+12D4x2+n=5n(n-1)Dnxn-2-2cn=1nDnxn+2-2γD1x-4γD2x2-2γn=3nDnxn+(ϵ-γ)D0+(ϵ-γ)D1x+(ϵ-γ)D2x2+(ϵ-γ)n=3Dnxn++γ2-a-3cD0x2+γ2-a-3c×n=1Dnxn+2=0.

Em cada somatório da forma n=αF(n) na expressão acima, redefinimos o índice como m=n-α, reescrevendo o somatório como m=0F(m+α). Deste modo, conseguimos reescrever todos os somatórios em um mesmo índice m com um mesmo valor mínimo para este índice (zero, no caso). Após coletarmos os termos de mesma potência em x, obtemos

(20)2D2+(ϵ-γ)D0+6D3+(ϵ-3γ)D1x+12D4+(ϵ-5γ)D2+(γ2-a-3c)D0x2+m=0(m+5)(m+4)Dm+5+ε-γ(2m+7)Dm+3+γ2-a-c(2m+5)Dm+1xm+3=0.

Para que a série no lado esquerdo de (20) seja identicamente nula (i.e., nula para todo x) todos os seus coeficientes devem ser nulos. Obtemos então o seguinte sistema algébrico de equações de recorrência para os coeficientes Dk de v(x):

(21){Ordem 0, 2D2+ (-γ)D0= 0;Ordem 1, 6D3+ (- 3γ)D1= 0;Ordemk, tal quek 2,(k+ 2)(k+ 1)Dk+2+ [-γ(2k+ 1)]Dk+ [γ2-a-c(2k- 1)]Dk-2= 0.

A última equação do sistema (21) foi obtida fazendo-se m=k-3 no somatório em (20). A condição m0 deste somatório traduzir-se-ia então como k3, mas podemos incorporar a esta equação o caso de ordem 2 fazendo k2.(Na verdade, podemos incorporar também as equações de ordem zero e de primeira ordem nesta única equação de recorrência; basta estabelecermos que Dk=0, se k<0).

Podemos obter soluções (aproximadas) com v(x)polinomial ao truncarmos a série (17) em um certo grau Ω — o que equivale a impormos

(22)DΩ0 e Dk=0,parak>Ω.

Observe ainda que as relações de recorrência (21) mantém os coeficientes com índices pares independentes dos coeficientes com índices ímpares.

Se Ω for par, v(x) será um polinômio par,

(23)v(x)=D0+D2x2+D4x4++DΩxΩ,
com Npar(Ω)=Ω2+1 coeficientes {D0,D2,D4,,DΩ}. As equações relevantes aqui são aquelas obtidas de (21) que envolvem estes coeficientes; portanto, todas aquelas equações de ordem par tal que kΩ+2.

Por outro lado, se Ω for ímpar, v(x) será um polinômio ímpar,

(24)v(x)=D1x+D3x3+D5x5++DΩxΩ,
com Nímpar(Ω)=12Ω+1 coeficientes {D1,D3,D5,,DΩ}. As equações relevantes aqui são aquelas obtidas de (21) que envolvem estes coeficientes; portanto, todas aquelas equações de ordem ímpar tal que kΩ+2.

Portanto, o número de coeficientes do polinômio v(x), de acordo com o seu grau Ω, será dado por

(25)N(Ω)={Npar(Ω)=Ω2+1, seΩfor par;Nímpar(Ω)=12(Ω+1), seΩfor ímpar.

Este é o número de variáveis envolvidas no problema de determinarmos v(x).

3.5. Determinação de c e γ.

Qualquer que seja a paridade de Ω, a equação de (21) correspondente a k=Ω+2, levando-se em conta a condição de truncamento (22), não fixa DΩ, mas estabelece uma nova restrição a γ,

(26)γ2-a-c(2Ω+3)DΩ=0γ2-a-c(2Ω+3)=0
que, combinada com (15) de modo a eliminar γ, nos leva a uma equação polinomial cúbica para a variável auxiliar c, em termos dos parãmetros a,b e do grau de truncamento Ω.
(27)42Ω+3c3+4ac2-b2=0.

O valor de c obtido deve ser real e positivo para que a solução (12) seja normalizável. (O valor de γ é determinado a partir do valor de c por (15), e ele também será positivo se c o for, pois b>0 por hipótese.) Podemos deduzir que sempre haverá uma raiz real e positiva colocando esta equação na forma

(28)c3+4a42Ω+3c2-b242Ω+3=0.

Uma das relações de Girard [14[14] G. Iezzi, Fundamentos de Matemática Elementar, vol. 6: Complexos, Polinõmios, Equações (Atual, São Paulo, 1977), 2ª ed.], que relacionam os coeficientes de uma equação polinomial com as suas raízes, estabelece que o produto das raízes de (28) deve ser igual a b242Ω+3, sendo portanto uma quantidade positiva. Os coeficientes de (28) são reais; logo, ou todas as suas raízes são reais, ou então temos uma raiz real e duas raízes complexas conjugadas, não havendo possibilidade de haver três raízes complexas. Se todas as raizes forem reais, o produto delas só será positivo se não houver nenhuma raiz nula e houver três raizes positivas ou então uma raiz positiva e duas negativas. Se houver, por outro lado, duas raízes complexas conjugadas, o produto delas é sempre positivo, então a raiz real restante deverá ser positiva. No caso de haver mais de uma raiz positiva, considera-se a de menor valor, de tal modo que o potencial auxiliar U(x) (10) seja o mais próximo possível do potencial original V(x) (9); esta raiz diminui monotonicamente, à medida em que Ω cresce [15[15] S.C. Chhajlany and V.N. Malnev, Phys. Rev. A 42, 3111 (1990).]. Como consequência, U(x)V(x), ou seja, o espectro de soluções relacionado a U(x) se aproxima cada vez mais do espectro de soluções para V(x).

Na solução do problema para um certo grau de truncamento Ω, portanto, utilizaremos as equações (21) até kΩ-2, tendo k sempre a mesma paridade de Ω, e além disso utilizaremos a equação (27), o que corresponde a k=Ω em (21).

3.6. Obtendo autovalores e autofunções de energia

Vamos a seguir mostrar, à guisa de exemplo, como determinar os autovalores e autovetores de energia aproximados do sistema, para dois valores do grau de truncamento, Ω=2 e Ω=7.

Considere primeiramente Ω=2. Estamos portanto considerando soluções da forma (12), com um polinômio par v(x) de segundo grau, da forma

(29) v ( x ) = D 0 + D 2 x 2 ,

de acordo com (17). De (21) e (22) para Ω=2 obtemos um sistema de equações lineares para as variáveis {D0,D2},

(30)2D2+(ϵ-γ)D0=0(ordem 0),
(31)(ϵ-5γ)D2+(γ2-a-3c)D0=0(ordem 2),
e também (27), que fornece
(32)28c3+4ac2-b2=0.

O sistema formado por (30) e (31) pode ser expresso matricialmente como

(33)ϵ-γ2γ2-a-3cϵ-5γD0D2=0.

Temos uma matriz de coeficientes de dimensão N(2)×N(2)=2×2, Para que este sistema admita soluções não-triviais (i.e., outras soluções além de D0=D2=0), a matriz dos coeficientes deve ter determinante nulo,

(34)ϵ-γ2γ2-a-3cϵ-5γ=0.

Esta é uma equação polinomial de segundo grau para ϵ.

Dados valores para {a,b}, determinamos c e γ através de (32) e (15), respectivamente. Considere o potencial (9) em x com a=1 e b=1, por exemplo; obtemos c=0.2879313419735601 e γ=1.736525091616854. Resolvemos então (34) para ϵ, determinando os autovalores de energia do sistema. Obtemos os autovalores ϵ0=1.419386847997858 e ϵ1=8.999763701703262. Substituímos então cada ϵ obtido no sistema matricial (33) e então determinamos a relação entre D0 e D2 a partir da primeira linha da matriz, por exemplo. (A segunda linha é linearmente dependente da primeira, dado que o determinante do sistema é nulo.) Para ϵ=ϵ0 obtemos D2=0.1585691218094980D0; portanto, a autofunção correspondente a este autovalor é, de acordo com o ansatz (12),

(35)η0(x)=D0(1)(1+0.1585691218094980x2)×e-0.07198283549339002x4-0.8682625458084270x2.

Para ϵ=ϵ1 obtemos D2=-3.631619305043204D0; portanto, a autofunção correspondente a este autovalor é

(36)η1(x)=D0(2)(1-3.631619305043204x2)×e-0.07198283549339002x4-0.8682625458084270x2.

Os superíndices (1) e (2) foram acrescentados apenas para distinguir as constantes arbitrárias remanescentes, em cada caso. Estas constantes serão determinadas pela condição de normalização da solução,

(37)-dxη(x)2=-dxe-c2x4-γx2v(x)2=1.

(O intervalo de integração, no caso, é x(-,+); se estivéssemos lidando com um problema na semi-reta, seria x[0,+).) Aplicando esta condição, obtemos D0(1)=0.8401401848348102 e D0(2)=0.6937218650838422.

Apesar da relativa simplicidade, obtemos assim aproximações ao estado fundamental do sistema e ao seu primeiro estado excitado, e informações sobre o gap ou intervalo de energia entre eles, que é informação relevante, por exemplo, no estudo de excitações do estado fundamental em sistemas de matéria condensada.

Consideremos agora o grau de truncamento Ω=7. O polinômio v(x) é da forma

(38)v(x)=D1x+D3x3+D5x5+D7x7.

Obtemos das equações de ordem ímpar de (21), o sistema de dimensão N(7)×N(7)=4×4

(39)ϵ-3γ600(γ2-a-5c)ϵ-7γ2000γ2-a-9cϵ-11γ4200γ2-a-13cϵ-15γD1D3D5D7=0.
e, da equação (27),
(40)68c3+4ac2-b2=0.

O valor de γ é obtido do valor de c pela equação (15). Os autovalores ϵ de energia serão determinados pela equação polinomial

(41)ϵ-3γ600(γ2-a-5c)ϵ-7γ2000γ2-a-9cϵ-11γ4200γ2-a-13cϵ-15γ=0.

Considere, por exemplo, a=36 e b=1.2 para o potencial (9) na semi-reta (x[0,)). Resolvendo (40), obtemos c=0.09776855216831178; e de (15) obtemos γ=6.13694267423620. Resolvendo a seguir (41), Obtemos as raízes ϵ0=18.12370382580685, ϵ1=42.61257182437350, ϵ2=67.47968149727072 e ϵ3=92.71397912505204. As correspondentes autofunções normalizadas (na semi-reta) são listadas a seguir.

(42) η 0( x ) = (5.803441504090680 x + 0.2777180801880438 x 3 + 0.004418522510481575 x 5 +0 .2337292030734449 × 10 4 x 7 ) e −0.02444213804207794 x 4 −3.068471337118098 x 2
(43) η 1( x ) = (7.178089385895051 x + 28.95371338381371 x 3 + .9220130214074282 x 5 −0 .7292962687093990 × 10 2 x 7 ) e −0.02444213804207794 x 4 −3.068471337118098 x 2
(44) η 2( x ) = (8.101476728663191 x − 66.25502908772484 x 3 + 80.75701084456534 x 5 +1.285154826970363 x 7 ) e −0.02444213804207794 x 4 −3.068471337118098 x 2
(45) η 3( x ) = (8.830319787947538 x − 109.3534309143132 x 3 + 271.5280811403464 x 5 −160.9296019369097 x 7 ) e −0.02444213804207794 x 4 −3.068471337118098 x 2

Os gráficos destas autofunções são exibidos na Fig. 2.

Figura 2
Autoestados para o caso Ω=7, para a=36 e b=1.2 no potencial (9).

Generalizando para um polinômio v(x) com grau Ω ímpar, as relações de recorrência (21), juntamente com as condições (22), resultam no sistema linear de equações

(46) ϵ - 3 γ 6 0 0 0 0 0 0 0 γ 2 - α - 5 c ϵ - 7 γ 20 0 0 0 0 0 0 0 0 0 A n B n C n 0 0 0 0 0 0 0 0 0 0 γ 2 - a - ( 2 Ω - 1 ) c ( ϵ - ( 2 Ω + 1 ) γ ) D 1 D 3 D 5 D 2 n + 1 D Ω = 0 0 0 0 0

em que n{0,1,,Ω-12} e

(47) A n = γ 2 - a - ( 4 n + 1 ) c ,
(48) B n = ϵ - ( 4 n + 3 ) γ ,
(49) C n = ( 2 n + 2 ) ( 2 n + 3 ) .

De modo análogo, para um polinômio v(x) com grau Ω par obtemos o sistema

(50)ϵ-γ20000000γ2-α-3cϵ-5γ12000000000FnGnHn0000000000γ2-a-(2Ω-1)c(ϵ-(2Ω+1)γ)D0D2D4D2nDΩ=00000
em que n{0,1,,Ω2} e
(51)Fn=γ2-a-(4n-1)c,
(52)Gn=ϵ-(4n+1)γ,
(53)Hn=(2n+1)(2n+2).

Em ambos os casos, utilizamos também a equação (27) para determinação do c e a equação (15) para o cálculo de γ.

4. Exemplos

Nesta seção aplicaremos o MCM para obtermos os autovalores e autofunções de energia em dois osciladores, o harmônico e o anarmônico quártico, na semi-reta, ou seja, x[0,) (potencial radial). A subseção 4.1. mostra os resultados da aplicação do MCM ao oscilador harmônico, enquanto que a subseção 4.2. faz o mesmo para o oscilador anarmônico quártico.

4.1. O Oscilador Harmônico na Semi-Reta

Nesta subseção mostramos os resultados da aplicação do MCM, desenvolvido na seção 3., a um oscilador harmônico simples na semi-reta (x0). com potencial da forma

(54)V(x)=12kx2,
e sob a condição Ψ(x=0,t)=0. As soluções analíticas (autofunções e espectro de energia) para este sistema físico são conhecidas, correspondendo às autofunções ímpares e correspondentes energias do oscilador harmônico na reta real.

A Tabela 1 compara os resultados fornecidos pelo programa CHHAJLANY, que implementa o MCM (vide apêndices 5. e 5.), com aqueles obtidos utilizando-se o método de diferenças finitas [16[16] G.A. Monerat, L.G. Ferreira Filho, E.V. Corrêa Silva, G. Oliveira-Neto, P.H.A.S. Nogueira e A.R.P. de Assumpção, Rev. Bras. de Ens. de Fis. 32, 1304 (2010).] e o método espectral de Galerkin [17[17] G.A. Monerat, E.V. Corrêa Silva, L.B. Leal e G. Oliveira-Neto, Rev. Bras. de Ens. de Fis. 37, 4301 (2015).]. Utilizamos um polinômio v(x) de ordem Ω=100, uma partícula de massa m=6 e uma constante elástica k=6, de tal modo que a frequência seja ω=1. Encontramos c=0. O MCM consegue, portanto, reproduzir com perfeição o espectro do oscilador harmônico na semi-reta. A Fig. 3 mostra algumas autofunções calculadas.

Tabela 1
Espectro de energia para os 5 primeiros níveis de energia do oscilador harmônico na semi-reta (com =1 e ω=1). Aqui, En=ωn+12 representam os níveis de energia exatos do oscilador harmônico na reta real (utilizamos apenas os valores ímpares de n); En(CM) são os níveis de energia calculados via MCM; En(df) os níveis obtidos via método das diferenças finitas; e En(esp) os níveis fornecidos pelo método espectral. Os erros absolutos dos níveis de energia calculados via método das diferenças finitas e via método espectral em relação aos níveis exatos são dados por δEn(df) e δEn(esp), respectivamente. Observe que os erros absolutos associados ao método espectral são muito menores do que os do método de diferenças finitas.
Figura 3
Os oito primeiros autoestados aproximados obtidos para o oscilador harmônico através do MCM, para um polinômio de grau Ω=100.

4.2. O Oscilador Anarmônico Quártico na semi-reta

Estudaremos a seguir o oscilador anarmônico quártico para uma partícula de massa m=6, caracterizado por um potencial na forma

(55)V(x)=ax2+12Λx4,
em que Λ é um parâmetro do sistema.

Para estimarmos a precisão do MCM, quando aplicado ao oscilador anarmônico quártico, fizemos o seguinte experimento numérico: aumentamos gradualmente o grau do polinômio característico até o valor máximo Ω=14001, com objetivo de verificarmos a convergência no espectro de energia. A Tabela 2 exibe nossos resultados do experimento. De acordo com esta tabela, o erro absoluto do espectro de energia é da ordem de 10-3, para polinômios de ordem até 14001. A Fig. 6 mostra o descréscimo da variável auxiliar c à medida em que Ω cresce. As autofunções para o caso do oscilador anarmônico quártico também foram obtidas pelo método e os seis primeiros autoestados são mostrados na Fig. 4. Analisamos também como o valor da energia do estado fundamental depende do valor da constante Λ que aqui é negativa. A Fig. 5 mostra essa dependência.

Tabela 2
Espectro de energia para os 5 primeiros níveis de energia do oscilador anarmônico quártico (com =1, a=36, Λ=-0.1 temos b=-1.2 e w=1). Aqui, os En representam os níveis de energia calculados via MCM, para polinômios característicos de diferentes ordens. O valor da variável auxiliar c encontrado em cada caso foi: para Ω=53 obtivemos c=0.08.877763929333429, para Ω=8001 obtivemos c=0.02749793712783993 e para Ω=14001 obtivemos c=0.02300475316347050. Observe que dependendo da ordem do polinômio característico obtemos uma precisão de 4 casas decimais.
Figura 4
Os seis primeiros autoestados aproximados e seus respectivos autovalores para o oscilador anarmônico quártico. Aqui utilizamos k=1 e Λ=-0.1, parâmetros estes mostrados no potencial em (55).
Figure 5
diferença de energia ΔEn=EnΛ-En0, para os cinco primeiros estados do oscilador harmônico na semi-reta (Λ=0), mostra como os autovalores da energia variam com Λ.
Figure 6
Comportamento da variável auxiliar c em função da ordem do polinômio. Em todos os casos observamos que c0 quando Ω.

5. Conclusão

Aplicamos o método semi-analítico e não perturbativo de Chhajlany e Malnev (MCM) a dois sistemas físicos: um oscilador harmônico e um oscilador anarmônico quártico. Em ambos os casos o método forneceu o espectro de energia e seus autoestados. No caso do oscilador harmônico, conseguimos com excelente aproximação o espectro de energia conhecido, usando um polinômio característico de ordem 100. Comparados os resultados àqueles fornecidos pelo métodos de diferenças finitas e pelo método espectral de Galerkin, o MCM mostrou-se mais eficaz. No caso do oscilador anarmônico quártico, verificou-se uma convergência no espectro com uma precisão entre 10-4 a 10-3, nos primeiros 10 níveis de energia, tendo sido usado para isso um polinômio característico de ordem 7000. Os resultados se tornam mais precisos a medida que polinômios de ordens mais elevadas são considerados. O método é eficaz, mas possui grande custo computacional, quando comparado a outros métodos; por outro lado, ele é bem simples, podendo ser apresentado a alunos de graduação em Física.

  • 1
    Podemos definir uma constante elástica efetiva k tal que k=2V2.
  • 2
    Analogamente, nas disciplinas de Mecânica clássica, a 2a Lei de Newton, que governa a dinâmica de uma partícula unidimensional, é a base para a modelagem de situações mais complexas, envolvendo maiores dimensões e conjuntos de partículas.

Agradecimentos

G. A. Monerat agradece à UERJ/FAPERJ pela bolsa Prociência, P. A. P. Simas e Silva e L. B. Leal agradecem ao programa PIBIC/UERJ pelas bolsas de Iniciação Científica. Este trabalho foi parcialmente realizado no Laboratório de Computação Avançada (LCA) do Depto. de Matemática, Física e Computação da Faculdade de Tecnologia da UERJ.

Referências

  • [1]
    W. Greiner, Quantum Mechanics: An Introduction (Springer-Verlag, Berlin, 2001), 4th ed.
  • [2]
    L.D. Landau and E.M. Lifshitz, Quantum Mechanics: Non-Relativistic Theory (Pergamon, Oxford, 1977), 3rd ed.
  • [3]
    C. Cohen-Tannoudji, B. Diu and F. Laloe, Quantum Mechanics, Vols. 1 and 2 (John Wiley & Sons, New York, 1977).
  • [4]
    David J. Griffiths, Introduction to Quantum Mechanics (Prentice-Hall, New Jersey, 1995).
  • [5]
    J.W. Thomas, Numerical Partial Differential Equations: Finite Difference Methods (Springer, New York, 1995).
  • [6]
    S.C. Chhajlany, D.A. Letov and V.N. Malnev, J. Phys. A 24, 2731 (1991).
  • [7]
    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).
  • [8]
    C.M. Bender and T.T. Wu, Phys. Rev. 184, 1231 (1969).
  • [9]
    J.E. Drummond, J. Phys. A: Math. Gen. 14, 1651 (1981).
  • [10]
    E.V. Corrêa Silva a, G.A. Monerat, G. de Oliveira Neto and L.G. Ferreira Filho, Computer Physics Communications 185, 380 (2014).
  • [11]
    J.J. Sakurai and J. Napolitano, Mecânica Quântica Moderna (Bookman, Porto Alegre, 2013), 2ª ed.
  • [12]
    E. Merzbacher, Quantum Mechanics (Wyley, New York, 1998), 3rd ed.
  • [13]
    P.A.M. Dirac, The Principles of Quantum Mechanics (Clarendon Press, Oxford, 1967), 4th ed.
  • [14]
    G. Iezzi, Fundamentos de Matemática Elementar, vol. 6: Complexos, Polinõmios, Equações (Atual, São Paulo, 1977), 2ª ed.
  • [15]
    S.C. Chhajlany and V.N. Malnev, Phys. Rev. A 42, 3111 (1990).
  • [16]
    G.A. Monerat, L.G. Ferreira Filho, E.V. Corrêa Silva, G. Oliveira-Neto, P.H.A.S. Nogueira e A.R.P. de Assumpção, Rev. Bras. de Ens. de Fis. 32, 1304 (2010).
  • [17]
    G.A. Monerat, E.V. Corrêa Silva, L.B. Leal e G. Oliveira-Neto, Rev. Bras. de Ens. de Fis. 37, 4301 (2015).

Apêndice A O Programa CHHAJLANY

Implementamos o MCM, descrito na seção 3., em um programa em linguagem Fortan, aqui denominado CHHAJLANY.FOR, cujo código fonte encontra-se disponível no apêndice deste artigo. Ele permite também a construção de pacotes de onda e do cálculo do valor esperado da posição.

Em relação à avaliação da precisão dos resultados, foram criadas duas ferramentas, tendo por base a verificação da conservação da norma da função de onda e a inspeção de como o potencial U(x) (10) se aproxima do potencial original V(x) (9), comportamento este determinado pela variável auxiliar c.

Em relação a c, observamos que seu valor diminui à medida que aumentamos o grau de truncamento Ω do polinômio v(x) (17), garantindo assim que U(x)V(x).

Para calcularmos o valor de c, usamos a rotina zplrc(coef, zero, ndeg) da biblioteca IMSL embarcada no compilador FORTRAN Profissional. O parâmetro coef é um vetor com dimensão ndeg+1, com ndeg=3, ou seja, coef contém os quatro coeficientes da Eq.(27). O parâmetro zero é um vetor de dimensão ndeg contendo as raízes da referida equação. Na sequência o programa seleciona a raiz real e positiva resultante e atribui este valor à variável auxiliar c. Esse resultado é mostrado na Fig. 5.

Os auto-valores e auto-vetores da matriz dada pela Eq.(21) são obtidos através da subrotina evcrg(a, eval, evec), onde eval é um vetor de dimensão (Ω+1)2 contendo os auto-valores e evec é uma matriz de dimensão (Ω+1)2×(Ω+1)2 contendo na n-esima coluna os coeficientes D correspondentes à expansão Eq.(17) do auto-estado n. O parâmetro SSE é a soma dos quadrados dos erros.

Além disso, o programa também faz a evolução de um pacote de onda em particular. Não comentaremos aqui esta sua funcionalidade.

Apêndice B Código fonte CHHAJLANY.FOR


Datas de Publicação

  • Publicação nesta coleção
    2017

Histórico

  • Recebido
    20 Maio 2016
  • Revisado
    19 Ago 2016
  • Aceito
    22 Ago 2016
Sociedade Brasileira de Física Caixa Postal 66328, 05389-970 São Paulo SP - Brazil - São Paulo - SP - Brazil
E-mail: marcio@sbfisica.org.br