{'fraction_non_alphanumeric': 0.0366100172194355, 'fraction_numerical': 0.010930598188215917, 'mean_word_length': 4.27567140600316, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Shell structure in the density profiles for noninteracting fermions in anisotropic harmonic confinement

We develop a Green's function method to evaluate the exact equilibrium particle-density profiles of noninteracting Fermi gases in external harmonic confinement in any spatial dimension and for arbitrary trap anisotropy. While in a spherically symmetric configuration the shell effects are negligible in the case of large number of particles, we find that for very anisotropic traps the quantum effects due to single-level occupancy and the deviations from the Thomas-Fermi approximation are visible also for mesoscopic clouds.

We develop a Green's function method to evaluate the exact equilibrium particle-density profiles of noninteracting Fermi gases in external harmonic confinement in any spatial dimension and for arbitrary trap anisotropy. While in a spherically symmetric configuration the shell effects are negligible in the case of large number of particles, we find that for very anisotropic traps the quantum effects due to singlelevel occupancy and the deviations from the Thomas-Fermi approximation are visible also for mesoscopic clouds. Ultra-cold Fermi gases of alkali atoms are novel quantum systems which are becoming available by the techniques of atom trapping and cooling. The quantum degeneracy regime has been reached with various atomic species [1] and a great effort is now devoted to study the collective properties [2] and to the possibility of observing spatial phase separation [3] and the transition to a superfluid phase [4]. Atomic Fermi gases are produced in magnetic or optical traps, providing an inhomogeneous external confinement which can be well approximated as harmonic. If a single Zeeman sublevel is populated, interparticle interactions can be neglected at the very low temperatures of the experiments since the Pauli principle forbids collisions in the s-wave channel and dipole-dipole interactions are extremely weak for alkali-metal atoms.
The high purity of the samples and the high resolution of the detection techniques make these systems ideal candidates for the study of single-level quantum properties on a mesoscopic scale, such as the shell structure in the equilibrium density profiles. In the experiments the strength of the external confinement can be tuned to vary the anisotropy of the harmonic trap in order to reach quasi-onedimensional and quasi-twodimensional configurations where, as we will show below, the quantum effects on the equilibrium profiles are greatly enhanced.
The exact equilibrium profiles for the particle and kinetic-energy density of non-interacting fermions under isotropic harmonic confinement at zero temperature have been obtained by Brack and van Zyl [5] in D dimensions (D=1,2,3) in terms of sums of Laguerre polynomials. Schneider and Wallis [6] have computed the three-dimensional density profiles of a cigar-shaped Fermi cloud using Hermite polynomials.
An alternative method which does not rely on the evaluation of high-order polynomials has been proposed for the calculation of 1D density profiles [7]. This Green's function method replaces the use of single-particle orbitals in favour of matrix elements of simple operators. The basic idea of the method is to exploit a formal analogy between the expression of the density of singleparticle states in the energy domain and the 1D particledensity profiles in order to rewrite the latter in terms of a suitably defined Green's function operator in coordinate space. As originally formulated the method is intrinsically onedimensional, since it maps the energy variable into a 1D spatial coordinate.
In this Letter we extend the Green's function method to any spatial dimension by making use of a systematic reduction of dimensionality in the evaluation of the density profiles. This method allows us to treat fully anisotropic Fermi clouds with a large number of particles and to investigate how the shell effects are enhanced with respect to the isotropic case. We find that very anisotropic clouds show the shell structure best in the tight direction of the confinement and that in this case the Thomas-Fermi approximation fails to describe the density profiles even at large numbers of particles.
Green's function method in dimension D > 1. The particle density profile for noninteracting fermions under external confinement at zero temperature is given by the zeroth moment of the Dirac density matrix [9] in terms of single-particle orbitals ψ {i} (r): (1) Here {i} is a complete set of quantum numbers which univocally identify the single-particle energy levels E {i} and the sum runs up to the highest occupied level corresponding to the Fermi energy E F . Let us first consider for simplicity the 2D case. For a noninteracting system the Hamiltonian is separable,Ĥ = H x +Ĥ y , and hence the eigenvalues can be found independently in the two cartesian directions, E {i} = E ix +E iy , taking {i} = (i x , i y ) as a set of quantum numbers. Since the wavefunction factorizes we can write the density profile as where the upper indices I x ,Ĩ y of summation are fixed by the value of the Fermi energy from the implicit relations E F = E (Ix,1) and E F = E (ix,Ĩy) , which we assume to be able to solve; by definition,Ĩ y depends on the value of the quantum number i x and its highest value I y is obtained by setting i x = 1. In Eq. (2) we have chosen to index the levels (i x , i y ) starting from (1, 1). We are thus led to the final result that the 2D density profile can be rewritten in the form of 1D density profile for I x +1 particles in the x direction, where the i x -th term is weighted by the 1D density profile ofĨ y + 1 particles in the y direction, This expression allows us to apply recursively the Green's function method in coordinate space: we first obtain the 1D density profiles [7] in the y direction, whereĜ(y) = (y + iε −ŷ) −1 is the Green's function and y is the position operator in the y direction expressed in the basis of the single-particle energy eigenstates, and the trace runs over the firstĨ y + 1 elements of the matrix. We can then obtain the full 2D density profile by taking the trace HereĜ(x) = (x+ iε −x) −1 is the Green's function operator in the x direction, and N 1D Iy (y) is a diagonal "weight" matrix whose i x -th element is given by n 1D Iy (y). Again, the trace is taken only over the first I x + 1 levels.
We have applied this general method to the case of 2D harmonic confinement. The Hamiltonian is given bŷ wheren x,y =â † x,yâ x,y are the number operators witĥ a x,y | ψ ix,y = i x,y − 1 | ψ ix,y−1 andâ † x,y | ψ ix,y = i x,y | ψ ix,y+1 , ψ ix,y being respectively the wave functions of the HamiltoniansĤ x,y of the i x,y -th energy levels, and we have considered the case of a generic anisotropy k = ω y /ω x . For this system the energy levels E ix,iy are known and the upper limits for the sums appearing in Eq. (2) can be evaluated explicitly: I x is defined by the relation and analogouslyĨ y reads Eq. (5) allows us to evaluate efficiently the 2D density profile by exploiting the relation between the trace of a matrix and the determinant of its inverse [7]: where I Ix+1 is a matrix with the first I x + 1 diagonal elements equal to 1 and null elsewhere. The calculation of the 2D density profile thus reduces to evaluating the determinant of a tridiagonal matrix, which can be expressed with a recursion relation as the product of infinite terms: where the factors can be written as a 1 = λn 1D Iy , for 1 < j ≤ I x + 1, a j = λn 1D Iy (y) + j/(2(x − a j−1 )) and for j > I x + 1, a j = j/(2(x − a j−1 )). In the latter expression we have scaled the coordinate x in units of the harmonic oscillator length l x in the x direction, n 1D (y) in units of l −1 y and the resulting 2D profile in units of (l x l y ) −1 , with l x,y = ( /mω x,y ) 1/2 .
The 3D system. The procedure outlined above can be applied recursively to describe the physical case of three spatial dimensions. In the case of noninteracting fermions the Hamiltonian is separable in the three directions and the density profile can be written as n 3D EF (x, y, z) = n 3D Ix,Iy,Iz (x, y, z) = where the highest indices of the sum are fixed by the implicit relations E F = E (Iz,1,1) , E F = E (iz,Ĩx,1) and E F = E (iz,ix,Ĩy) and I x and I y are the highest value ofĨ x andĨ y . Analogously to the case of the 2D confinement, the 3D particle density (11) can be reduced to an expression for an effective 1D density profile for I z + 1 particles, with the i z -th term weighted by the 2D density profile n 2D Ix,Ĩy (x, y) Here,Ĝ(z) is the Green's function for theẑ operator, N 2D Ix,Ĩy (x, y) is a diagonal matrix whose i z -th element is given by n 2D Ix,Ĩy (x, y) andĨ y =Ĩ y (i x = 1; i z ) is the highest value ofĨ y at fixed i z .
In the case of harmonic confinement, the indices of the occupied levels are functions of the anisotropy parameters k = ω y /ω x and l = ω x /ω z of the trap, and the density profile, in units of (l x l y l z ) −1 takes the simple form where the effect of the presence of the confinement in the {x, y} plane is taken into account by the a j factors, with a 1 = λn 2D Ix,Iy (x, y), for 1 < j ≤ I z +1, a j = λn 2D Ix,Ĩy (x, y)+ j/(2(z − a j−1 )) and for j > I z + 1, a j = j/(2(z − a j−1 )). In Eq. (13) the coordinates are in units of the harmonic oscillator lengths l x , l y and l z .
Numerical results. We show here some results for the density profiles in the case of a 2D harmonic trap for various values of the anisotropy parameter. Our 2D results correspond to column densities of 3D systems in which only one level is occupied in the integrated z direction. Our numerical procedure is the following. The first step is to calculate and store the 1D density profiles from 1 tõ I y + 1 particles. Then, the scheme given in Eqs. (9)-(10) is easily implemented, performing the calculation of the determinant (10) up to the product of M terms. This approximation, which corresponds to neglect all the excited states larger than M , is the same as the one used for the evaluation of the diagonal elements of the matrix N 1D Iy (y). Our calculations have been performed using M = 10 6 and ε = 0.01. In Figs. 1-3 we report the density profiles of three different closed-shell systems with a thousand of particles. Figs. 1 and 2 refer to the case a trap with a large anisotropy, having only two and three levels occupied in the y direction, respectively. The prominent shell structure which comes out in the direction of tight confinement is completely lost in the case of an isotropic trap for such a large number of fermions (Fig. 3). Another peculiarity of the system with a large anisotropy is that the tails of the profiles in the direction of weak confinement have a 1D behavior: for large values of x only the Hermite polynomials of high degree have a significant weight, but these terms (the last I x + 1 − kI y ones) are weighted by the same factor -the density of a single fermion in the y direction -giving rise to an essentially 1D profile in the tails as is shown in the inset of Fig. 4(a). However, in the case of a mesoscopic system where a large number of levels are occupied in the longitudinal direction, the 1D shell structure is not well visible.
We have also compared the exact profiles with those given in the Thomas-Fermi (TF) approximation for N fermions at anisotropy parameter k with E F = (2kN ) 1/2 ω x evaluated in the TF limit. While for the symmetric case (k = 1) and N ≃ 1000 the TF profile is practically indistinguishable from the exact one, the local density approach completely fails in describing not only the narrow profile in the y direction (see Fig. 4(b)) but also the profile in the x direction ( Fig. 4(a)) even for a large number of occupied levels.
In conclusion in this Letter we have given a general formula for the exact particle density of 2D and 3D fermions in external confinement in terms of a Green's function in coordinate space. Our approach allows us to treat systems with arbitrary anisotropy since we can deal separately with each cartesian direction by reducing recursively the problem to the evaluation of onedimensional profiles, where the other degrees of freedom are taken into account through renormalization factors. The same idea can be exploited to evaluate the higher-order moments of the density matrix. We have used this method, which is particularly suited to the case of harmonic confinement, to evaluate 2D density profiles of mesoscopic systems for various values of the anisotropy parameter. In the case of large anisotropy we obtain a prominent shell structure, which should be experimentally observable in mesoscopic clouds at temperatures lower than the harmonic-oscillator energy spacing in the tight direction [10]. For such anisotropic systems we have found that the local density approximation completely fails in reproducing the exact profiles even for a large number of atoms (N ≃ 1000). As a consequence, we expect that the TF density functional, which has been shown to work well for isotropic traps [5] should be inadequate to deduce the kinetic energy density in the case of large anisotropy.

ACKNOWLEDGMENTS
We are indebted to Prof. M.P. Tosi for many useful discussions and encouragement. We acknowledge support from INFM through the PRA2001 program.