Equilibrium and dynamics of a trapped superfluid Fermi gas with unequal masses

G. Orso,1 L. P. Pitaevskii,2,3 and S. Stringari21Laboratoire Physique Théorique et Modèles Statistiques, Université Paris Sud, Bâtiment 100, 91405 Orsay Cedex, France

2Dipartimento di Fisica, Università di Trento and BEC-INFM, 1-38050 Povo, Italy3Kapitza Institute for Physical Problems, 117334 Moscow, Russia

�Received 11 September 2007; revised manuscript received 19 December 2007; published 13 March 2008�

Interacting Fermi gases with equal populations but unequal masses are investigated at zero temperatureusing local density approximation and the hydrodynamic theory of superfluids in the presence of harmonictrapping. We derive the conditions of energetic stability of the superfluid configuration with respect to phaseseparation and the frequencies of the collective oscillations in terms of the mass ratio and the trappingfrequencies of the two components. We discuss the behavior of the gas after the trapping potential of a singlecomponent is switched off and show that, near a Feshbach resonance, the released component can still remaintrapped due to many-body interaction effects. Explicit predictions are presented for a mixture of 6Li and 40Kwith resonant interaction.

DOI: 10.1103/PhysRevA.77.033611 PACS number�s�: 03.75.Kk, 05.30.Fk, 03.75.Lm, 03.75.Ss

The superfluid behavior of dilute interacting Fermi gasesat very low temperature is now rather well-understood bothfrom the experimental and theoretical point of view �1�. Inparticular the attractive nature of the interaction between thetwo different spin components of the gas is known to play acrucial role along the whole BCS-Bose Einstein condensa-tion �BEC� crossover. This includes the case of small andnegative values of the scattering length, where the ordinaryBCS regime of superfluidity holds, the BEC regime charac-terized by the formation of molecules in the Bose-Einsteincondensed state and the unitary regime where the scatteringlength takes a divergent value and the attraction results inpeculiar many-body effects.

A more recent and intriguing direction is the search forsuperfluidity in mixtures of Fermi gases belonging to differ-ent species, and hence having different masses. Experimen-tally, the most promising candidates are ultracold mixtures of40K and 6Li, where the mass ratio is 6.7, near heteronuclears-wave Feshbach resonances �2�. Equilibrium configurationsof the uniform superfluid phase where the atom densities ofthe two species are equal but the masses are different havebeen theoretically investigated in Ref. �3� within BCS meanfield theory and, more recently, by quantum Monte Carlomethods �4�. Other recent works have explored the interplaybetween different masses and different populations of thetwo components �5–7�.

In this paper we investigate the equilibrium and the dy-namic properties of a trapped Fermi gas with unequal massesusing local density approximation and developing the hydro-dynamic theory of superfluids at zero temperature. We as-sume that the external potentials for the two components,hereafter called ↑ and ↓, are harmonic and given by Vho

� �r�=m���x�

2 x2+�y�2 y2+�z�

2 z2� /2, where �= ↑ ,↓ and m↑ and m↓are the atomic masses. Since the two species have differentmagnetic and optical properties, the trapping frequencies �i�can be tuned separately. We consider a mixture of two dif-ferent fermionic species with equal populations N↑=N↓=N /2, corresponding to the most favorable condition forCooper pairing. We assume that the gas is superfluid at zerotemperature and we explore configurations where the densi-

ties of the two components are equal and move in phase:n↑=n↓�n /2, v↑=v↓�v. We will not consider here exoticpolarized phases, like the FFLO phase, whose relevance fortrapped Fermi gases is still unclear at present.

At equilibrium, where v=0, the atomic density profilen0�r� of the gas is given by the local density �also calledThomas-Fermi� approximation for the chemical potential

�0 = �„n0�r�… + Vho�r� , �1�

where �0 is fixed by the normalization condition �n0�r�dr=N and ��n�=�e /�n is the chemical potential of uniformmatter, e�n� being the energy per unit volume of the homo-geneous phase. In Eq. �1� we have introduced the effectivetrapping potential

Vho�r� =1

2�Vho

↑ + Vho↓ � =

m

2��x�

2 x2 + �y�2 y2 + �z�

2 z2� , �2�

given by the average of the two potentials and

�i2 =

m↑�i↑2 + m↓�i↓

2

m↑ + m↓�3�

with m= �m↑+m↓� /2. The averaging of the effective potentialfollows from the assumption that the density profiles of thetwo components are the same. Actually in the superfluidphase the densities of the two components are equal, even ifthe trapping potentials, in the absence of interactions, wouldgive rise to different equilibrium profiles at zero temperature,i.e., even if the oscillator lengths � / �m����1/2 of the twocomponents do not coincide. If the oscillator lengths areequal, the effective frequencies �3� simply reduce to the geo-metrical averages �i=��i↑�i↓.

Let us discuss the behavior of the equation of state alongthe BCS-BCS crossover and the corresponding shape of thedensity profiles. At unitarity, where the scattering length adiverges, the equation of state takes the universal form �1��= �3�2�2/3�2�1+��n2/3 /4mr, where mr=m↑m↓ / �m↑+m↓� isthe reduced mass and � is a dimensionless parameter de-pending on the mass ratio m↑ /m↓ and accounting for the

PHYSICAL REVIEW A 77, 033611 �2008�

1050-2947/2008/77�3�/033611�4� ©2008 The American Physical Society033611-1

interacting effects in this strongly interacting regime. By in-serting �=�n2/3 in Eq. �1�, with �= �3�2�2/3�1+�� /4mr, wefind that the density distribution of the resonant gas takes the

usual form n�r�= ��0− Vho�r��3/2 /�3/2, the Thomas-Fermi ra-dii Ri, where the density vanishes, being given by

Ri = aho�24N�1/6�1 + ��1/4 �ho

�i

, �4�

where �ho= ��x�y�z�1/3 is the geometrical average of thethree effective oscillator frequencies and aho

2

=� / �ho�m↑m↓�1/2. The value of � is known for the specialcase of equal masses m↑=m↓, where �=−0.58 �8,9�. Prelimi-nary Monte Carlo calculations suggest that � depends veryweakly on the mass ratio �10�. In the deep BCS limit, corre-sponding to a weakly attractive interaction �n�a�3�1�, theequation of state and the Thomas-Fermi radii are given bythe same expressions holding at unitarity by simply setting�=0.

In the limit of small and positive scattering length, corre-sponding to na3�1, the gas instead corresponds to a BEC ofdiatomic heteronuclear molecules of mass m↑+m↓=2m andmolecular density n /2. In this regime the equation of state isfixed by the repulsive interaction between molecules andtakes the usual bosonic form �m=gmn /2, where the couplingconstant gm is related to the molecule-molecule scatteringlength am by gm=2��2am /m. The exact value of the molecu-lar scattering length has been calculated by Petrov et al. �11�as a function of the atom scattering length a and the massratio m↓ /m↑. The bosonic chemical potential is related to thefermionic one by �m=−Eb+2�, where Eb=−�2 /2mra

2 is thetwo-body binding energy. From Eq. �1�, the density profile is

then given by n�r�= ��0− Vho�r��2m /��2am corresponding tothe Thomas Fermi radii

RiBEC = aho

BEC15Nam

2ahoBEC 1/5 �i

�ho

, �5�

where �ho= ��x�y�z�1/3 and ahoBEC=� /�2m�ho. Notice that

ahoBEC differs from the oscillator length aho defined above for

the resonant case.We now take advantage of the fact that the trapping po-

tentials of the two different species can be tuned separatelyto suggest an experiment pointing out in a direct way theattractive role of the interactions in the presence of a Fesh-bach resonance. After generating the equilibrium configura-tion discussed above we switch off the confining potential ofa single species, say Vho

↓ =0, corresponding to a change of thetrapping frequencies �3� into the new values

�inew =� m↑m↑ + m↓

�i↑. �6�

The potential can be switched off either adiabatically, bring-ing the system into a new equilibrium configuration, or sud-denly, giving rise to the excitation of collective oscillations�see discussion in the second part of this work�.

In the absence of interactions, the ↓ atoms would fly awayleaving an ideal gas of ↑ fermions trapped in the harmonicpotential Vho

↑ . This will be also the case in the deep BCS

superfluid regime where interactions are too weak to keepthe ↓ fermions confined. In the other �more robust� superfluidregimes, however, the released atoms do not necessarily es-cape to infinity but can remain trapped due to the attractiveinteraction with the other species. This statement is obviousin the BEC regime, where each ↓ fermion forms a boundmolecule with a corresponding ↑ particle. At unitarity, how-ever, the two-body binding energy Eb=−�2 /2mra

2 vanishesmeaning that no molecules can exist in vacuum. In this casethe trapping of the ↓ component is a pure many-body effectreflecting the attractive nature of the interatomic force.

It is not difficult to derive explicit conditions for the en-ergetic stability of the new superfluid configuration. The sta-bility is ensured if the energy ES of the configuration wherethe two components remain trapped and fully overlapped issmaller than the energy EN of the normal state where the gasis phase separated and only the ↑ atoms are trapped. At uni-tarity the energy of the trapped superfluid state is given by

ES = ��ho�3N�4/3

8��m↑ + m↓�

mr

�1 + � �7�

and, in the absence of trapping for the ↓ atoms, one has�ho=�m↑ / �m↑+m↓���x↑�y↑�z↑�1/3. Since the energy of atrapped gas of noninteracting ↑ fermions is given by EN=���x↑�y↑�z↑�1/3�3N�4/3 /8, we find that the superfluid gas,where both components remain trapped, is energeticallystable if the condition

�1 + ���m↑ + m↓�

m↓� 1 �8�

is satisfied. Taking into account that � barely depends on themass ratio �10� and hence remains close to the equal massvalue �=−0.58, the above condition is always satisfied ifm↑m↓, showing that the superfluid remains energeticallystable if we release the potential of the heavy species. Con-versely, if the condition �8� is violated, the superfluid con-figuration corresponds to a metastable state which is ener-getically unstable toward phase separation. It is alsointeresting to compare the Thomas Fermi radii RiS and RiN ofthe trapped cloud in the �new� superfluid and in the separatednormal phase, respectively. A simple calculation yieldsRiS /RiN= �1+��1/4��m↑+m↓� /m↓�1/4, showing that the super-fluid phase corresponds to the configuration with smaller ra-dii if and only if Eq. �8� is satisfied.

Let us now discuss the macroscopic dynamic behavior ofthe superfluid. This is obtained by deriving the hydrody-namic equations of motion in terms of the atom densityn�r , t� and the velocity field v�r , t�. The Lagrangian L of thesystem, in the local density approximation �LDA�, is givenby

L =� dr�e�n� + �Vho↑ + Vho

↓ �n

2+

n

2

�

�t+

1

2�m↑ + m↓�v2n

2 ,

�9�

where is the phase of the order parameter��↑�r , t��↓�r , t������↑�r��↓�r���ei�r,t�, ���r , t� being thefermionic field operators. The superfluid velocity v and the

ORSO, PITAEVSKII, AND STRINGARI PHYSICAL REVIEW A 77, 033611 �2008�

033611-2

phase entering the Lagrangian �9� are related by the mostimportant condition

v =�

m↑ + m↓� , �10�

where m↑+m↓ is the mass of the pair. Equation �10� can bederived microscopically by noticing the state moving withvelocity v is obtained from the steady state by applying thegauge transformation ���r�→���r�eim�v·r/� to the twoFermi field operators �12�.

Taking Eq. �10� into account, the equations of motion ofLagrangian �9� yield the hydrodynamic equations

�n

�t+ ��nv� = 0, �11�

m�v

�t+ ���n� + Vho�r� +

1

2mv2 = 0. �12�

Equations �11� and �12� apply to weakly as well as tostrongly interacting superfluids and permit one to calculatethe macroscopic dynamics of the system �expansion and col-lective oscillations� once the equation of state is known.

In uniform matter �Vho=0� the linearized solutions arephonons with sound velocity fixed by c2= �n /m��� /�n.At unitarity, the sound velocity is given by c

=�kF�1 / �3m↑m↓��1+�, where kF= �3�2n�1/3. In the BCSlimit, the same formula holds with �=0 �13�.

For trapped configurations, the frequencies of the collec-tive oscillations are obtained �14� by linearizing the hydro-dynamic equations around the equilibrium distribution n0�r�and depend on the applied external potentials only throughthe effective oscillator frequencies �3�. They also depend onthe equation of state and in particular on its density depen-dence �15�. A special case is the center of mass oscillationsof the cloud whose frequency, in the superfluid phase, isindependent of the equation of state and given, for each di-rection, by �i. We emphasize that in general �i differ from

either �i↑ and �i↓, pointing out that in the superfluid phasethe two components always oscillate in phase as a conse-quence of the pairing mechanism.

Equations �11� and �12� remain valid even if the effectivetrapping frequencies depend explicitly on time, i.e., if �i= �i�t�. Remarkably, if the chemical potential has the powerlaw dependence ��n on the density and the confining po-tential is harmonic, the hydrodynamic equations admit aclass of exact solutions given by v�t ,r�=� ·r and n�t ,r�=n0�x /bx ,y /by ,z /bz� / �bxbybz�, where �i�t� ,�i�t� are time-dependent parameters. From the scaling form of the density,we see that the parameters �i�t� are related to the ThomasFermi radii Ri�t� of the evolving cloud according to Ri�t�=Ri�0�bi�t�. Inserting the scaling ansatz in Eq. �11� yields the

relationships �i= bi /bi and one then obtains a set of coupledordinary differential equations �16�.

bi + �i2�t�bi =

�i2�0�bi

1

�bxbybz� . �13�

Equation �13� applies both to resonance, where =2 /3, andto the BEC regime where =1.

If one suddenly switches off the trapping potential for the↓ component, corresponding to suddenly setting the effectivefrequencies to the new values �6�, the superfluid gas will startoscillating around the new equilibrium configuration. Thelatter is described by the Thomas Fermi radii Rinew=Ri�0� /�1/4, with �=m↓ / �m↑+m↓�, as follows from the staticsolution of Eq. �13�. The dynamics of the oscillating cloudcan be calculated from Eq. �13� by substituting �i�t�0�= �inew and employing the initial conditions bi�0�=1 and

bi�0�=0.For simplicity, let us assume that the two trapping poten-

tials are axis-symmetric ��x�=�y������. We further as-sume that the two components have initially equal oscillatorlengths �and hence, in the degenerate limit, equal densityprofiles also in the absence of interaction�, corresponding to��

2 �0�=��↑��↓ and �z2�0�=�z↑�z↓.

We have solved Eq. �13� for a mixture of 40K and 6Lifermions in an elongated trap with aspect ratio �= �z�0� / ���0�=0.2. In Fig. 1 we plot the calculated timeevolution of the radius of the cloud in the radial direction,after we suddenly switch off the trapping potential of 40Kfermions, corresponding to �=0.869. Since � is close to 1,the initial configuration Rinew=1.036Ri�0� is close to equilib-rium and the resulting oscillation is a linear superpositionof the two breathing modes, in the radial and longitudinaldirections. For elongated cloud ���1�, the correspondingfrequencies are given by �rad=�10 /3�1/2���0� and �axial

=�12 /5�1/2��z�0� �1�.In Fig. 2 we instead plot the time dependence of the trans-

verse radii of the cloud, after we suddenly switch off thetrapping potential of the 6Li fermions, corresponding to �=0.131 and Rinew=1.66Ri�0�. In this case, the initial configu-ration is far from equilibrium and nonlinear effects play animportant role. In particular we see that the breathing of thecloud around the equilibrium configuration is no longer sym-metric. We emphasize that in this second case the superfluid

0 20 40 60τ

1

1.04

1.08R

⊥(τ

)/R

⊥(0

)

FIG. 1. Time evolution of the size of the cloud along the radialdirection after switching off suddenly the trapping potential for 40Kfermions. Here �= ���0�t. The dashed line corresponds to the newequilibrium value �see text�.

EQUILIBRIUM AND DYNAMICS OF A TRAPPED … PHYSICAL REVIEW A 77, 033611 �2008�

033611-3

configuration is energetically unstable. A major question inthis case is to understand the decay mechanisms and the roleplayed by the sudden excitation of the collective modes. Infact the energy E calculated after the sudden release of V↓can be significantly larger than the energy EN of the phaseseparated configuration, the corresponding ratio being given,in the general case, by

E

EN= �1 + ��m↑ + m↓

m↓

2 + �

2�1 + �, �14�

where �=m↓�↓2 /m↑�↑

2. As a consequence, differently fromthe case of Fig. 1 where m↓�m↑ and hence E /EN�1, in thecase of Fig. 2 the gas can oscillate in the superfluid phaseonly in the presence of a sufficiently high energy barrierpreventing it from phase separation. This experiment wouldconsequently open new perspectives in the study of the de-cay mechanism of metastable configurations.

In conclusion, we have derived the equilibrium conditionsand hydrodynamic equations of a superfluid Fermi gas withunequal masses and investigated the behavior of a gas atunitarity after the release of the trapping potential of a singlecomponent. We have shown that, under appropriate condi-tions, the superfluid phase is stable against phase separationof the two components. As a result, the released fermionsremain confined in the trap due to the pairing with the othercomponent, pointing out in a direct and remarkable way theattractive nature of the interatomic forces near a Feshbachresonance.

This work was supported by the Marie Curie programunder Contract No. EDUG-038970 and by the Ministerodell’Istruzione, dell’Università e della Ricerca �MIUR�.

�1� S. Giorgini, L. Pitaevskii, and S. Stringari, e-print arXiv:cond-mat/0706.3360, Rev. Mod. Phys. �to be published�.

�2� E. Wille et al., Phys. Rev. Lett. 100, 053201 �2008�.�3� H. Caldas, C. W. Morais, and A. L. Mota, Phys. Rev. D 72,

045008 �2005�.�4� J. Von Stecher, C. H. Greene, and D. Blume, Phys. Rev. A 76,

053613 �2007�.�5� W. V. Liu and F. Wilczek, Phys. Rev. Lett. 90, 047002 �2003�.�6� S. T. Wu, C. H. Pao, and S. K. Yip, Phys. Rev. B 74, 224504

�2006�; C. H. Pao, shin-Tza Wu, and S. K. Yip, Phys. Rev. A76, 053621 �2007�.

�7� M. Iskin and C. A. R. Sa de Melo, Phys. Rev. Lett. 97, 100404�2006�; Phys. Rev. A 76, 013601 �2007�.

�8� J. Carlson, S.-Y. Chang, V. R. Pandharipande, and K. E.Schmidt, Phys. Rev. Lett. 91, 050401 �2003�.

�9� G. E. Astrakharchik, J. Boronat, J. Casulleras, and S. Giorgini,Phys. Rev. Lett. 93, 200404 �2004�.

�10� G. E. Astrakharchick, D. Blume, and S. Giorgini �private com-munication�.

�11� D. S. Petrov, C. Salomon, and G. V. Shlyapnikov, J. Phys B38, S645 �2005�.

�12� See, for instance, L. P. Pitaevskii, S. Stringari, and G. Orso,Phys. Rev. A 71, 053602 �2005�.

�13� L. He, M. Jin, and P. Zhuang, Phys. Rev. B 73, 220504�R��2006�.

�14� S. Stringari, Phys. Rev. Lett. 77, 2360 �1996�.�15� S. Stringari, Europhys. Lett. 65, 749 �2004�; G. E. Astra-

kharchik, R. Combescot, X. Leyronas, and S. Stringari, Phys.Rev. Lett. 95, 030404 �2005�.

�16� C. Menotti, P. Pedri, and S. Stringari, Phys. Rev. Lett. 89,250402 �2002�. These equations have been first derived forBEC � =1� by Y. Castin and R. Dum, ibid. 77, 5315 �1996�and by Y. Kagan, E. L. Surkov, and G. V. Shlyapnikov, Phys.Rev. A 54, R1753 �1996�.

0 40 80 120τ

1

2

3R

⊥(τ

)/R

⊥(0

)

FIG. 2. Time evolution of the size of the cloud along the radialdirection after switching off suddenly the trapping potential for 6Lifermions. Here �= ���0�t. The dashed line corresponds to the newequilibrium value �see text�.

ORSO, PITAEVSKII, AND STRINGARI PHYSICAL REVIEW A 77, 033611 �2008�

033611-4