A Manybody Formalism for Fermions, Enforcing the Pauli Principle on Paper
aa r X i v : . [ qu a n t - ph ] J un A Manybody Formalism for Fermions, Enforcing the Pauli Principle on Paper
D.K. WatsonUniversity of OklahomaHomer L. Dodge Department of Physics and AstronomyNorman, OK 73019 (Dated: September 24, 2018)Confined quantum systems involving N identical interacting fermions are found in many areas ofphysics, including condensed matter, atomic, nuclear and chemical physics. In a previous series ofpapers, a manybody perturbation method that is applicable to both weakly and strongly-interactingsystems of bosons has been set forth by the author and coworkers. A symmetry invariant pertur-bation theory was developed which uses group theory coupled with the dimension of space as theperturbation parameter to obtain an analytic correlated wave function through first order for asystem under spherical confinement with a general two-body interaction. In the present paper, weextend this formalism to large systems of fermions, circumventing the numerical demands of ap-plying the Pauli principle by enforcing the Pauli principle on paper. The method does not scale incomplexity with N and has minimal numerical cost. We apply the method to a unitary Fermi gasand compare to recent Monte Carlo values. PACS numbers: 03.65.Ge,03.75.Ss,31.15.xh,31.15.xp
Introduction. – Confined quantum systems of fermionsare widespread across physics. They include, for ex-ample, atoms, atomic nuclei, neutron stars, quantumdots, and cold Fermi gases. These systems possess froma few tens to millions of particles and span an enor-mous range of interparticle interaction strength present-ing a challenge for N -body methods when mean-field ap-proaches fail. In the last decade, ultracold Fermi gaseshave emerged as a testing ground for manybody methodsdue to their precise controllability in experiments. Usinga magnetic field to tune the scattering length of atomsallows the exploration of the physics over many lengthscales, including the “unitary gas” defined by an infinitescattering length. Without a defining length scale otherthan the interparticle distance, the unitary gas does notyield to conventional perturbation treatments, typicallyrequiring a full manybody treatment for an accurate de-scription.Fermi systems in the unitary regime are currently ofgreat interest. This strongly interacting regime, stabi-lized by the Pauli exclusion principle, exists on the cuspof the BCS-BEC crossover and exhibits universal ther-modynamic behavior which has been verified in the lab-oratory to within a few percent[1]. Verifying this univer-sal behavior theoretically requires access to the partitionfunction, i.e. the energy spectrum with degeneracies. Re-quired resources on a classical computer for an exact so-lution of even the ground state of the N -body problemscale exponentially with N . Traditional full configurationinteraction can cope with at most some ten particles[2]and methods such as coupled-cluster[3] with a compu-tational time polynomial in N , O ( N ), are extremelyexpensive. Quantum Monte Carlo (QMC) methods forfermions suffer from the so-called “sign” problem[4–6]resulting in an exponential growth in simulation times.Other methods that have been applied to strongly inter-acting Fermi systems include the method of correlated basis functions[7, 8], density functional theory[9], dia-grammatic approaches[10], and the stochastic variationalmethod [11, 12]. To date, the determination of the energyspectrum of systems containing four or more particles re-mains a challenge.[12]In this paper, we investigate an alternative approach tostudy large systems of fermions. Our symmetry-invariantperturbation theory, SPT, offers a perturbation approachfor a systematic study of correlation including the uni-tary regime. The perturbation parameter depends on theinverse dimensionality of space ( δ = 1 /D ), rather thanthe strength of the interaction, so SPT is equally ap-plicable to weakly or strongly interacting systems. Themethod is essentially analytic [13], with N , the numberof particles, entering as a parameter, allowing results forany N to be obtained from a single calculation[14]. Thelowest-order result includes correlation and in theory, canbe systematically improved by going to higher order[15].This method couples group theory with conventional di-mensional perturbation theory (DPT)[16–18] to take ad-vantage of the high degree of symmetry possible amongidentical particles in higher dimensions. The terms in theperturbation expansion are invariant under N ! symmetryoperations, greatly reducing the number and complex-ity of the building blocks of these terms. This methodcurrently includes full manybody effects that are exactthrough first order. Excited states are obtained from thesame analytic calculation[14], differing only in the num-ber of quanta in the different normal modes. In principle,the full energy spectrum is accessible.In past work, DPT has been applied to fermion sys-tems including single atoms[17, 19], quantum dots[20],and small molecules[21]. The symmetry invariant ap-proach, SPT, has been applied thus far only to bosonswith spherical confinement, determining first-order ener-gies, normal mode frequencies[13, 14], the lowest-orderSPT wave function[22] and density profile[23]. We ap-plied these results to a BEC for which the density profileis a directly observable manifestation of the quantizedbehavior. In a later series of papers, we extended thiswork for bosons to first-order wave functions and densityprofiles[24, 25]. We demonstrated that this method effec-tively rearranges the numerical work for this manybodyproblem into analytic building blocks at each order giv-ing the exact result order by order in the perturbationseries[25, 26]. The complexity of the rearranged problemscales with the order of the perturbation series, not withthe number of particles[26].In this extension of SPT to fermions, we tackle thechallenge of applying the Pauli principle, typically an ex-pensive numerical task. We describe how this is achieved“on paper” for any value of N , thus circumventing heavynumerical effort. We apply the method in this initialstudy to large systems of cold fermions in the unitaryregime for which a number of very accurate calculationsare available. Results are obtained using analytic build-ing blocks which have been calculated and stored pre-viously and which have been extensively checked usingan independent solution of a model system of harmon-ically confined, harmonically interacting particles[25].Our results compare well to very accurate Monte Carloresults[8, 27] including some recent benchmark calcula-tions using the auxilliary field Monte Carlo method[28].We begin the perturbation analysis by defining di-mensionally scaled quantities: ¯ E = κ ( D ) E, and ¯ H = κ ( D ) H where κ ( D ) is a scale factor which regularizesthe large-dimension limit[22]. The scaled version of theSchrodinger equation becomes¯ H Φ = (cid:18) κ ( D ) ¯ T + ¯ V eff (cid:19) Φ = ¯ E Φ , (1)where barred quantities indicate variables in scaled units( κ ( D ) = D / ( ~ ¯ ω ho ) for this work. See Ref. [22]). Theterm ¯ T contains the derivative terms of the kinetic energyand ¯ V eff includes centrifugal, two particle, and confine-ment potentials[22].We assume a totally symmetric, large-dimension con-figuration at which the effective potential is a minimum.The N particles are arranged on a hypersphere, each par-ticle with a radius, ¯ r ∞ , from the center of the confiningpotential. Furthermore, the angle cosines between eachpair of particles take on the same value, γ ∞ , i.e.lim D →∞ ¯ r i = ¯ r ∞ (1 ≤ i ≤ N ) , lim D →∞ γ ij = γ ∞ (1 ≤ i < j ≤ N ) . (2)(This symmetric high-dimensional structure is not unlikethe localized structure found in a hyperspherical treat-ment of the confined two-component Fermi gas in the N → ∞ limit[29].) In scaled units the δ → D → ∞ )approximation for the energy is simply the effective po-tential minimum, i.e. ¯ E ∞ = ¯ V eff (¯ r ∞ , γ ∞ ; δ = 0) .In this δ → V eff , which is nonzero even for the ground state, is a zero-point energy contribution satis-fying the minimum uncertainty principle[30]. The valueof γ ∞ , which is zero in the mean-field approximation forthe L = 0 angular momentum states considered here, is,in fact, not zero, an indication that beyond-mean-fieldeffects are included in the δ → δ → N identical objects[31], S N ,allowing a largely analytic solution. The δ → S N symmetry greatly simplifies this task since the inter-action terms individually have to transform as a scalarunder the S N point group.The perturbation series has the form:¯ E = ¯ E ∞ + δ ∞ X j =0 (cid:16) δ (cid:17) j ¯ E j Φ = ∞ X j =0 (cid:16) δ (cid:17) j Φ j . (3)In practice ¯ E j = 0 ∀ j odd. The j = 0 terms are ob-tained from a harmonic equation, and referred to as theenergy and wave function at harmonic order. To obtainthis harmonic correction for small values of δ , we expandabout the minimum of the δ → ω µ . The number of roots, λ µ , ( λ µ = ¯ ω µ ) , ofthe secular equation, N ( N + 1) / S N symmetry of the problemthere is a reduction to five distinct roots,The F G matrix is invariant under S N , so it does notconnect subspaces belonging to different irreducible rep-resentations (irreps.) of S N [33]. Thus the normal co-ordinates must transform under irreps. of S N . Thenormal coordinates are linear combinations of the ele-ments of the internal displacement vectors which trans-form under reducible representations of S N . One canshow that these reduce to two 1-dimensional [ N ] irreps.denoted by + , − , two ( N − N − , + , − , and one angular N ( N − / N − ,
2] irrep. denoted by [13].The energy through harmonic order in δ is [13] E = E ∞ + δ " X µ = { ± , ± , } ( n µ + 12 d µ )¯ ω µ + v o , (4)where n µ is the total number of quanta in the normalmode with the frequency ¯ ω µ ; µ is a label which runsover − , + , − , + , and , regardless of the numberof particles in the system (see Refs. [13]and Ref.[15] in[22]), and v o is a constant. The multiplicities of the fiveroots are: d + = 1 , d − = 1 , d + = N − , d − = N − , d = N ( N − / normal modes are phonon, i.e. compressionalmodes; the ± modes show single-particle character, andthe ± modes describe center-of-mass and breathing mo-tions. Enforcing the Pauli Principle. – To generalize SPTfrom quantum systems of bosons to quantum systemsof fermions, we must enforce the Pauli principle, thus re-quiring the N -body wave function to be totally antisym-metric. This is enforced by placing certain restrictions onthe occupancies of the normal modes, i.e. on the values ofthe normal mode quantum numbers, n µ , µ = ± , ± , inEq. (4)[17]. The possible assignments can be found by re-lating the normal mode states | n + , n − , n + , n − , n > to the states of the confining potential which is a spher-ically symmetric three dimensional harmonic oscillator( V conf ( r i ) = mω ho r i ) for which the restrictions im-posed by antisymmetry are known. These two seriesof states can be related in the double limit D → ∞ , ω ho → ∞ where both representations are valid.For large D , the normal mode description given byEq. (4) is exact. Applying the large ω ho limit results in: E = N D ~ ω ho + (2 n + + 2 n − + 2 n + + 2 n − + 2 n ) ~ ω ho (5)Now consider ω ho → ∞ first and then D → ∞ . Theharmonic oscillator levels are exact: E = N X i =1 (cid:20) (2 ν i + l i ) + D (cid:21) ~ ω ho = N D ~ ω ho + N X i =1 (2 ν i + l i ) ~ ω ho (6)where ν i is a radial quantum number and l i the orbitalangular momentum quantum number. Equating thesetwo expressions which are equal in the double limit, thequantum numbers in the two representations can now berelated to show the restrictions on normal mode statesimposed by antisymmetry. Because of the clean separa-tion of radial and angular motions, two conditions result:2 n − + 2 n − = N X i =1 ν i n + + 2 n + + 2 n = N X i =1 l i . (7)These equations determine a set of possible normal modestates | n + , n − , n + , n − , n > from the known set ofpermissible L = 0 harmonic oscillator configurations. Application: The Unitary Gas. – The Schr¨odingerequation for an N -body system of fermions, N = N + N with N spin up and N spin down fermions, confined by FIG. 1: Ground state energies of the harmonically trappedunitary Fermi gas (units ~ ω ho = 1). Our first-order per-turbation results (filled diamonds) are compared to GFMC(+’s)[27], fixed-node DMC (open circles) from Ref. [8], andAFMC results (filled circles, N =6, 8, 14, 20, 30 only)[28]. a spherically symmetric potential is H Ψ = N X i =1 h i + N X i =1 N X j =1 g ij Ψ = E Ψ . (8)where h i and g ij are the single-particle Hamiltonian andthe two-body interaction potential, respectively. We as-sume a T = 0 K condensate with N = N confined byan isotropic, harmonic trap with frequency ω ho .To study the unitary regime, we replace the actualatom-atom potential by an attractive square well poten-tial of radius R : V int ( r ij ) = (cid:26) − V o , r ij < R , r ij ≥ R . , (9)For fixed range R , the potential depth V is adjusted sothe s-wave scattering length, a s is infinite. The rangeis selected so R << a ho ( a ho = p ~ / ( mω ho )) and canbe systematically reduced to extrapolate to zero-rangeinteraction. We dimensionally continue the square wellpotential so that it is differentiable away from D = 3 ,allowing us to perform the dimensional perturbationanalysis[13, 14]). Thus, we take the interaction to be V int ( r ij ) = V ( δ ) (cid:20) − tanh (cid:20) − δ ( r ij − δR ) (cid:21)(cid:21) , (10)where V o ( δ ) = − bδ . The potential depth V o is adjustedby adjusting the value of b so the scattering length is infi-nite when δ = . This interaction becomes a square wellof radius R in the physical D = 3 limit. The functionalform of the potential at D = 3 is not unique. Other formscould be chosen with equal success as long as the formis differentiable and reduces to a square well potential at D = 3 . We simply choose a form that allows a gradualsoftening of the square well.In Fig. 1 we plot the ground state energies from N = 6 to N = 30 and compare to Green’s function FIG. 2: Ground state energies of the harmonically trappedunitary Fermi gas (units ~ ω ho = 1). Our first-order perturba-tion results (filled diamonds) are compared to DMC resultsRef. [28] (open squares) Monte Carlo (GFMC) energies[27], to fixed-node diffu-sion Monte Carlo (DMC) energies[8] which provide ac-curate upper bounds to the ground state energy, and torecent benchmark auxiliary field Monte Carlo calcula-tions (AFMC)[28] which are exact, but subject to finitelattice size errors for which corrections have been made.The AFMC results are currently the most accurate re-sults available. Our energies which include full many-body effects through first order compare well to theseAFMC results and except for N = 8 they are slightlycloser to the AFMC results than the GFMC and DMCresults shown in Figure 1. Our numbers, as well as theGFMC and DMC numbers from Refs. [8, 27], show a dis-tinct odd/even oscillation, but no obvious shell effects.Obtained from analytic building blocks that have beencalculated and stored previously, our numerical require-ments take a few seconds on a work station.No parameters are used to obtain these results otherthan those used to produce a square well potential withinfinite scattering length. For the GFMC and DMC stud- ies shown in Figure 1, no results exist for values of N higher than shown. Accurate results for higher N are,of course, increasingly difficult to achieve. In Figure 2we compare our results at first order above N = 30 withavailable DMC results from Ref. [28]. Our first-order re-sults show the expected increase in error as N increasessuggesting that for larger N, higher order terms may benecessary. Conclusions – In this paper we have extended thesymmetry-invariant perturbation method from bosonsto fermions, applying the Pauli principle “on paper” toavoid heavy numerical expense. The method has beentested in the unitary regime, which is of particular in-terest for manybody methods since its infinite scatteringlength and the lack of a natural scale typically requireintensive numerical simulation for an accurate descrip-tion. Our analytic results through first order yield ener-gies that are comparable in accuracy with recent MonteCarlo results. As N increases, our error increases sug-gesting the need for higher order terms. It may also bepossible to rearrange the perturbation series to minimizethe importance of higher order terms.The theory applied in this paper is applicable to L = 0states of spherically confined systems with general at-tractive or repulsive interparticle interactions and is alsoapplicable to both weakly and strongly correlated sys-tems. The fact that γ ∞ is not zero is an indication thatbeyond-mean-field effects are included in this result evenin the D → ∞ limit. This theory is readily generalizableto systems with a cylindrical confining potential.Higher-order calculations may be required for larger N . A detailed program for calculating higher-order DPTcorrections to N -body systems has been laid out[15] andapplied to high order for small- N systems[34]. For large- N systems the S N point-group symmetry greatly simpli-fies the calculation of these higher-order terms.We thank Doerte Blume and Stefano Gandolfi for dis-cussions and Monte Carlo results. [1] H. Hu, P.D. Drummond, and X.-J. Liu, Nat. Phys. ,469, (2007); M. Horikoshi et al. Science , 442(2010);S. Nascimbene et al. Nature(London( , 1057(2010);M.J.H. Ku et al., Science , 563(2012); C. Cao et al., ibid , 58(2011).[2] G.H. Booth et al., Nature . 365 (2013).[3] See for example L.S. Cederbaum, O.E. Alon, and A.I.Streltsov, Phys. Rev. A , 043609 (2006); R.J. Bartlettand M. Musial, Rev. Mod. Phys. , 291 (2007).[4] See for example D. Landau and K. Binder, A Guide toMonte-Carlo Simulations in Statistical Physics , (Cam-bridge University Press, Cambridge, 2001); M. Holz-mann, W. Krauth, and M. Naraschewski, Phys. Rev. A , 2956 (1999).[5] K.E. Schmidt and M.H. Kalos, in K. Binder(ed.): Appli-cations of the Monte Carlo Method in Statisitcal Physics (Springer, 1984).[6] M. M. Forbes, S. Gandolfi and A. Gezerlis, Phys. Rev. A , 053603 (2012).[7] See for example S. Fantoni and A. Fabrocini, in it Micro-scopic Quantum Many Body Theories and Their Appli-cations, J. Navarro and A. Polls eds., Lecture Notes inPhysics Vol. 150 (Springer-Verlag, Berlin, 1998) p. 119.[8] D. Blume,J. von Stecher, and C. H. Greene, Phys. Rev.Lett. , 233201(2007).[9] W. Kohn, Rev. Mod. Phys. , 1253 (1999).[10] H. Hu, X.-J. Liu, and P.D. Drummond, Europhys. Lett. , 574(2006).[11] K. Varga, Y. Suzuki, Phys. Rev. C. , 2885(1995).[12] K.M. Daily and D. Blume, Phys. Rev. A , 053615(2010).[13] B.A. McKinney et al., Ann. Phys. (NY) , 56 (2003).[14] B.A. McKinney, M. Dunn, D.K. Watson, Phys. Rev. A , 053611 (2004).[15] M. Dunn et al., J. Chem. Phys. Dimensional Scaling in Chemical Physics , edited by D.R.
Herschbach, J. Avery, and O. Goscinski (Kluwer, Dor-drecht, 1992).[17] J.G. Loeser, J. Chem. Phys. , 5635 (1987).[18] A. Chatterjee, Phys. Rep. , 249 (1990).[19] L.D. Mlodinow and N. Papanicolaou, Ann. Phys.(N.Y.) ,2299(1981).[20] B.A. McKinney and D.K. Watson, Phys. Rev. B , 4958(2000).[21] S.-W. Huang et al., Phys. Rev. A , 250 (1998); A. A.Svidzinsky, M. O. Scully, and D. R. Herschbach, Phys.Rev. Lett. , 080401, 2005.[22] M. Dunn, D.K. Watson, and J.G. Loeser, Ann. Phys.(NY), , 1939 (2006).[23] W.B. Laing, M. Dunn, and D.K. Watson, Phys. Rev. A , 063605 (2006).[24] W.B. Laing, M. Dunn, and D.K. Watson, J. of Math.Phys. , 062105 (2009).[25] W.B. Laing et al., J Phys A , 205307 (2009); M. Dunnet al., Phys. Rev A , 062108 (2009). [26] D.K. Watson and M. Dunn, Phys. Rev. Lett. , 020402(2010).[27] S.Y. Chang and G.F. Bertsch, Phys. Rev. A , 021603(2007).[28] J. Carlson and S. Gandolfi, Phys. Rev. A , 735 (1985).[31] See for example M. Hamermesh, Group Theory andits Application to Physical Problems , (Addison-Wesley,Reading, MA, 1962).[32] E.B. Wilson, Jr., J.C. Decius, P.C. Cross,
Molecular vi-brations: The theory of infrared and raman vibrationalspectra . McGraw- Hill, New York, 1955.[33] See for example reference 32, Appendix XII, p. 347.[34] J.R. Walkup et al., Phys. Rev. A, , 4668, (1998); J.C.Carzoli, M. Dunn and D.K. Watson, Phys. Rev. A ,182 (1999); J.R. Walkup, M. Dunn, and D.K. Watson, J.Math. Phys. , 218 (2000); J.R. Walkup, M. Dunn, andD.K. Watson, Phys. Rev. A63