Resonantly Interacting Fermions In a Box
Michael McNeil Forbes, Stefano Gandolfi, Alexandros Gezerlis
lla - ur 10 - , int - pub - - , nt @ uw - - Resonantly Interacting Fermions In a Box
Michael McNeil Forbes, , , Stefano Gandolfi, and Alexandros Gezerlis , Institute for Nuclear Theory, University of Washington, Seattle, Washington – USA Department of Physics, University of Washington, Seattle, Washington – USA and Theoretical Division, Los Alamos National Laboratory, Los Alamos, New Mexico , USA (Dated: October , )We use two fundamental theoretical frameworks to study the finite-size (shell) properties of the unitarygas in a periodic box: ) an ab initio Quantum Monte Carlo ( qmc ) calculation for boxes containing to particles provides a precise and complete characterization of the finite-size behavior, and ) a newDensity Functional Theory ( dft ) fully encapsulates these effects. The dft predicts vanishing shellstructure for systems comprising more than particles, and allows us to extrapolate the qmc resultsto the thermodynamic limit, providing the tightest bound to date on the ground-state energy of theunitary gas: ξ S (cid:54) . ( ) . We also apply the new functional to few-particle harmonically trappedsystems, comparing with previous calculations. PACS numbers: . .-d, . .Mb, . .E-, . .Ss, . .Cn, . .Hh, . .-n T he fermion many - body problem plays a fundamen-tal role in a vast array of physical systems, from dilutegases of cold atoms to nuclear physics in nuclei andneutron stars. The universal character of this problem– each system is governed by a similar microscopic the-ory – coupled with direct experimental access in coldatoms, has led to an explosion of recent interest (seeRefs. [ ] for a review). Despite this broad applicability,we are far from fully understanding even the simplestsystem: the “unitary gas” comprising equal numbers oftwo fermionic species interacting with a resonant s -waveinteraction of infinite scattering length a s → ∞ . Lack-ing any scale beyond the total density n + = n a + n b ,the unitary gas eschews perturbative expansion and re-quires experimental measurement or accurate numericalsimulation for a quantitative description – the latter ispresently more precise. Typical Quantum Monte Carlo( qmc ) calculations, however, can access at most a fewhundred particles. Density Functional Theory ( dft ) pro-vides a complementary approach through which onemay extrapolate these results to large systems beyondthe reach of direct simulation.In this Letter, we present the most precise qmc cal-culations to date of the unitary gas in a periodic box,studying from to particles, thereby providing abenchmark for many-body theories. We use this to cali-brate a local dft , then use this dft to study the finite-sizeeffects (“shell” effects in nuclear physics) and extrapolatethe qmc results to the thermodynamic limit. We providethe most precise bound to date of the universal Bertschparameter [ ] ξ S = E / E FG (cid:54) . ( ) . ( E FG = + E F is the energy density of a free Fermi gas with the sametotal density n + = k / ( ) , and E F = (cid:32) h k / ( ) is theFermi energy.) We also explore the finite-size propertiesof the dft – a crucial element in the program to calculateproperties of finite nuclei with a universal dft [ ]. Wefind that a local dft can capture the finite-size effects in these systems without the need for particle number pro-jection. We limit our discussion to symmetric systems( n a = n b ), leaving odd-even staggering to future work,as the dft then requires an additional dimensionlessparameter to characterize the asymmetry n a (cid:54) = n b .The qmc results presented here are directly applicableto cold Li or K atoms, and constrain dilute neutronmatter in neutron stars [ ]; likewise, the general dft approach has myriads of applications throughout cold-atom and nuclear physics (see Ref. [ ] for a review).Our calculation of ξ S is consistent with previous results,but an order of magnitude more precise. Continuum qmc bounds ξ S (cid:46) . – . with an uncertainty nobetter than the last digit [ , – ]. Lattice qmc resultsrange from ξ S ≈ . – . [ – ], comparable to ana-lytic results [ ]. Experimental groups found qualitativeagreement [ ], which led to precision measurements:notably with Duke [ ] and Paris [ ] quoting . ( ) and . ( ) respectively.D ft is an in principle exact approach, widely usedin in nuclear physics [ ], and in quantum chemistryto describe normal (i.e., non-superfluid) systems. It hasrecently been extended to describe the unitary gas [ , – ]. We build upon one approach – the Superfluid LocalDensity Approximation ( slda ) – which was originallyconstrained by qmc calculations of the continuum state,and then validated with qmc calculations in a harmonictrap [ , ] (see also Fig. ). We focus on translationallyinvariant systems in a periodic box to isolate the finite-size effects from the gradient effects. We find that theinclusion of an anomalous density is crucial: functionalsattempting to model the superfluid by adding only gra-dient or kinetic corrections [ , ] are unable to evenqualitatively characterize the finite-size effects.Our qmc results are based on a fixed-node DiffusionMonte Carlo approach that projects out the state of low-est energy from the space of all wave functions with a r X i v : . [ c ond - m a t . qu a n t - g a s ] J un (= N + ) ( ξ =) 0 0.1 0.2 (= k F r e )0.380.400.42 ξ ( k F r e ) for N + = Figure . (color online) Ground-state energy-density ξ = E / E FG of N + fermions in a periodic cubic box at the uni-tary limit. The circles with error bars are the result of using aquadratic least-squares extrapolation to zero effective range ofour new qmc results. The solid curve is the best fit slda dft .The light dotted curve is the functional considered in [ ] with α = . . For comparison, we have plotted the previous bestestimate ξ S = . ( ) (red square) and the current estimate ξ S = . ( ) below it to the far right of the figure. Inset: weshow the typical effective-range dependence ξ ( k F r e ) with thebest fit σ error bounds for all-point cubic (solid dark green)and five-point quadratic (hatched light yellow) polynomialfits. Note that: a) the five-point quadratic model is consistentwith the full cubic model and has a comparable extrapolationerror, and b) the inflection point near k F r e ≈ . necessitatesa higher-order fit for larger ranges (cubic is sufficient for theranges shown here). Results for N + = show the samequalitative behaviour; hence, for the other points we use thefive-point quadratic extrapolation. fixed nodal structure as defined by an initial many-bodywave function (ansatz). By varying the ansatz, we obtaina variational upper bound on the ground-state energy.In this work, we use the trial function introduced in [ ]: Ψ T = A [ φ ( r (cid:48) ) φ ( r (cid:48) ) · · · φ ( r nn (cid:48) )] (cid:89) ij (cid:48) f ( r ij (cid:48) ) , ( )where A antisymmetrizes over particles of the same spin(either primed or unprimed) and f ( r ) is a nodeless Jas-trow function introduced to reduce the statistical error.The antisymmetrized product of s -wave pairing func-tions φ ( r ij (cid:48) ) defines the nodal structure: φ ( r ) = (cid:88) n α (cid:107) n (cid:107) e i k n · r + ˜ β ( r ) . ( )The sum is truncated (we include ten coefficients) andthe omitted short-range tail is modelled by the phe-nomenological function ˜ β ( r ) chosen to ensure smoothbehavior near zero separation. We use the same formfor ˜ β ( r ) as in [ ] with the values b = . and c = , andvary the coefficients α (cid:107) n (cid:107) for each N + to minimize theenergy as described in Ref. [ ]. Representative nodal N + ξ ( N + ) ξ box ( N + ) ξ N /ξ box . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . . ( ) . ( ) . Table I . Values of ξ ( N + ) = E ( N + ) /V/ E FG shown in Fig. . Tofacilitate comparison with other normalizations in the litera-ture, we include the values ξ box ( N + ) = E ( N + ) /E FG ( N + ) where E FG ( N + ) is the energy of N + non-interacting particles in a box.The conversion factor is shown in the last column. structures are defined by the coefficients in Table II. Thesame ansatz suffices for different effective ranges, butan independent optimization is required for each N + .We simulate the Hamiltonian: H = (cid:32) h (cid:32) − N + (cid:88) k = ∇ − µ (cid:88) i , j (cid:48) sech ( µr ij (cid:48) ) (cid:33) , ( )with an interspecies interaction of the modified P ¨oschl-Teller type (off-resonance intraspecies interactions areneglected). We tune to infinite s -wave scattering lengthby setting v = : the effective range becomes r e = .To extrapolate to the zero-range limit, we simulate at µ/k F ∈ { . , , , , } for which .
03 < k F r e < 0 . .A careful examination of additional ranges up to k F r e ∼ . for N + = and N + = (see the inset in Fig. ) N + a a a a a a a a a a
10 1600 350 49 16 12 14 14 11 9 . .
740 160 91 27 0 . - . - .
086 2 . . . . -
24 13 12 8 . . . . . . . - -
17 0 .
51 7 . . . . . . . Table II . Sample coefficients of the pairing function ( ) α (cid:107) n (cid:107) = − a I where I = (cid:107) n (cid:107) = n + n + n = k L /4π . Higher-order coefficients are set to zero. reveals that a three-parameter quadratic model in r e is necessary and sufficient to extrapolate our resultswithout a systematic bias; the results are shown in Fig. .The energies exhibit definite finite-size effects for N + (cid:46) , but are essentially featureless for larger N + .This lack of structure is confirmed by the best fit dft (dis-cussed below) and disagrees with the results presentedin Ref. [ ]. The values of ξ for N + > 50 are distributedabout the best fit value ξ S ≈ . ( ) , and represent thelowest variational bounds to date. Part of the decreasefrom previous results is due to the careful extrapolationto zero effective range. The remainder is due to theimproved optimization of the variational wave function.To model the finite-size effects we turn to a local dft for the unitary Fermi gas that generalizes the slda originally presented in Ref. [ ]. In addition tothe total density n + = (cid:80) n | v n | , the slda includesboth kinetic τ + = (cid:80) n | ∇ v n | and anomalous densities ν = (cid:80) n u n v ∗ n . (The + index signifies the sum of thecontributions coming from the two components a and b ; u n ( r ) and v n ( r ) are the Bogoliubov quasiparticle wavefunctions.) The original three-parameter slda is E slda = (cid:32) h m (cid:18) α2 τ + + β 310 ( ) n + (cid:19) + gν † ν , ( )where α is the inverse effective mass; β is the self-energy;and γ controls the pairing through the regularized cou-pling g = ( n + /γ − Λ/α ) where Λ → ∞ is a momen-tum cutoff that we take to infinity (see Ref. [ ] for details).One can use numerically the equations for homogeneousmatter in the thermodynamic to replace the parameters β and γ with the more physically relevant quantities ξ S and η = ∆/E F , where ∆ is the pairing gap.In principle, the dft can be expressed in terms ofonly the density n + and its gradients. References [ ]consider local formulations of this type (called Ex-tended Thomas-Fermi ( etf ) functionals). Since gradientsvanish in the periodic box, etf functionals reduce to E etf ( n + ) ≡ ξ S E FG and exhibit no finite-size structure,contrary to the qmc results. Reference [ ] adds ατ + ,but without ν † ν the finite-size effects do not correlatewith the qmc behavior (see Fig. ) and the best fit toour results is also flat ( α → ). Furthermore, such func-tionals cannot qualitatively reproduce the quasiparticledispersion relationship, an attractive feature of the slda (see also Ref. [ ]).The best fit three-parameter slda functional ( ) – α = . ( ) , ξ S = . ( ) , and η = . ( ) – is shown inFig. . It fits the
23 qmc points from N + = to N + = with a reduced chi squared χ r = . , indicatingcomplete consistency. Although remarkable, the fit isnot completely satisfactory: ) Fitting the exact two-particle energy ξ = − . · · · raises χ r = . , and ) the best fit gap parameter η and inverse effective ( = x = ( + ) − ) E / ( N + ) /
50 30 20 10 6 4 ( N + =) Figure . (color online) Ground-state energy of the harmon-ically trapped unitary Fermi gas (in units where (cid:32) hω = )scaled to demonstrate the asymptotic form / ( + ) = ξ S (cid:0) + cx + O ( x ) (cid:1) predicted by the low-energy effectivetheory of Ref. [ ]. The best fit slda (solid blue line) is com-pared with zero-range results for N + ∈ { , } from Ref. [ ],and finite-range qmc results from Ref. [ ] (upper red dots)and Ref. [ ] (green pluses). The latter have significantly lowerenergy, despite having a slightly large effective range, sug-gesting that the wave functions in Ref. [ ] were not fullyoptimized. We expect careful optimization and zero-rangeextrapolation to bring the qmc results for large N + in line withthe dft as discussed in the text. mass α are inconsistent with the values η = . ( ) and α = . ( ) obtained from the N + = qmc quasiparticledispersion relation [ , ], and the values η = . ( ) [ ]and η = . ( ) [ ] extracted from experimental data.These deficiencies might be remedied by generalizingthe slda . As noted in Ref. [ ], the following combinationof divergent kinetic and anomalous densities is finite: K = (cid:32) h τ + + gα ν † ν = (cid:32) h τ + + ν † ναn + /γ − Λ . ( )The lack of scales thus dictates the functional form: E ( K , n + ) = ξ ( Q ) E FG ( n + ) , Q = K/ E FG ( n + ) , ( )where Q is dimensionless, and the regularization con-dition depends on Q through the function γ ( Q ) . Theoriginal slda is linear ξ ( Q ) = αQ + β with constant γ ( Q ) = γ . This generalized functional can fit any mono-tonic ξ ( N + ) , including the exact N + = point. For N + > 6 , ξ ( N + ) is not monotonic and the functional isin principle constrained. For example, requiring that ξ = ξ S at both N + ≈ . ( ) and N + = ∞ fixes the ratio η/α = . ( ) . (As an aside, we note that the momentumdistribution n k in the dft relates this to the “contact” C : η/α = √ ≈ . – . ; see Refs. [ ] and referencestherein, though it is not clear that this property shouldbe trusted.) In practice, the errors and the discretenessin N + leave room for flexibility in the functional form,and we have found several generalized functional formswith χ ≈ . while constraining η = . . We may haveto accept the discrepancy in α as a limitation of the dft .However, generalizing the slda may not be needed:analyzing the “symmetric heavy-light ansatz” [ ], (jus-tified by lattice qmc calculations [ ]), we find that the simple three-parameter slda suffices ( χ ≈ . ) withreasonable α = . ( ) , η = . ( ) , and ξ = . ( ) –slightly higher than the ξ = . ( ) extracted in [ ].It is not trivial that the simple dft ( ) captures allfinite-size effects above N + = to high precision in bothcalculations, indicating that the slda may be used toextrapolate to the thermodynamic limit. We note thatno particle-number projection is required – a quite ill-defined procedure often considered necessary in nuclearphysics [ ]: Perhaps improved nuclear functionals maysimilarly capture finite-size effects through local anoma-lous densities in the spirit of ν .To finish, we consider harmonically trapped systemsin Fig. . As discussed in [ ], the energy may be ex-pressed as E ( N + ) = (cid:32) hω √ ξ S ( + ) (cid:0) + cx + O ( x ) (cid:1) where x = ( + ) − and c is expressed in terms of low-energy coefficients. As demonstrated by the zero-range N + ∈ { , } results of [ ], the dft still over-estimates theenergy for small systems, most likely because we haveomitted the gradient terms in the functional that vanishin homogeneous systems.For large N + the dft has the expected asymptotic formwith intercept ξ S = . unlike the finite-range qmc results of Refs. [ , ]. This is qualitatively consistentwith the leading effective-range corrections which scaleasymptotically as x − ; the systematic overestimationof the energy by the variational qmc approach might alsocontribute. We defer further discussion until carefullyextrapolated zero-range results are published.To summarize, we present the most precise QuantumMonte Carlo calculations to date of a symmetric unitaryFermi gas in a periodic box comprising to particles.By carefully characterizing and extrapolating these re-sults to zero effective range, we have completely mappedout the finite-size effects. These results are used to ana-lyze the structure of a Density Functional Theory for thesymmetric unitary gas, and it is shown that the simplestthree-parameter form of Eq. ( ) fully accounts for allshell effects to within the statistical errors of the qmc results without the need for particle-number projection;a more complicated form, however, may be requiredto capture both the finite-size effects and the quasipar-ticle dispersions. The dft predicts no significant shellcorrections beyond particles, and the qmc calcula-tions confirm this, allowing us to extract a precise upperbound on the universal equation of state ξ S (cid:54) . ( ) ,an order of magnitude improvement in precision overprevious bounds and the lowest bound of any varia-tional method to date. The functional in its latest formis well constrained, but leads to slight disagreementswith qmc predictions for harmonic traps. Convergingboth qmc and dft approaches promises to be a fruitfuldirection of future research.We thank Aurel Bulgac, Joe Carlson, and Dean Lee for useful discussions. This work is supported, inpart, by us Department of Energy ( d o e ) grants de - fg02 - , de - fg02 - , & de - ac52 - , d o e contracts de - fc02 - ( unedf s ci dac ) & de - ac52 - , and by the ldrd program at LosAlamos National Laboratory ( lanl ). Computations forthis work were carried out through Open Supercom-puting at lanl , on the uw Athena cluster, and at theNational Energy Research Scientific Computing Center( nersc ). [ ] M. Inguscio, W. Ketterle, and C. Salomon, eds., Ultra-cold Fermi Gases , International School of Physics “EnricoFermi”, Vol. (I os Press, Amsterdam, ); S. Giorgini,L. P. Pitaevskii, and S. Stringari, Rev. Mod. Phys. , ( ).[ ] “The Many-Body Challenge Problem ( mbx ) formulatedby G. F. Bertsch in ”; G. A. Baker, Jr., Phys. Rev. C , ( ); Int. J. Mod. Phys. B, , ( ).[ ] Do e s ci dac unedf project, .[ ] A. Gezerlis and J. Carlson, Phys. Rev. C , ( ).[ ] A. Bulgac, M. M. Forbes, and P. Magierski in BCS-BECCrossover and the Unitary Fermi Gas , edited by W. Zwerger,Lecture Notes in Physics (Springer, ).[ ] J. Carlson, S. Y. Chang, V. R. Pandharipande, and K. E.Schmidt, Phys. Rev. Lett. , ( ).[ ] G. E. Astrakharchik, J. Boronat, J. Casulleras, andS. Giorgini, Phys. Rev. Lett. , ( ).[ ] J. Carlson and S. Reddy, Phys. Rev. Lett. , ( ).[ ] A. Gezerlis, S. Gandolfi, K. E. Schmidt, and J. Carlson,Phys. Rev. Lett. , ( ).[ ] A. J. Morris, P. L´opez R´ıos, and R. J. Needs, Phys. Rev. A , ( ).[ ] A. Bulgac, J. E. Drut, and P. Magierski, Phys. Rev. A , ( ).[ ] D. Lee, Phys. Rev. C , ( ); T. Abe and R. Seki,Phys. Rev. C , ( ).[ ] S. Zhang, K. E. Schmidt, and J. Carlson, private com-munication; J. Drut, A. Gezerlis, and T. A. L¨ahde, inpreparation ( ); M. Endres, D. Kaplan, J.-W. Lee, andA. Nicholson., in preparation ( ).[ ] R. Haussmann, W. Rantner, S. Cerrito, and W. Zwerger,Phys. Rev. A , ( ); Y. Nishida, , ( ).[ ] M. Bartenstein et al. , Phys. Rev. Lett. , ( );J. Kinast et al., Science , ( ); G. B. Partridge,W. Li, R. I. Kamar, Y. an Liao, and R. G. Hulet, , ( ).[ ] L. Luo and J. E. Thomas, J. Low Temp. Phys. , ( ).[ ] N. Navon, S. Nascimb`ene, F. Chevy, and C. Salomon,Science , ( ).[ ] J. E. Drut, R. J. Furnstahl, and L. Platter, Prog. Part. Nucl.Phys. , ( ); A. Gezerlis and G. F. Bertsch, Phys.Rev. Lett. , ( ); S. Gandolfi, J. Carlson, andS. C. Pieper, Phys. Rev. Lett. , ( ).[ ] T. Papenbrock, Phys. Rev. A , ( ).[ ] A. Bulgac, Phys. Rev. A , ( ).[ ] G. Rupak and T. Schaefer, Nucl. Phys. A , ( ); L. Salasnich and F. Toigo, Phys. Rev. A , ( ).[ ] S. Y. Chang and G. F. Bertsch, Phys. Rev. A , ( ).[ ] D. Blume, J. von Stecher, and C. H. Greene, Phys. Rev.Lett. , ( ).[ ] S. Sorella, Phys. Rev. B , ( ).[ ] A. Bhattacharyya and R. J. Furnstahl, Phys. Lett. B , ( ).[ ] A. Bulgac and M. M. Forbes, Phys. Rev. Lett. , ( ).[ ] J. Carlson and S. Reddy, Phys. Rev. Lett. , ( ). [ ] A. Schirotzek, Y. Shin, C. H. Schunck, and W. KetterlePhys. Rev. Lett. , ( ).[ ] S. Gandolfi, K. E. Schmidt, and J. Carlson, Phys. Rev. A , ( ); J. E. Drut, T. A. L¨ahde, and T. Ten,arXiv: . .[ ] D. Lee, Eur. Phys. J. A , ( ).[ ] T. Duguet, M. Bender, K. Bennaceur, D. Lacroix, andT. Lesinski, Phys. Rev. C , ( ).[ ] D. T. Son and M. B. Wingate, Ann. Phys. (NY) , ( ).[ ] D. Blume and K. M. Daily, C. R. Phys. , ( ).[ ] S. Gandolfi, J. Carlson, and K. E. Schmidt, private com-munication (2010