Particle Projection Using a Complex Langevin Method
PParticle Projection Using a Complex Langevin Method
Christopher R.
Shill and Joaquín E.
Drut Department of Physics and Astronomy, University of North Carolina, Chapel Hill, NC, 27599, USA
Abstract.
Using complex stochastic quantization, we implement a particle-number pro-jection technique on the partition function of spin-1 / ffi cients of one-dimensional, attractively interacting fermions. Calculating the thermodynamics of any quantum many-body system is a challenging but generally in-teresting problem for a variety of applications. By far the main roadblock in Monte Carlo calculationsof such systems is the sign problem, whose impact a ff ects fields ranging from ultracold atoms [1, 2] toQCD [3]. In the last few years, a concerted e ff ort to tackle this issue has been taking place, as severaltechniques have been put forward to that e ff ect (see e.g. [4]). In this work, we make use of one of thoseideas, namely that of complex stochastic quantization (i.e. we use the complex Langevin approach)to explore the feasibility of particle-number projection in the calculation of the thermal properties offinite clusters and in the numerical determination of virial coe ffi cients.We begin by recalling that the equilibrium thermodynamics is encoded in the grand canonicalpartition function, Z = tr (cid:104) e − β ( ˆ H − µ ˆ N ) (cid:105) , (1)where ˆ H is the Hamiltonian, ˆ N is the particle number operator, β is the inverse temperature, and µ isthe chemical potential. This in itself is nearly impossible to solve for essentially all interacting cases.To gain a better understanding, in particular in dilute regimes, we may expand about z = e βµ = z is the fugacity. Such a low-fugacity expansion is what we call the virial expansion. Thecorresponding coe ffi cients accompanying the powers of z in the expansion of the grand-canonicalpartition function Z are the N -particle partition functions; specifically, for spin-1 / Z = ∞ (cid:88) n ↑ , n ↓ = z n ↑ ↑ z n ↓ ↓ Q n ↑ , n ↓ (2)where z ↑ = e βµ ↑ and z ↓ = e βµ ↓ are the fugacity for the spin-up and spin-down particles, respectively,and Q n ↑ , n ↓ are the canonical partition functions for n ↑ + n ↓ particles. a r X i v : . [ h e p - l a t ] O c t o further analyze the thermodynamics of many-body systems we can study the grand-canonicalpotential, − β Ω = ln Z . Again, expanding in powers of z ↑ and z ↓ we achieve an expansion about thedilute limit. This yields, − β Ω = ln Z = Q ∞ (cid:88) n ↑ , n ↓ = b n ↑ , n ↓ z n ↑ ↑ z n ↓ ↓ (3)Where b n ↑ , n ↓ are the “virial coe ffi cients” of order N = n ↑ + n ↓ . These are of particular interest asthey give access to the equation of state (EoS), and therefore to observables for these many-bodysystems (i.e. density, pressure, compressibility, etc.). Moreover, current experiments with ultracoldatoms have gained increasing access to finely controllable parameters such as temperature, coupling,polarization, mass imbalance, and trap shape (harmonic, hard wall, lattices), which allows directcomparison of thermodynamic predictions in the dilute limit, which is the focus of this work.In this contribution, we show our progress towards accessing Q n ↑ , n ↓ for the universal regime calledthe unitary Fermi gas [5], with the aim of characterizing properties of finite systems, i.e. those derivedfrom the Helmholtz free energy − β F = ln Q n ↑ , n ↓ of finite clusters. Our second main goal is to use asimilar approach to compute virial coe ffi cients, b n ↑ , n ↓ , at various orders for the one-dimensional Fermigas with attractive, short-range interactions. To that end, we use a variant of the complex Langevinmethod as a means of stochastic quantization. Below, we will make use of the path integral representation of the grand-canonical partition function,namely Z = (cid:90) D σ det M [ σ ] , (4)where σ is a Hubbard-Stratonovich field representing a two-body interaction and the system is as-sumed to be non-relativistic. The specific dynamics of the system at hand is encoded in the matrix M [ σ ]. In particular for a two-species non-relativistic system, it is not di ffi cult to show [6] thatdet M [ σ ] = det( I + z ↑ U [ σ ]) det( I + z ↓ U [ σ ]) , (5)where U [ σ ] does not depend on the chemical potential. The above equation is the source of theso-called sign problem in polarized non-relativistic matter, as U [ σ ] is real for attractive interactionsbut z ↑ (cid:44) z ↓ leads to di ff erent determinants for each species (generally of di ff erent, fluctuating sign).In the sections below, we implement a method to compute canonical partition functions and virialcoe ffi cients, both based on the idea of particle projection, which is in essence a Fourier transform ona variable φ if the fugacity is replaced by taking z → ze i φ . Such a transformation, which makes z acomplex variable, introduces a phase problem in conventional lattice Monte Carlo approaches. The idea of stochastic quantization is that a random process that obeys the Langevin equation ∂σ ( θ ) ∂θ = − δ S [ σ ] δσ ( θ ) + η, (6)ill yield configurations σ that are distributed according to e − S [ σ ] . Here, θ is Langevin time and η isgaussian random noise (see e.g. [3]). We assume our Hubbard-Stratonovich field σ is real. However,when S [ σ ] is complex, the Langevin equation naturally implies that σ must be complexified into σ = σ R + i σ I . Thus, stochastic quantization requires complexified fields and the (spacetime local)Langevin equation now has both real and imaginary parts, σ R ( n + = σ R ( n ) + (cid:15) K R ( n ) + √ (cid:15)η ( n ) , (7) σ I ( n + = σ I ( n ) + (cid:15) K I ( n ) , (8)where the Langevin time has been discretized as θ = n (cid:15) , and (cid:15) is the time step. Again, η , is the realGaussian noise such that (cid:104) η ( n ) (cid:105) =
0. The drift terms are now specified as, K Ra , x = − Re (cid:34) δ S δσ (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) σ → σ R + i σ I (cid:35) , (9) K Ia , x = − Im (cid:34) δ S δσ (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) σ → σ R + i σ I (cid:35) . (10) As a starting point, we use the auxiliary-field representation of the grand-canonical partition functionon the lattice, which reads Z = (cid:90) D σ det( I + z ↑ U [ σ ]) det( I + z ↓ U [ σ ]) , (11)where σ is the Hubbard-Stratonovich auxiliary field, and U [ σ ] encodes all the properties of the Hamil-tonian of the system. Expanding in terms of the low-fugacity limit we rewrite Z = ∞ (cid:88) n , m z n ↑ z m ↓ Q n , m . (12)Particle projection consists in implementing delta functions that select n ↑ and n ↓ in the above sum,that is δ n ↑ , n and δ n ↓ , m . To that end, we first take Z [ z ↑ , z ↓ ] → Z [ z ↑ e i φ ↑ , z ↓ e i φ ↓ ], such that our fugacityexpansion of Z becomes, Z [ z ↑ e i φ ↑ , z ↓ e i φ ↓ ] = (cid:88) n , m z n ↑ e in φ ↑ z m ↓ e im φ ↓ Q n , m . (13)We may thus project out the desired individual partition functions Q n ↑ , n ↓ using a Fourier integral inthe usual way, namely (cid:90) π d φ ↑ π d φ ↓ π Z [ z ↑ e i φ ↑ , z ↓ e i φ ↓ ] z − n ↑ ↑ e − in ↑ φ ↑ z − n ↓ ↓ e − in ↓ φ ↓ = (14) = (cid:90) π d φ ↑ π d φ ↓ π (cid:88) n ↑ , n ↓ z n ↑ e in φ ↑ z m ↓ e im φ ↓ Q n , m z − n ↑ ↑ e − in ↑ φ ↑ z − n ↓ ↓ e − in ↓ φ ↓ (15) = (cid:88) n ↑ , n ↓ z n ↑ − n ↑ z n ↓ − m ↓ Q n , m δ n , n ↑ δ m , n ↓ = Q n ↑ , n ↓ . (16)n practice, we may define Φ = ( σ, φ ↑ , φ ↓ ) as a new auxiliary field variable and perform the Fourierintegral as part of stochastic process (see below). For this purpose, we include φ ↑ and φ ↓ as part ofour action and integral as follows, Q n , m = (cid:90) D Φ exp (cid:0) − S n , m [ Φ ] (cid:1) , (17)where S n , m [ Φ ] = − ln (cid:104) det( I + z ↑ U [ σ ]) det( I + z ↓ U [ σ ]) z − n ↑ e − in φ ↑ z − m ↓ e − im φ ↓ (cid:105) . (18)From this point on, thermodynamic properties of finite clusters can be obtained in the usual wayby formally taking derivatives of Q n , m with respect to the desired source and computing expectationvalues in the extended ensemble defined by the variable Φ . By definition, the virial expansion of the partition function Z can be written asln Z [ z ] = Q ∞ (cid:88) n = b n z n , (19)where z is the fugacity, b n are the virial coe ffi cients and Q is the single-particle partition function.Based on the above expression, we may use Fourier projection to obtain the virial coe ffi cients: b n = Q (cid:90) π d φ π e i φ n ln Z [ z → e − i φ ] . (20)More generally, we may define the function b n ( ξ ) = Q (cid:90) π d φ π e i φ n ln Z [ z → ξ e − i φ ] = b n ξ n , (21)where the last identity is easy to see from our first definition above.This last equation, however, does not help in computing any of the above: while we have a path-integral form for Z , the presence of the logarithm makes a direct evaluation impractical. Nevertheless,inserting the path-integral form of the partition function and di ff erentiating with respect to ξ we obtaina more useful expression: ∂ b n ( ξ ) ∂ξ = Q (cid:90) π d φ π e i φ n (cid:90) D σ P [ σ, ξ e − i φ ] tr (cid:104) M − ∂ M /∂ξ (cid:105) , (22)where P [ σ, z ] ≡ det M [ σ, z ] / Z [ z ]. Using shorthand notation for the expectation value with respect tothe above weight P , ∂ b n ( ξ ) ∂ξ = Q (cid:90) π d φ π e i φ n (cid:68) tr (cid:104) M − ∂ M /∂ξ (cid:105)(cid:69) φ,ξ = nb n ξ n − , (23)In practice, we do not need a continuous representation for the Fourier integral, particularly becausewe cannot include this as part of the stochastic process as in our particle projection method due tothe logarithm between these integrals (i.e. the integral over the projecting angles does not appear inhe normalization of P ). Rather, we may use an exact discrete representation via the associated anglevalues φ k = π k / N k , where k = , . . . , N k − N k is the number of discretization points, i.e. ∂ b n ( ξ ) ∂ξ = Q N k N k − (cid:88) k = e i φ k n (cid:68) tr (cid:104) M − ∂ M /∂ξ (cid:105)(cid:69) φ k ,ξ = nb n ξ n − , (24)What is interesting about Eq. (23) is that, by using the complex Langevin method to determine thecomplex expectation value as a function of φ k and ξ , followed by summing over φ k with appropriatephases, we may have access to b n itself.In this approach, varying ξ could be very advantageous, as we know the exact ξ dependence of thefinal result via the last identity of Eq. (23). To enhance the sensitivity of the approach, it is preferableto take ξ to be as large as possible, since the b n tend to be small numbers. However, we anticipatethat taking ξ ∼ . ξ makes the approach viable. Our knowledge of thefunctional dependence mentioned above should connect the data points for all ξ with a single curve.In fact, the expected dependence can be accounted for by plotting b n = n ξ n − ∂ b n ( ξ ) ∂ξ (25)as a function of ξ and fitting a constant. If the number of Fourier discretization points N k is largeenough, we expect to see such a constant behavior beyond the statistical precision. We pursue suchan approach below. Figure 1: Running average over Langevin time, τ , where τ is the total number of Langevin time steps, of the totalparticle numbers, (cid:104) ˆ n ↑ (cid:105) and (cid:104) ˆ n ↓ (cid:105) . Separate runs are seento converge to the expected particle number for the 2 + + V = N x , N x = β = .
0, tuned to the unitary limit.The simplest and most important verifi-cation of this method is computing theprojected particle numbers for each fla-vor. The observable in the grand canon-ical ensemble is given by (cid:104) ˆ N (cid:105) = (cid:104) tr( M − [ σ ] ∂ M [ σ ] /∂ z ) (cid:105) , (26)which is to be averaged over Langevintime. As Langevin time progresses, therunning average for the particle num-bers of each flavor, n ↑ and n ↓ , shouldconverge to the selected values. Thisis demonstrated in Figure. (1), wherewe indeed verify that the method con-verges to the correct solutions. Fig-ure. (1) is computed in a relativelysmall spatial volume, V = , but at β = .
0, which implies N τ =
60 lat-tice sites in the time direction. Computations for larger volumes (i.e. N x = , ,
12) are needed tolimit finite-size e ff ects. Likewise, larger values of β are needed to extrapolate to the continuum limit(see below). .2 Virial Projection Results As mentioned above, the goal of this method is to access the virial coe ffi cients b n by carrying outcalculations for varying values of ξ . A first step toward extracting b n is to show that the results overthe range of ξ yield a constant value, at least in some window. Fig. 2, shows results for b n ( ξ ) at β = . .
0, both at √ βg = .
1. We find consistent results for up to order n = ξ = . − .
1, but alsofind statistical noise at smaller ξ , which is particularly evident in the higher order virial coe ffi cients;the latter correspond to progressively higher frequencies in the Fourier projection. Nevertheless, theroughly constant behavior allows us to extract estimates of b n up to n =
5. Systematic e ff ects in thenumber of Fourier points N k and the Langevin time T are currently under study. The data shown herecorresponds to N k =
100 and T = −0.4−0.2 0 0.2 0.4 0.6 0.8 1 0 0.02 0.04 0.06 0.08 0.1 b n ξ b b b b −0.4−0.2 0 0.2 0.4 0.6 0.8 1 0 0.02 0.04 0.06 0.08 0.1 b n ξ b b b b b Figure 2: Virial coe ffi cients b n , for n = , . . . ,
5, as a function of ξ [see Eq. (25)], for N x = β = . . λ = . λ ∈ [0 . , . b n vs. λ . An extrapolation of our answers to λ = λ it may bepossible to obtain approximate virial coe ffi cients at strong coupling. These results for b n vs. λ at β = . β = . b ( λ ) (see Ref. [7]).Notice in Fig. (3) at low coupling our results match extremely well the analytical solution for b ( λ ), but that agreement deteriorates as λ increases. More statistics may be needed at larger couplingto obtain more consistent results. Fig. (3) shows continuous results for most orders of b n , but againwith inconsistencies at higher order as λ increases. Both are encouraging results that point to thepossibility of computing high-order virial coe ffi cients and extrapolating to larger coupling with fairaccuracy. This approach, in that sense, complements the conventional, high-precision way to obtainvirial coe ffi cients, based on solving the n -body problem and computing canonical partition functions. In this contribution we have presented two calculations that use the complex Langevin approach toextract information about finite systems from finite-temperature, grand-canonical calculations. Wefocused on particle-number projection to obtain thermodynamics via canonical partition functions, b n λ = √ β g b b b b b Non−Interactingb Analytic −0.4−0.2 0 0.2 0.4 0.6 0.8 1 0 0.25 0.5 0.75 1 b n λ = √ β g b b b b b Non−Interactingb Analytic
Figure 3: Virial coe ffi cients b n ( λ ) for n = , . . . ,
5, at β = . . N x =
31. Theblack crosses denote the non-interacting values at the continuum limit. Solid line shows the knowncontinuum-limit result for b (see Ref. [7]).and to obtain virial coe ffi cients. We have shown that both methods converge to reasonable resultswithin some constraints on various parameters.More data are needed in order to estimate the free energy for various systems and particle numbers,and to further understand systematic e ff ects. To extrapolate to the continuum limit, we must take1 (cid:28) λ T (cid:28) N x (where λ T = (cid:112) πβ is the thermal wavelength), which amplifies the computationtime for these calculations as both β and N x must be taken to be large. Due to this limitation, weare exploring ways to optimize these calculations. We have shown, however, that particle projectionusing complex Langevin does converge to the correct state, and with the correct number of particlesfor each flavor.Furthermore, our method for virial projection has yielded promising results for one-dimensionalsystems, demonstrating possibilities for computing higher order coe ffi cients (i.e. b , b , and above)as well as extrapolation to stronger values of the coupling. We are currently exploring finite-sizee ff ects of this method by varying the inverse temperature, β , and the volume. In addition, we areinvestigating the e ff ects of varying the total number of Fourier points N k , as decreasing this numberlikewise decreases the total run time of the simulation (linearly). Increased statistics (longer totalLangevin time) is needed for improved accuracy for larger virial coe ffi cients, and particularly forlarger values of the coupling.Both particle projection methods have produced promising results that, with more in-depth study,may provide access to various observables in finite systems that remain unknown and may be directlycompared with ultracold atom experiments. Looking beyond the latter, the ultimate goal is to performfull-fledged, finite-temperature calculations of finite nuclei with realistic interactions, for which thepresent work provides a proof of principle.This material is based upon work supported by the National Science Foundation under Grant No.PHY1452635 (Computational Physics Program). eferences [1] S. Giorgini, L. Pitaevskii, S. Stringari, Rev. Mod. Phys. , 1215 (2008)[2] I. Bloch, J. Dalibard, W. Zwerger, Rev. Mod. Phys. , 885 (2008)[3] G. Aarts, J. Phys. Conf. Ser. , 022004 (2016)[4] C. Gattringer, K. Langfeld, International Journal of Modern Physics A , 1643007 (2016)[5] W. Zwerger, BCS-BEC crossover and the Unitary Fermi Gas (Springer-Verlag, Berlin, 2011)[6] J.E. Drut, A.N. Nicholson, J.Phys.
G40 , 043101 (2013), [7] M.D. Ho ff man, P.D. Javernick, A.C. Loheac, W.J. Porter, E.R. Anderson, J.E. Drut, Phys. Rev. A91