Study of implosion in an attractive Bose-Einstein condensate

TL;DR

This paper investigates the collapse behavior of an attractive Bose-Einstein condensate using analytical and numerical solutions of the Gross-Pitaevskii equation, focusing on critical atom numbers and stability dynamics.

Contribution

It provides a combined analytical and numerical analysis of implosion phenomena in attractive BECs, including a new calculation of the critical atom number using variational methods.

Findings

Determined the critical atom number for stability.

Analyzed collapse dynamics through numerical simulations.

Validated analytical predictions with numerical results.

Abstract

By solving the Gross-Pitaevskii equation analytically and numerically, we reexamine the implosion phenomena that occur beyond the critical value of the number of atoms of an attractive Bose-Einstein condensate (BEC) with cigar-shape trapping geometry. We theoretically calculate the critical number of atoms in the condensate by using Ritz's variational optimization technique and investigate the stability and collapse dynamics of the attractive BEC by numerically solving the time dependent Gross-Pitavskii equation.

Equations

$$\left[-i\hslash \frac{\partial }{\partial t}-\frac{{\hslash }^2}{2m}{\nabla }^2+V(\tilde{\text{r}})+gN|\Psi (\tilde{\text{r}},\tau ){|}^2\right]\Psi (\tilde{\text{r}},\tau )=0.$$

$$\left[-i\frac{\partial }{\partial t}-{\nabla }^2+\frac{1}{4}({\rho }^2+{\lambda }^2{z}^2)+\frac{8\pi N{a}_s}{l}{\left|\varphi \right|}^2\right]\varphi =0,$$

$$\varphi (\rho ,z,t)=u(\rho )\psi (z,t).$$

$$i\frac{\partial \psi }{\partial t}+\frac{{\partial }^2\psi }{\partial z^2}-{\lambda }^2{z}^2\psi -\frac{4\pi N{a}_s}{l}|\psi {|}^2\psi =0,$$

$${\int }_{-\infty }^{\infty }|\psi {|}^{2}dz=N/\pi .$$

$$\delta \int \mathcal{L}(z,t,\psi ,{\psi }^{\ast },{\psi }_z,{\psi }_z^{\ast },{\psi }_t,{\psi }_t^{\ast })dt=0.$$

$$\mathcal{L}(t)=$$

$$\psi (z,t)=A(t)\exp \left[-\frac{z^2}{2R(t{)}^2}+\frac{i}{2}\beta (t)z^2\right],$$

$$<\mathcal{L}(t)>$$

$$\delta \int <\mathcal{L}>dt=0.$$

$$\frac{\partial <\mathcal{L}>}{\partial A^{\ast }}=$$

$$\frac{\partial <\mathcal{L}>}{\partial A}=$$

$$\frac{\partial <\mathcal{L}>}{\partial R}=$$

$$\frac{\partial <\mathcal{L}>}{\partial \beta }=$$

$$\frac{\partial }{\partial t}(RA{A}^{\ast })$$

$$R|A{|}^2$$

$$\beta$$

$$\frac{1}{2}{\left(\frac{dR}{dt}\right)}^2+2{\lambda }^2{R}^2+\frac{\sqrt{2}\pi N{a}_s}{l}\frac{1}{R}+\frac{2}{R^2}$$

$$V(R)$$

$$\frac{dV(R)}{dR}$$

$$4{\lambda }^2{R}^4+|P|R-4$$

$$4{\lambda }^2{R}^4-2|P|R+12$$

$$R$$

$$|P|$$

$$N$$

$$N_c$$

Full text

1 Introduction

Since the experimental realization of Bose-Einstein condensates (BECs) two decades ago Anderson1995 ; Davis1995 ; Bradley1995 , the field of quantum gases has attracted much theoretical and experimental interest. At absolute zero temperature the properties of a condensate are well described by the nonlinear, mean-field Gross-Pitaevskii equation (GPE) for both repulsive and attractive inter-atomic interactions, see the early comprehensive reviews Dalfovo1999 ; Leggett2001 ; Kevrekidis2004 ; Abdullaev2005 ; for more recent overviews in the area of ultracold quantum Bose and Fermi gases and related phenomena, see, for instance, Refs. Giorgini2008 . The nonlinear terms in the GPE arise from the interactions between atoms in the condensate, characterized by their two-body s-wave scattering length, $a_{s}$. The static and dynamical properties of trapped BEC crucially depend on the sign and strength of the interatomic interactions. Both the strength and the sign of the interaction can be controlled by varying the scattering length $a_{s}$ of atoms near the Feshbach resonance FBR . The stability of the condensate under magnetic traps has been studied both numerically and analytically Ruprecht ; Dalfovo ; attractive . If the interaction is repulsive, the BEC is stable and its size and number of particles have no fundamental limit. If it is attractive, only a limited number of atoms can form a condensate. Moreover, the BEC was predicted to be metastable, which means that it remains stable for some finite time, only when the number of atoms is below some critical number ($N_{c}$). The estimated critical number for 7Li atoms is about $10^{3}$ Bradley1997 . Beyond $N_{c}$, we expect the BEC to ‘collapse’ or ‘implode’ due to the strong attactive interaction. In particular, in a trapped attractive gas whose number of particles or the strength of the inter-particle interaction exceeds a critical value, the kinetic energy cannot balance the (negative) interaction energy PRA83 ; Shuryak . The attraction brings the bosons so tightly close such that the spatial extension of the wave function of the system shrinks to a point and the condensate eventually implodes PRA83 ; PhysRevA82 . A BEC can avoid implosion only as long as the number of atoms in it, $N$, is less than a critical value $N_{c}$. Beyond $N_{c}$, we show in our theoretical study that if we increase the number of atoms, the BEC first explodes, then implodes as expected, and subsequently a sudden explosion of atoms occurs. The explosion followed by an implosion has been already studied experimentally. This ‘explosion’ actually corresponds to a small amount of energy by normal standards. The explosion of a collapsing BEC was termed as ‘bosenova’ Donley2001 .

The study of collapsing phenomena of attractive BECs has been the subject of much interest for the past two decades; see, for instance, Refs. Ruprecht ; Dalfovo ; attractive ; Bradley1997 ; PRA83 ; Shuryak ; PhysRevA82 ; Donley2001 ; Melo ; talukdar ; HU2000 . It refers to the situation when strong self-focusing of a beam leads to a catastrophic increase (blow-up) of its intensity in a finite time or after a finite propagation distance. In this context, the collapse process is a fast, collective phenomenon consisting of the destruction of a multi-particle system happening abruptly on the timescale that governs the usual dynamics. One example is the collapse of the gravitational core, initiating a supernova. Tailoring both the external confining potential and the interaction between the atoms allows us to control the properties of the condensate. It is an ideal system to study not only problems of condensed matter physics, but also the dynamics of a collapse as well. The collapse of an attractive BEC of ${}^{7}{Li}$ or ${}^{85}Rb$ atoms has been observed by various experimental groups Bradley1997 ; roberts (collapsing condensates were first observed in 7Li NJP ) and theoretically analyzed by many authors for various external trapping potentials such as, single-well potential Ruprecht ; Sackett , potential without axial confinement Garcia ; Carr , toroidal confinement Parola , double-well potential Sakellari , and periodic potential Adhikari . In recent past, Ghosh Ghosh studied the duality-symmetry, which are abundant in different branches of physics and astrophysics (see Ghosh and refs. therein) and showed that invariance under the duality-symmetry leads to explosion-implosion duality in one- and two-dimensional Bose-Einstein condensation without the harmonic trap.

In this work, we study the sequence of implosion and explosion phenomena of an attractive BEC beyond the $N_{c}$, by using the time-dependent, nonlinear, mean-field GPE. First, we calculate theoretically the $N_{c}$ for the BEC by solving the quasi one-dimensional (Q1D) GPE taking recourse to the use of Ritz’s variational optimization technique Anderson1983 . Then, we study the collapsing phenomenon of an attractive BEC beyond $N_{c}$ through extensive numerical simulations of the GPE. The organization of the present paper is as follows. In Sec. 2, we present a brief overview of the nonlinear mean-field model. In Sec. 3, we derive the equation of motion of the Q1D system to find out the critical number of atoms through variational approximation (VA) method and discuss the stability of the condensate. Then, in Sec. 4 we numerically study the collapsing dynamics of the BEC by analyzing the GPE through split-step Crank-Nicholson (SSCN) method Muruganandam2009 . Finally, we give the concluding remarks in Sec. 5.

2 Nonlinear mean-field model

The GPE can be used at low temperatures, to explore the macroscopic behavior of the system. The time-dependent BEC wave function $\Psi(\tilde{\textbf{r}},\tau)$ can be described by the following mean-field nonlinear GP equation Dalfovo1999 ,

[TABLE]

Here $\nabla^{2}=\frac{\partial^{2}}{\partial x^{2}}+\frac{\partial^{2}}{\partial y^{2}}+\frac{\partial^{2}}{\partial z^{2}}$, $g=4\pi\hbar^{2}a_{s}/m$, $a_{s}$ is the atomic s-wave scattering length, which is negative or positive for attractive or repulsive interaction between atoms in the condensate, $m$ is the mass of a single bosonic atom and $N$ is the number of atoms in the condensate. The external magnetic trap potential may be written as $V(r)=\frac{1}{2}m\omega^{2}(\nu^{2}x^{2}+\kappa^{2}y^{2}+\lambda^{2}z^{2})$, where $\omega_{x}\equiv\nu/\omega$, $\omega_{y}\equiv\kappa/\omega$ and $\omega_{z}\equiv\lambda/\omega$ are the angular frequencies of the trap in the $x,y$, and $z$ directions, respectively. The normalization condition is $\int|\Psi(\tilde{\textbf{r}},\tau)|^{2}=1$.

To study the collapse dynamics beyond the critical number of atoms, we are interested to work with quasi one-dimensional (Q1D) cigar-shaped BEC dispersing along the $z$ direction. In the following, we consider equation (1) in a geometry in which the trapping potential in $z$ is much weaker than the corresponding potential in $\tilde{\rho}=\sqrt{(\tilde{x}^{2}+\tilde{y}^{2})}$. Further, we make the transformation of variables as $\rho=\sqrt{2}\tilde{\tilde{\rho}}/l$, $z=\sqrt{2}\tilde{z}/l$, $t=\tau\omega$, $l=\sqrt{\hbar/(m\omega)}$ and $\phi(\rho,z,t)=\Psi(\tilde{\rho},\tilde{z},\tau)(l^{3}/2)^{1/2}$. Then, the GP equation (1) becomes Muruganandam2009 ,

[TABLE]

We assume a separable ansatz for the solution of equation (2) such that Garcia

[TABLE]

Substituting Eq. (3) in Eq. (2), after simplification one can get the following Q1D equation

[TABLE]

with

[TABLE]

Equation (4) represents a desired form of the evolution equation in which the atom-atom interaction is characterized by a negative s-wave scattering length. The realistic 1D limit in (4) is not a true 1D system because this equation involves the effect of transverse degrees of freedom through $z$ and $l$.

3 Variational calculation of critical number $N_{c}$

In this section by using Ritz’s variational optimization technique Anderson1983 , we shall calculate the critical number of atoms ($N_{c}$) to analyze the collapsing dynamics of the matter waves beyond $N_{c}$. The Eq. (4) can be restated as the following variational problem

[TABLE]

with the Lagrangian density written as

[TABLE]

Now we use the variational approach with the Gaussian trial wave function for the solution of Eq. (4) Garcia ; Sabari2010 :

[TABLE]

where $A(t)$, $R(t)$, and $\beta(t)$ are the time dependent amplitude, width, and chirp, respectively. The initial condensate at rest will have $dR(t)/dt=0$. The trial wave function equation (8) is substituted in the Lagrangian density and the averaged effective Lagrangian is calculated by integrating the Lagrangian density as $L_{eff}=\int^{\infty}_{\infty}\mathcal{L}\,dz$ to write

[TABLE]

The reduced variational principle

[TABLE]

gives a set of coupled ordinary differential equations for the parameters of our trial solution. Let us now obtain the variational equations for the Gaussian parameters $A(t)$, $A(t)^{*}$, $R(t)$, and $\beta(t)$, which follow from the vanishing conditions of $\frac{\partial<\mathcal{L}>}{\partial A^{*}}$, $\frac{\partial<\mathcal{L}>}{\partial A}$, $\frac{\partial<\mathcal{L}>}{\partial R}$, and $\frac{\partial<\mathcal{L}>}{\partial\beta}$. These equations are given by

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

From equations (11) and (12) we get

[TABLE]

and consequently

[TABLE]

where the constant $Q$ is related to the number of particles in the condensate since the value of the integral (5) is $\sqrt{\pi}R|A|^{2}$. Combining the equations (14) and (16), we get

[TABLE]

Equations (16) and (17) clearly show that if we can derive a method to calculate the values of $R$, the other parameters of the condensate will be automatically determined. It is fortunate that we are able to write a second-order ordinary differential equation from (11), (12), (13), and (17). This gives a first integral of the form

[TABLE]

where $E$ is the the constant of integration. The equation for $R$ in Eq. (18) is related to the motion of a particle of unit mass in a potential field $V(R)$ of the form

[TABLE]

The constant of the motion $E$, i.e. the total energy of the particle, can be determined by the initial conditions of the second-order differential equation from which Eq. (18) has been extracted. Now it is easy to solve Eq. (18) and to look for the dynamics of the condensate. However, the analysis of the equilibrium point obtained from the extremum of $V(R)$ written as

[TABLE]

can give some illuminating results. For bright solitons the nonlinear interaction is attractive and the scattering length $a_{s}<0$. In this case we shall use $P=-|P|$ and carry out the subsequent analysis by using only the numerical values of $a_{s}$. We shall make use of Eq. (20) to derive a simple physical picture for the collapse dynamics of bright solitons when the trap of the BEC is relaxed in one direction. From Eqs. (19) and (20) with $P=-|P|$ we get

[TABLE]

The equilibrium point determined by Eq. (21) should be a minimum for our system to support a soliton solution and the condition for minimum ($\frac{d^{2}V(R)}{dR^{2}}>0$) gives

[TABLE]

where $\mu>0$. Eliminating $|P|$ from Eqs. (21) and (22) we find that

[TABLE]

is a particular solution of Eqs. (21) and (22). From Eq. (23) and Eq. (21) or Eq. (22) we get

[TABLE]

The form of equation (24) imposes a further restriction on the values of $\mu$ than that given in equation (23) and sets a lower bound for it. Using $\mu=\gamma\lambda^{2}$ we write Eq. (24) in the form

[TABLE]

Thus non-zero values of $P$ will be obtained for $\gamma>16$ only.

For $\gamma>16$ the interaction term vanishes and the corresponding GPE becomes linear and the soliton formation becomes impossible. From Eqs. (19) and (25) we obtain the expression

[TABLE]

for the number of atoms that are present in the system. In Fig. 1 we plot the potential $V(R)$ in (19) as a function of $R$; here we have used the value of $N$ from Eq. (26). In this figure we show three curves represented by $V_{23}(R)$, $V_{24}(R)$, and $V_{25}(R)$ corresponding to $\gamma=23,\,24,$ and $25$, respectively. A common feature of all these potentials is that each of them exhibits a minimum. The curve for $V_{25}$ represents a potential well. The minimum of the well is negative. A mechanical analogy suggests a solution that oscillates between the zeros of $V_{25}$. In this situation the BEC soliton will become unstable and leads to a collapse.

A similar situation arises for other values of $\gamma>24$. For $\gamma=24$ the potential well degenerates into a single point such that $V_{24}$ touches the axis at a particular point, where the potential has a stable minimum. Understandably, a particle released at this point will stay there. In the present context this implies that for our chosen value $\gamma=24$ the BEC soliton will be critically stable. The attractive BEC is stable, only when the number of atoms is below some critical number of atoms calculated from Eq. (26)

[TABLE]

For the ${}^{7}Li$ atom, mass $m=1.1524\times 10^{-26}$ Kg, the scattering length $a_{s}=(-27.6\pm 0.5)a_{0}$ Bradley1997 and $=(-27.3\pm 0.8)a_{0}$ Bradley1995 , with Bohr radius $a_{0}=0.529\times 10^{-10}$ m, and using the frequencies of the trap in the radial direction $\omega_{\rho}=433$ Hz and axial direction $\omega_{z}=39$ Hz for the cigar-shaped geometry, and correspondingly $\bar{\omega}=194.095$ and $\lambda=0.20093$, we have calculated the critical number of atoms $(N_{c})$= 1667 from Eq. (27). Our calculated result is closely matched with the critical value $N_{c}=1400$ for this geometry reported earlier by the experiment Bradley1995 ; Bradley1997 ; Bradley2 and by the theory wadati . The $N_{c}$ can change with respect to experimental parameters i.e., different trapping frequencies and interaction between the atoms in the condensate are related to the s-wave scattering length tuned by Feshbach resonance technique FBR . Due to these factors the initial wavefunction of the condensate is more spread or confined, and such a configuration imposes a more or less severe restriction on the collapse and also change the $N_{c}$ wadati .

4 Numerical results

We have numerically simulated the stability and collapse dynamics of Q1D BEC that supports our analytical study of the critical number $N_{c}$. The experimental realization of cigar-shaped attractive BECs has been possible under the strong transverse binding which, in the case of weak or no axial binding, imply one dimensional BECs. We study the stability and the evolution of a collapsing cigar-shaped BEC by solving numerically the time-dependent three-dimensional (3D) GPE (1) with strong transverse binding ($\omega_{\rho}=433$ Hz) and less axial trapping ( $\omega_{z}=39$ Hz), through the SSCN method Muruganandam2009 . The advantage of the present split-step procedure is that the nonlinear and other linear non-derivative terms can be treated very precisely and this improves the accuracy and stability of the method compared to other methods. We have used spatial grid points $40\times 40\times 80$; the typical space and time steps for discretization are 0.01 and 0.0001.

Figure 2 depicts the stability of the attractive condensate for the number of atoms present in the system both below $N_{c}$ and around $N_{c}$. Here, in the upper plot (see Fig. 2), we have shown that the initial (dashed line) and final (solid line) profiles illustrate the stability of the condensate below the critical number $N(<N_{c})$ = $1600$. In the lower left-hand-side subplot (see Fig. 2), we have presented the time evolution of the amplitude $|\psi|$ of attractive BEC, which shows the gradual shrinking of the width of the condensate with respect to time, but still it is stable for more than $t=20$ ms. In this case, the system is spatially confined and when the number of BEC atoms is below a certain critical value $N_{c}$, the zero-point motion of the atoms serves as a kinetic obstacle against collapse, which allows to form a metastable BEC.

But, in Fig. 3, we have shown the sudden explosion of attractive BEC followed by an implosion during the evolution after $t=20$ ms for $N=1660$. Due to high attractive interaction of the BEC atoms the condensate implodes and then a sudden explosion occurs because locally near the center of BEC where the atomic density exceeds a certain critical value, the centripetal force then weakens, and the atoms that gathered in this narrow central region are ejected due to the ’quantum pressure’ arising from the uncertainty principle. Here, interestingly, we have shown an explosion followed by an implosion. Further, when we increase the number of atoms in the condensates, say, $N=1700$, the condensate collapses earlier at $t=14ms$ and explodes; this phenomenon is clearly shown in Fig. 4. Then we have considered up to $N=2000$ atoms and we have shown in Fig. 5 another type of explosion pattern.

We have shown in Fig. 6 the amplitude/contour plots corresponding to the final condensate for $N=$ (left$\rightarrow$ $1600$), (middle$\rightarrow$ $1660$ and $1700$) and (right$\rightarrow$ $2000$), respectively. Below the critical number of atoms we can see that the condensates are in the metastable states. As we increase the number of the atoms to critical value, the condensate is shrinking with time and explodes into four major peaks followed by its implosion. If we increase the number of atoms even more, another type of explosion phenomenon occurs giving rise to a different pattern containing a central high density peak surrounded by small density peaks (see the right panel in Fig. 6).

5 Conclusion

In conclusion, we have shown theoretically the possibility of a subsequent explosion after implosion of an attractive BEC by solving numerically the time-dependent GPE (1) through the SSCN method. We have calculated the critical number of atoms through the VA method and investigated the stability of the condensate by numerical techniques. The condensate collapses when the number of atoms exceeds the critical value. Beyond $N_{c}$, the attractive BEC first exhibits explosion followed by implosion and then it exhibits a subsequent explosion. Our results are in good agreement to those obtained in the Bosenova experiment at JILA Donley2001 and at Rice University Bradley1997 ; HU2000 that the condensates with a particle number exceeding the critical value are unstable against collapse.

Acknowledgement

S.S. thanks University Grants Commission (UGC) for offering support through a Post-Doc under the Dr. D.S. Kothari Post Doctoral Fellowship Scheme. K.P. thanks to the DST-DFG, DST, CSIR and UGC, Government of India, for the financial support through major projects.