Ground state of the dipolar Bose-Einstein condensate

Ground state of the dipolar Bose-Einstein condensate

T. F. Jiang*

Institute of Physics, National Chiao Tung University, Hsinchu, 300, Taiwan W. C. Su

Department of Civil Engineering, National Chiao Tung University, Hsinchu, 300, Taiwan

共Received 10 June 2006; published 4 December 2006兲

We present a variational method to study the ground state of the newly realized52Cr dipolar Bose-Einstein condensate. Besides the usual contact potential term in the mean-field equation, there is an additional long-range and anisotropic dipole-dipole interaction potential. We develop an efficient Newton-Raphson’s scheme to solve the condensate state. The solution shows a double-peak feature in the loosely confined dimension. Compared to the existing single-peak calculations, the double-peak solution has lower energy and reflects the distinct property of the dipole-dipole interaction. Our method is easy and efficient to use for future investiga-tions of the dipolar systems.

DOI:10.1103/PhysRevA.74.063602 PACS number共s兲: 03.75.Hh, 05.30.Jp

I. INTRODUCTION

Since the realization of Bose-Einstein condensation 共BEC兲 of alkali metal atoms in 1995 关1兴, its explosive

progress has made an important impact on science. The con-tact atomic interaction has played a significant role in sys-tems studied so far. Syssys-tems with long-ranged interaction were proposed a few years later关2兴. The properties of

dipole-dipole forces in an ultracold atomic systems were studied. Yi and You studied the properties of an electric field induced dipolar system关3兴. The stability of a dipolar BEC system and

the possible species for realization were discussed by Santos et al. and by Fischer关4兴. Góral and Santos also studied the

ground state and excitations of the dipolar BEC 关5兴. The

collapse of dipolar BEC was discussed by Lushnikov 关6兴.

Góral et al. studied the dipolar BEC in optical lattice 关7兴.

Dipolar spinor BEC was also investigated recently关8兴.

Cal-culations through the quantum hydrodynamics model and by Thomas-Fermi approximations were also carried out关9兴.

Experimentally, Bose-Einstein condensation for systems with long-ranged interatomic interactions was first realized in 2005 by Stuttgart group with aligned chromium atoms 关10兴. The 52Cr atom used has permanent magnetic dipole moment M of 6 bohr magneton. In their recent experiment 关11兴, there are 105_{condensate atoms in an anisotropic trap of} frequencies ␻x:␻y:␻z= 942: 712.5: 116.5共Hz兲. The atomic

magnetic dipoles are aligned along the y direction by a dc magnetic field, and the s wave scattering wavelength as was

measured to be 105 bohr radius. Under the mean-field theory, the Gross-Pitaevskii equation 共GPE兲 that describes the BEC is written as

冉

− ប 2 2mⵜ 2_{+ V} ext共rជ兲 + Ngcn共rជ兲 + NVdd共rជ兲

冊

⌿共rជ兲 =␮⌿共rជ兲. 共1兲 In the GPE, N is the number of condensate atoms, and m is the atomic mass. The trap potential is

Vext= m

2共␻x2x2+␻y2y2+␻z2z2兲, gc=

4␲ប2_a

s

m . The order parameter

is normalized as兰兩⌿兩2_d3_rជ_{= 1, and n共r}ជ_{兲=兩⌿共r}ជ_兲兩2_{. The} dipole-dipole potential Vdd with dipoles aligned along the y axis is

written as follows: Vdd共rជ兲 = gd

冕

冋

1 兩rជ− rជ

⬘

兩3− 3共y − y

⬘

兲2 兩rជ− rជ

⬘

兩5

册

兩⌿共rជ

⬘

兲兩 2_d3_rជ

_⬘

_{, 共2兲} where gd= ␮0M2

4␲ . We can see that the interaction potential among atoms becomes nonlocal instead of just the pseudo-potential of zero-ranged interaction.

In addition to the nonlocal property in the dipolar poten-tial, the potential is also anisotropic. Along the polarization direction, the forces between dipoles are attractive and are repulsive in the directions perpendicular to the polarization axis. The effect was discussed by Pfau et al. 关2兴 and the

anisotropic nature was shown for a spherical trap case. The dipole-dipole interaction leads to the result of double-peak in the order parameter along the direction of the repulsive force. To our knowledge, no other theoretical calculations treated this property further because it is not easy to implement double-peak trial functions. In practical experiment关11兴, the

trap is not spherical but is anisotropic. It is tightly confined in the x- and y-directions and loosely confined in the z-direction. Since the dipole-dipole force is much weaker than the confined forces in the x- and y-directions, the dipole-dipole effect is negligible. However, along the z-direction confinement is much weaker and the additional repulsive dipole-dipole force may exhibit double-peak fea-ture in the order parameter under experimental conditions.

We present in this paper a variational study for this prop-erty. We found that the double-peak in order parameter will lower down the total energy of the system in comparison with the single Gaussian wave function which has been popularly used in most previous related calculations. Due to the nonlocality, in addition to the nonlinearity in the GPE, the solution of the equation is not straightforward. We present a robust but simple way to obtain the variational solutions. We compare in detail the total energy, chemical potential, order parameter of single peak and double-peak *Email address: [email protected]

trial solutions in the z-direction with experimental param-eters. We show that the double-peak order parameter corre-sponds to the real experiment and it is quite different from the single-peak picture. The method is easy to apply for fu-ture dipolar BEC systems.

The layout of the paper is as follows. In Sec. II, we de-scribe the energy variational formulation and introduce our method of calculations. In Sec. III we present our results. The discussion and conclusions are followed in Sec. IV.

II. FORMULATION AND METHOD OF CALCULATION

The GPE of Eq.共1兲 can be obtained by minimization of

the following energy functional per particle with constraint 兰兩⌿兩2_d3_rជ_{= 1:} E关⌿兴 =

冕

冉

ប 2 2m兩ⵜ⌿兩 2_{+ V} ext共rជ兲兩⌿兩2+ Ngc 2 兩⌿兩 4 +N 2Vdd共rជ兲兩⌿兩 2

冊

_d3_{rជ. 共3兲}

In the following, we use ប␻y as the energy unit, where

␻y= 2␲⫻712.5 Hz, and L=

冑

ប/m␻y= 0.52261␮m as the

length unit. To include the double-peak property in the z-direction, we employ the following trial functions:

⌿共x,y,z兲 =␾1共x兲␾2共y兲␾3共z兲, ␾1共x兲 = 共␲␣12兲−1/4exp

冉

− x2 2␣₁2

冊

, ␾2共y兲 = 共␲␣2 2_兲−1/4_exp

冉

₋ y 2 2␣₂2

冊

, ␾3共z兲 = c exp

冉

−共z −␣4兲 2 2␣₃2

冊

+ c exp

冉

− 共z +␣4兲2 2␣₃2

冊

, 共4兲 where␣iwith i = 1 , 2 , 3 , 4 are the four variational parameters

and c is used to normalize␾3共z兲. Here two Gaussian func-tions with peaks displaced at ␣4 with the origin in the z-direction are assumed. To compare with the commonly used variational results, we also perform calculations with

␣4= 0 and then there will be three variational parameters. It corresponds to the usual single-peak Gaussian trial solution. We denote the former calculation by four-parameter type and the latter by three-parameter type. In the variational cal-culations, we need to calculate the total energy functional E关␣i兴 and its derivatives. The minimization conditions are

⳵E

⳵␣i

= 0, i = 1,2,3,4. 共5兲 To find out the optimum parameters, we need to calculate the energy functional. In the calculation of the kinetic energy term, we transform the trial function in Eq.共4兲 into

momen-tum space and perform the integration by numerical quadra-ture. The calculations of trap potential energy and contact

energy terms are numerically straightforward, we use the Gauss-Legendre quadratures关12兴 for integrations. The

calcu-lation of the dipolar energy term in coordinate space is a six-dimensional integral and hence is not easy. To make it feasible, the dipolar energy term is transformed into the mo-mentum space, and the integral becomes three-dimensional. With transformations

␳共kជ兲 =

冕

兩⌿共rជ兲兩2_eikជ·rជ_d3_rជ,

Vdd共kជ兲 =

冕

Vdd共rជ兲eikជ·rជd3rជ, 共6兲

then the dipolar energy term is equal to

Edd= N 2 1 8␲3

冕

兩␳共kជ兲兩 2_V dd共kជ兲d3kជ. 共7兲

We can calculate the momentum space density function␳共kជ兲 directly by quadratures. The momentum representation of the dipole-dipole potential has been studied关2,13兴. With the

di-poles aligned along the y-axis and under the mean-field theory, Vdd共kជ兲 = 4␲gd

冋

1 − 3ky2 k2

册

冋

cos共ka兲 共ka兲2 − sin共ka兲 共ka兲3

册

. 共8兲 With these efforts, the dipolar energy term can be calculated efficiently by Gauss-Legendre quadratures too.

The next problem is to find out the variational parameters for a given number of trapped atoms N such that the total energy is minimum. The Newton-Raphson’s scheme is useful for the purpose. Let␣_i共n兲be the value of parameter␣iin the

nth iteration, the next iteration value of␣iwill be␣i

共n+1兲_,

␣i共n+1兲=␣i共n兲−关J−1兴ij

冋

⳵E

⳵␣j

册

n

; i, j = 1,2,3,4. 共9兲 In the expression, repeated index means summation from 1 to 4.关J兴−1is the inverse of the Jacobian matrix关J兴, where in our problem, the matrix elements of关J兴 are

关J兴ij=

⳵2_E关⌿兴

⳵␣i⳵␣j

. 共10兲

Even with the simple Gaussian trial functions with single-peak in the z-direction, the analytic expressions for the prob-lem and the Newton-Raphson’s scheme is still very hard. An example of an analytic expression of total energy in a cylin-drical trap with a two-parameter Gaussian trial function can be found in Ref.关14兴. With the double-peak Gaussian in the

z-direction, the analytic expression for Newton-Raphson’s method becomes even much harder. Our method to the prob-lem is simply using a numerical central difference for the derivatives. For example,

⳵E关⌿兴

⳵␣i

=E共␣i+␦␣i兲 − E共␣i−␦␣i兲 2␦␣i

. 共11兲

We have calibrated our quadrature results of the total energy with the three-parameter trial functions so that the analytic

results can easily be calculated. The accuracy of our numeri-cal quadratures is to the fifth decimal place with 200 grids in the x- and y-axes, 400 grids in the z-axis. So the numerical derivatives work well. In Fig.1, we show the examples of iterations for a case of four-parameter and another for three-parameter. Within a few iterations, the convergence is achieved.

III. RESULTS

We present in Table I the results of double-peak, four-parameter calculations, and results of single-peak, three-parameter calculations. In the calculations, we vary the num-ber of atoms from 1000 to 105_{with experimental parameters} of trap frequencies, gcand gd. As expected, the tabulated␣1 and ␣2 show not much differences in widths of ␾1共x兲 and

␾2共y兲 for both three-parameter and four-parameter calcula-tions because the dipole-dipole effect is relatively smaller than the confined potential and the contact potential. But the widths in␾3共z兲 do show notable differences; and through the wide range of particle numbers, the total energies of double-peak calculations are lower than the single-double-peak calculations. However, the difference is not drastic. This is because the energy of trap potential and contact energy terms are an or-der of magnitude larger than the dipole-dipole term. Here we just perform simulations according to the real experimental parameters. We can expect that if we lower down the trap frequencies and tune the scattering length to a smaller value by the Feshbach resonance method关15兴, the more dominant

TABLE I. Results of variational calculations.␣₁,␣₂are widths of Gaussian density in the x- and in the y-directions, respectively.␣₃in three-parameter calculation is the width of density in the z-direction, while in four-parameter calculations, it is the width of the single Gaussian in a sum of double Gaussians separated at a distance of 2␣₄. The numbers of atoms are from 1000 to 100 000.

N / 1000 Four-parameter Three-parameter Etotal ␣1 ␣2 ␣3 ␣4 Etotal ␣1 ␣2 ␣3 1 2.044 1.031 1.306 3.826 3.995 2.062 1.035 1.313 6.253 5 3.297 1.236 1.679 5.544 5.881 3.338 1.243 1.691 9.180 10 4.202 1.369 1.905 6.452 6.858 4.260 1.378 1.919 10.699 15 4.874 1.462 2.057 7.035 7.483 4.943 1.471 2.073 11.672 20 5.428 1.535 2.174 7.475 7.954 5.505 1.545 2.191 12.405 25 5.905 1.595 2.270 7.832 8.335 5.991 1.606 2.287 12.999 30 6.330 1.648 2.352 8.135 8.659 6.423 1.659 2.370 13.503 35 6.716 1.694 2.424 8.399 8.941 6.814 1.706 2.443 13.942 40 7.070 1.736 2.489 8.633 9.191 7.174 1.748 2.507 14.332 45 7.400 1.774 2.547 8.845 9.417 7.509 1.786 2.566 14.684 50 7.708 1.808 2.600 9.038 9.623 7.822 1.821 2.620 15.006 55 7.999 1.841 2.650 9.216 9.814 8.118 1.854 2.670 15.302 60 8.274 1.871 2.695 9.382 9.991 8.397 1.884 2.716 15.577 65 8.537 1.899 2.738 9.537 10.156 8.664 1.913 2.759 15.835 70 8.787 1.925 2.779 9.683 10.311 8.919 1.939 2.800 16.076 75 9.028 1.950 2.817 9.820 10.457 9.163 1.965 2.838 16.304 80 9.259 1.975 2.853 9.949 10.595 9.397 1.989 2.875 16.520 85 9.482 1.997 2.887 10.075 10.729 9.623 2.012 2.910 16.725 90 9.697 2.019 2.920 10.192 10.854 9.842 2.034 2.943 16.921 95 9.905 2.040 2.952 10.303 10.973 10.053 2.055 2.974 17.109 100 10.106 2.060 2.982 10.411 11.088 10.258 2.075 3.005 17.288 (a) (b)

FIG. 1. Examples of convergence.共a兲 N=1000, four-parameter case, and共b兲 N=40 000, three-parameter case.

FIG. 2. 共Color online兲 The comparisons of ␾₃共z兲 for three-parameter 共solid line兲 and four-parameter 共dashed line兲 results with N = 1000, 10 000, 50 000, and 100 000. The horizontal axis is z in units of L = 0.522 61␮m, vertical axis is␾₃共z兲.

FIG. 3. 共Color online兲 The comparisons of momentum space ␾₃共k_z兲 for three-parameter 共solid line兲 and four-parameter 共dashed line兲 re-sults with N = 1000, 10 000, 50 000, and 100 000. The horizontal axis is kz, and the vertical axis is

role of dipole-dipole effect will emerge from the adjustment. In Fig.2, we plot the comparison of order parameter␾3共z兲 for four- and three-parameter calculations with N = 1000, 10 000, 50 000, and 100 000. The four-parameter results are combinations of two Gaussian functions. The difference in order parameters at z = 0 is mainly due to the result of the dipole-dipole effect. Also, Fig. 3 shows the corresponding order parameter in momentum space. There are negative val-ues of␾3共kz兲 for four-parameter results which are very

dif-ferent from three-parameter results. Corresponding to the ex-perimental case with 100 00052Cr atoms trapped and aligned in the y direction关11兴, we plot the density as a function of y

and z共see Fig. 4兲. The trap is far more loosely confined in

the z direction and the dipole-dipole repulsions among atoms exhibit the notable double-peak structure in density. This is a

very unique feature of dipolar BEC that does not appear in BEC systems without the long-ranged interaction.

We show in Fig.5the chemical potentials from both cal-culations. The chemical potential can be derived as

␮=

冕

冉

ប 2 2m兩ⵜ⌿兩 2_{+ V} ext共rជ兲兩⌿兩2+ Ngc兩⌿兩4+ NVdd共rជ兲兩⌿兩2

冊

d3rជ. 共12兲 We can see that the there are visible differences for double-peak and single-double-peak results for all N. The double-double-peak chemical potential is smaller than single-peak for each N. It means variationally that the four-parameter results are better than the three-parameter ones.

One question that naturally happens is how reliable is the simple double-peak trail function? To shed some light on this, we perform many-Gaussian trail functions for the cases of N = 1000 and N = 105. We make an expansion with three sets of four-parameter Gaussian trial functions:

⌿共x,y,z兲 =␾1共x兲␾2共y兲␾3共z兲, ␾1共x兲 =

兺

i=1 3 aiexp

冉

− x2 2␣i2

冊

, ␾2共y兲 =

兺

i=1 3 biexp

冉

− y2 2␤i2

冊

, ␾3共z兲 =

兺

i=1 3 ci

冋

exp

冉

− 共z −␥i兲2 2␦i2

冊

+ exp

冉

−共z +␥i兲 2 2␦i2

冊

册

. 共13兲 We obtain the value of total energy as 2.043 83 compared to the single set of four-parameter 2.0443 for N = 1000, and

to-FIG. 4.共Color online兲 The three-dimensional density plot for the experimental case. The number of atoms is 100 000 and trap fre-quencies and scattering length were described in Ref. 关11兴. The

dipoles are aligned along the y direction. The coordinates are in units of

冑

ប/m␻y= 0.522 61␮m.

FIG. 5. The chemical potential as a function of number of con-densate atoms. Solid dots are results of four-parameter calculations while empty dots are three-parameter results.

FIG. 6. The interaction energy of the dipolar BEC vs the ratio of

gc/ gd. For a ratio less than 3.2, the interaction energy becomes

negative and instability may happen. The parameter g_dis fixed with the atomic magnetic dipole moment of a chromium atom. gc is

tal energy of 10.0708 vs 10.1063 for N = 100 000. Thus we conclude that the four-parameter model is satisfactory for dipolar BEC.

With the method described above, we now turn to the instability problem of the dipolar BEC system. As shown in 关1兴, the instability will occur when the mutual atomic

inter-action becomes attractive. Then due to the increase of three-body collisional loss, the BEC may collapse. In the 52Cr dipolar system, the contact interaction energy Ec

=1₂Ngc兰兩⌿兩4d3rជ is repulsive with positive gc, while the

dipole-dipole energy term Eddis attractive. A predicted

crite-rion of the instability, gc/ gd⬍4␲/ 3 can be found in关2,6兴 for

instance. Since the scattering length can be tuned by the Feshbach resonance technique, we perform our calculation with experimental gdand adjust the scattering length near the

criterion region. In Fig.6, we show the results with 100 000 atoms共which is in the Thomas-Fermi regime兲. The total in-teraction energy Ec+ Edd becomes negative for gc/ gd⬍3.2.

This is quite close to the estimated criterion condition while the calculation here is in realistic experimental parameters.

IV. DISCUSSIONS AND CONCLUSIONS

In the past few years, simple single-peak Gaussian trial functions are used to model the order parameter of dipolar BEC. From our study, we find that the double-peak order parameter is a special property of dipolar BEC. We present a simple but efficient method to perform the variational calcu-lation of its ground state. Our simucalcu-lation results agree with the real experimental parameters. Due to the relatively larger trap potential and contact potential in comparison to the dipole-dipole potential, the double-peak property does emerge but is not drastically dominant. For the future dipolar systems, the dipole-dipole effect may become dominant if the physical conditions are adjusted. The method can be ap-plied to future systems straightforwardly.

ACKNOWLEDGMENT

T.F.J. acknowledges the support from the National Sci-ence Council of Taiwan under Contract No. NSC94-2112-M009-018 and from the MOE/ATU Program.

关1兴 For a general review, see, F. Dalfovo, S. Giorgini, L. P. Pitae-vskii, and S. Stringari, Rev. Mod. Phys. 71, 463共1999兲. 关2兴 K. Gorál, K. Rzążewski, and T. Pfau, Phys. Rev. A 61,

051601共R兲 共2000兲; J.-P. Martikainen, M. Mackie, and K. A. Suominen, ibid. 64, 037601共2001兲; P. O. Schmidt, S. Hensler, J. Werner, A. Griesmaier, A. Gorlitz, T. pfau, and A. Simoni, Phys. Rev. Lett. 91, 193201共2003兲.

关3兴 S. Yi and L. You, Phys. Rev. A 61, 041604共R兲 共2000兲. 关4兴 L. Santos, G. V. Shlyapnikov, P. Zoller, and M. Lewenstein,

Phys. Rev. Lett. 85, 1791共2000兲; U. R. Fischer, Phys. Rev. A 73, 031602共R兲 共2006兲.

关5兴 K. Gorál and L. Santos, Phys. Rev. A 66, 023613 共2002兲. 关6兴 P. M. Lushnikov, Phys. Rev. A 66, 051601共R兲 共2002兲. 关7兴 K. Gorál, L. Santos, and M. Lewenstein, Phys. Rev. Lett. 88,

170406共2002兲.

关8兴 S. Yi, L. You, and H. Pu, Phys. Rev. Lett. 93, 040403 共2004兲; L. Santos and T. Pfau, e-print cond-mat/0510634.

关9兴 D. H. J. O’Dell, S. Giovanazzi, and C. Eberlein, Phys. Rev.

Lett. 92, 250401共2004兲; C. Eberlein, S. Giovanazzi, and D. H. J. O’Dell, Phys. Rev. A 71, 033618共2005兲.

关10兴 A. Griesmaier, J. Werner, S. Hensler, J. Stuhler, and T. Pfau, Phys. Rev. Lett. 94, 160401共2005兲; J. Werner, A. Griesmaier, S. Hensler, J. Stuhler, T. Pfau, A. Simani, and E. Tiesinga,

ibid. 94, 183201共2005兲; J. Stuhler, A. Griesmaier, T. Koch,

M. Fattori, T. Pfau, S. Giovanazzi, P. Pedri, and L. Santos,

ibid. 95, 150406共2005兲.

关11兴 A. Griesmaier, J. Stuhler, and T. Pfau, Appl. Phys. B: Lasers Opt. 82, 211共2006兲.

关12兴 W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flan-nery, Numerical Recipes in FORTRAN, 2nd ed. 共Cambridge Univ. Press, New York, 1992兲.

关13兴 Kerson Huang, Statistical Mechanics, 2nd ed. 共Wiley, New York, 1987兲, Sec. 10.5.

关14兴 S. Yi and L. You, Phys. Rev. Lett. 92, 193201 共2004兲. 关15兴 S. Giovanazzi, A. Görlitz, and T. Pfau, Phys. Rev. Lett. 89,