Matter-wave solutions in the Bose-Einstein condensates with the harmonic and Gaussian potentials

TL;DR

This paper derives exact matter-wave solutions for the quasi-one-dimensional Gross-Pitaevskii equation with combined harmonic and Gaussian potentials, analyzing their stability and potential experimental relevance.

Contribution

It introduces new exact solutions for the GP equation with modulated potentials and nonlinearity, including stability analysis and parameter regimes.

Findings

Several families of exact solutions are identified.

Some solutions are found to be stable through numerical analysis.

Parameter regimes for stability are characterized.

Abstract

We study exact solutions of the quasi-one-dimensional Gross-Pitaevskii (GP) equation with the (space, time)-modulated potential and nonlinearity and the time-dependent gain or loss term in Bose-Einstein condensates. In particular, based on the similarity transformation, we report several families of exact solutions of the GP equation in the combination of the harmonic and Gaussian potentials, in which some physically relevant solutions are described. The stability of the obtained matter-wave solutions is addressed numerically such that some stable solutions are found. Moreover, we also analyze the parameter regimes for the stable solutions. These results may raise the possibility of relative experiments and potential applications.

Equations

$$i\hslash \frac{\partial \Psi }{\partial t}=\left[-\frac{{\hslash }^2}{2m}\frac{{\partial }^2}{\partial x^2}+{\mathcal{V}}_{\mathrm{ext}}(x,t)+{\mathcal{G}}_{1D}(x,t)|\Psi {|}^2+i\vartheta (t)\right]\Psi ,$$

$$a_s(x,t)=\hat{a}\left(1+\frac{\Delta }{B_0-B(x,t)}\right),$$

$$i\frac{\partial \psi }{\partial t}=\left[-\frac{1}{2}\frac{{\partial }^2}{\partial x^2}+V(x,t)+g(x,t)|\psi {|}^2+i\Gamma (t)\right]\psi ,$$

$$\mathcal{L}=i\left(\psi \frac{\partial {\psi }^{\ast }}{\partial t}-{\psi }^{\ast }\frac{\partial \psi }{\partial t}\right)+{\left|\frac{\partial \psi }{\partial x}\right|}^2$$

$$+g(x,t)|\psi {|}^4+2[V(x,t)+i\Gamma (t)]|\psi {|}^2,$$

$$\psi (x,t)=\rho (x,t)e^{i\phi (x,t)}\Phi [X(x,t)]$$

$$\epsilon \Phi (X)=-\frac{d^2\Phi (X)}{d{X}^2}+G|\Phi (X){|}^2\Phi (X),$$

$$({\rho }^2{)}_t+({\rho }^2{\phi }_x{)}_x-2\Gamma (t){\rho }^2=0,({\rho }^2{X}_x{)}_x=0,$$

$$g(x,t)=\frac{G{X}_x^2}{2{\rho }^2},X_t+{\phi }_x{X}_x=0,$$

$$V(x,t)=\frac{1}{2}\left(\frac{{\rho }_{xx}}{\rho }-{\phi }_x^2-\epsilon X_x^2\right)-{\phi }_t.$$

$$\rho (x,t)=e^{\int \Gamma (t)dt}\sqrt{\frac{\gamma }{F^{\prime }(\xi )}},$$

$$\phi (x,t)=-\frac{{\gamma }_t}{2\gamma }{x}^2-\frac{{\delta }_t}{\gamma }x+\alpha (t),$$

$$g(x,t)=\frac{\gamma G}{e^{2\int \Gamma (t)dt}}{F}^{\prime 3}(\xi ),$$

$$V(x,t)=\frac{1}{2}\left(\frac{{\rho }_{xx}}{\rho }-{\phi }_x^2-\epsilon {\gamma }^2{F}^{\prime 2}(\xi )\right)-{\phi }_t,$$

$$g(x,t)=\frac{\gamma G}{2e^{2\int \Gamma (t)dt}}{e}^{-3{\xi }^2}$$

$$V(x,t)=\frac{p(t)}{2}{x}^2+q(t)x+h(t)-\frac{\epsilon }{2}{\gamma }^2(t)e^{-2{\xi }^2},$$

$$\displaystyle\begin{array}[]{l}p(t)=\gamma^{4}e^{8\int\Gamma(t)dt}+(\gamma\gamma_{tt}-2\gamma_{t}^{2})/\gamma^{2},\vspace{0.1in}\cr q(t)=\gamma^{3}\delta e^{6\int\Gamma(t)dt}+(\gamma\delta_{tt}-2\gamma_{t}\delta_{t})/\gamma^{2},\vspace{0.1in}\cr h(t)=\gamma^{2}e^{4\int\Gamma(t)dt}(1+\delta^{2})/2-\delta_{t}^{2}/2\gamma^{2}-\alpha_{t}.\end{array}$$

$$V(x,t)=\frac{p(t)}{2}{x}^2-\frac{\epsilon }{2}{\gamma }^2(t)e^{-2{\gamma }^2(t)x^2}$$

$$X(x,t)={\int }_{-\infty }^{\gamma (t)x}{e}^{-s^2}ds,$$

$$\rho (x,t)=\sqrt{\gamma (t)}\exp \left[{\gamma }^2(t)x^2/2+\int \Gamma (t)dt\right],$$

$$g(x,t)=\frac{G}{2}\gamma (t)\exp \left[-3{\gamma }^2(t)x^2-2\int \Gamma (t)dt\right],$$

$$\phi (x,t)=-\frac{{\gamma }_t}{2\gamma }{x}^2+\frac{1}{2}\int {\gamma }^2(t)dt.$$

$$i\frac{\partial \psi }{\partial t}=\left[-\frac{1}{2}\frac{{\partial }^2}{\partial x^2}+V(x)+g(x)|\psi {|}^2\right]\psi .$$

$${\psi }_{sta}(x,t)=\rho (x)e^{i\phi (t)}\Phi [X(x)]$$

$$X(x)={\int }_{-\infty }^{\sqrt{\omega }x}{e}^{-s^2}ds\equiv \frac{\sqrt{\pi }}{2}\left[1+\mathrm{erf}(\sqrt{\omega }x)\right],$$

$$\rho (x)={\omega }^{\frac{1}{4}}{e}^{\frac{\omega x^2}{2}},\phi (t)=\frac{\omega t}{2}$$

$$g(x)=\frac{G\sqrt{\omega }}{2}{e}^{-3\omega x^2},$$

$$V(x)=\frac{{\omega }^2}{2}{x}^2-\frac{\epsilon \omega }{2}{e}^{-2\omega x^2},$$

$$-\frac{d^2\Phi (X)}{d{X}^2}+G|\Phi (X){|}^2\Phi (X)=0,$$

$${\Phi }_{sn}[X(x)]=\frac{\sqrt{2}\nu k}{\sqrt{G}}\mathrm{sn}[\nu X(x),k],$$

$$u={\int }_0^{\sigma }\frac{ds}{\sqrt{1-k^2{\sin }^{2}s}},k\in (0,1)$$

$$\mathrm{sn}(u,k)=\sin \sigma ,$$

$$K(k)={\int }_0^{\pi /2}\frac{ds}{\sqrt{1-k^2{\sin }^{2}s}}.$$

$$\mathrm{sn}[2nK(k),k]\equiv 0,k\in (0,1),n=0,1,2,...$$

$$\nu =\frac{2n}{\sqrt{\pi }}K(k),$$

$${\psi }_{sn,n}(x,t)=\frac{2\sqrt{2\sqrt{\omega }}nkK(k)}{\sqrt{\pi G}}{e}^{\frac{\omega x^2}{2}}\mathrm{sn}[{\theta }_n(x),k]e^{\frac{i\omega t}{2}},$$

Full text

I Introduction

The nonlinear Schrödinger (NLS) equation (alias the Gross-Pitaevskii (GP) equation in Bose-Einstein condensates) is one of the most useful physical models, which actually arises in many physics fields, such as nonlinear quantum field theory, condensed matter, plasma physics, nonlinear optics, photonics, fluid mechanics, semiconductor electronics, phase transitions, biophysics, star formation, econophysics, and so on (see, e.g., b1 ; b2 ; b3 ; b4 ; b5 ; b6 ). The first exact solutions of the one-dimensional (1D) NLS equation were obtained by using the inverse scattering method bb01 , and over these years there have been many significant contributions to the studies of exact solutions for some extensions of the NLS equation bb02 ; bb03 ; bb04 ; j8 .

Recently, exact localized solutions of the GP equations with varying coefficients have drawn much attention due to their potential applications in the theory of Bose-Einstein condensates (BECs) j5 ; j55 ; j6 ; j7 ; fr . It is in general difficult to analytically solve these GP equations due to their nonlinear nature, but there have been some works investigating exact solutions of the GP equations with (time, space)-modulated potential and nonlinearity, such as the harmonic potential $x^{2}$ liang ; bb07 ; bb08 ; yan11 , the van der Waals law potential $1/x^{3}$ yu9 , the harmonic potentials in the three-dimensional space yan06 , in the $d$-dimensional space vvk , as well as in the bounded two-dimensional domains yan12 .

For the usual 3D external potential $V_{\rm ext}(x,y,z)=\frac{1}{2}m(\omega_{x}x^{2}+\omega_{y}y^{2}+\omega_{z}z^{2})$ in BECs, where $m$ is the atomic mass and the trap frequencies $\omega_{x,y,z}$ along the three directions are in general different and may be used to control the shape of the condensates. The isotropic case $\omega_{x}=\omega_{y}=\omega_{\bot}\approx\omega_{z}$ denotes that the 3D BEC is almost spherical, but the other two weakly anisotropic cases $\{\omega_{\bot}<\omega_{z}\,\,{\rm or}\,\,\omega_{\bot}>\omega_{z}\}$ imply that the 3D BEC is cigar shaped or disc-shaped, respectively. Moreover the strongly anisotropic cases $\{\omega_{\bot}\gg\omega_{z}\,\,{\rm or}\,\,\omega_{\bot}\ll\omega_{z}\}$ are related to the lower-dimensional BECs, i.e., the effectively quasi-1D or quasi-2D BECs, respectively b5 . The lower-dimensional BECs in highly asymmetric traps have been verified in the study of the Bosonic and Fermionic gases low .

Here we consider a cigar-shaped BEC of a relatively low density, when the energy of two body interactions is much less than the kinetic energy in the transverse direction, i.e., $|a_{s}|\mathcal{N}\ll a_{\bot}$ with $\mathcal{N}$ being a total number of atoms and $a_{\perp}=\sqrt{\hbar/(m\omega_{\bot})}$ j55 , the macroscopic wave function of a quasi-one-dimensional condensate in BECs obeys the effective 1D GP equation j5 ; b5 ; j55 ; j6

[TABLE]

where $\Psi\equiv\Psi(x,t)$ denotes the condensate wave function, $m$ is the atomic mass, $\mathcal{V}_{\rm ext}(x,\!t)$ denotes the time- and space-modulated external confining potential, which can be usually chosen in the form of the harmonic well or the optical lattice j5 , the 1D effective coupling coefficient $\mathcal{G}_{\rm 1D}(x,t)=2\hbar\omega_{\bot}a_{s}(x,t)$ stands for a measure of the nonlinear two-body interactions between the atoms in quasi-1D geometries with strong transverse confinement with $\omega_{\bot}$ being the transverse trapping frequency and $a_{s}(x,t)$ being the time-dependent $s$-wave scattering length for elastic atom-atom collisions modulated by a Feshbach resonance j5 ; fr , in which the scattering length $a_{s}(x,t)$ is a function of the varying magnetic field $B(x,t)$:

[TABLE]

where $\hat{a}$ is the value of the scattering length far from resonance, $\Delta$ represents the width of the resonance, $B_{0}$ is the resonant value of the magnetic field, and $B(x,t)$ is external (space, time)-modulated magnetic field fr . The scattering length $a_{s}$ can be chosen as either positive or negative values corresponding to the attractive interactions (e.g., for 7Li or 85Rb in the BECs) or repulsive interactions (e.g., for 87Rb or 23Na in the BECs) between the atoms b5 ; j5 ; j55 ; j6 . The time-dependent gain or loss term $\vartheta(t)$ corresponds to the mechanism of continuously loading external atoms into the BEC (gain case) by optical pumping or continuous depletion of atoms from the BEC (loss case) gl .

Normalizing the density $|\Psi|^{2}$, length, time and energy in Eq. (1) in units of $2a_{s},\,a_{\bot}=\sqrt{\hbar/(m\omega_{\bot})},\,\omega_{\bot}^{-1}$, and $\hbar\omega_{\bot}$, respectively, we obtain a 1D effective GP equation with the space- and time-modulated potential and nonlinearity and the time-dependent gain or loss term in the dimensionless form b5 ; j5 ; j55 ; j6 ; bb08

[TABLE]

where the external potential $V(x,t)$, the nonlinearity $g(x,t)$, and the gain or loss distribution $\Gamma(t)$ are related to the functions $\mathcal{V}_{\rm ext}(x,t)$,  $\mathcal{G}_{\rm 1D}(x,t)$, and $\vartheta(t)$ in Eq. (1), respectively. Eq. (2) is associated with the Euler-Lagrange equation $\delta\mathcal{L}/\delta\psi^{\ast}=0$ in which the Lagrangian density can be written as lan ; yan06

[TABLE]

where $\psi^{\ast}\equiv\psi^{\ast}(x,t)$ denotes the complex conjugate of the complex field $\psi(x,t)$. Eq. (2) with the harmonic potential $V(x,t)=a(t)x^{2}$ has been shown to possess exact solutions only for the attractive interaction $g(x,t)<0$ bb08 .

In this paper, to translate our results into units relevant to the experiments exp00 ; exp0 ; exp1 ; exp2 , we perform this protocol for the attractive case by preparing BECs of up to $\sim 10^{4}-10^{5}$ 7Li atoms with $m=11.65\times 10^{-27}{\rm Kg}$ and $\omega_{\bot}=2\pi\times 10^{3}$Hz, in which a unity of the dimensionless space corresponds to $a_{\bot}=\sqrt{\hbar/(m\omega_{\bot})}\simeq 1.19{\rm\mu m}$ and a unity of the dimensionless time corresponds to $1.6\times 10^{4}\,{\rm s}$. For the repulsive interaction, we prepare BECs of up to $\sim 10^{5}$ 87Rb atoms with $m=1.44\times 10^{-25}{\rm Kg}$ and a unity of the dimensionless space corresponds to $a_{\bot}=\sqrt{\hbar/(m\omega_{\bot})}\simeq 0.34{\rm\mu m}$. Here, we want to point out that in the following discussion, both the space and time extents of the obtained solutions are both in the range of quasi-one-dimensional validity j5 ; b5 ; j55 ; j6 . We consider spatially localized solutions of Eq. (2) with two kinds of interactions (the repulsive interaction $g(x,t)>0$ and the attractive interaction $g(x,t)<0$) and the combination of the harmonic and Gaussian potentials $V(x,t)$. To do this, we connect exact solutions of the GP equation with space- and time-modulated potential and nonlinearity given by Eq. (2) with those of the stationary GP equation (5) by using the proper similarity transformation and some constraints about the varying coefficients. In addition, the stability of the solutions and the effects of amplitude and phase noises for the solutions are analyzed via the numerical simulations.

The paper is organized as follows. In Section II, we describe the similarity transformation and choose proper variables in order to obtain the combination of the harmonic and Gaussian potentials. In Sections III and IV, we present several families of exact solutions of the GP equation with both the stationary and the (space, time)-modulated potentials and nonlinearities. Moreover, we also study the stability of the obtained solutions numerically to find some stable solutions. In particular, we analyze the effects of both the amplitude noise and the phase noise for the stable solutions. We also give the parameter regimes for the stable solutions. These results may raise the possibility of relative experiments and potential applications in BECs and other related fields. Finally, the results are summarized in Section V.

II Similarity solutions

II.1 Similarity transformation and constraints

To seek spatially localized solutions of Eq. (2) satisfying the condition $\lim_{|x|\rightarrow\infty}|\psi(x,t)|=0$, we choose the complex field $\psi(x,t)$ in the form

[TABLE]

and require that the real-valued function $\Phi[X(x,t)]$ solves the 1D stationary GP equation with constant coefficients

[TABLE]

whose solutions can be obtained by using some transformations (see, e.g., yan06 ; yan ). Here $X\equiv X(x,t)$ denotes the similarity variable and is a real function of space $x$ and time $t$ to be determined below, $\varepsilon$ is the eigenvalue of the nonlinear stationary equation, and $G$ is a constant interaction ($G>0$ corresponds to the repulsive interaction and $G<0$ corresponds to the attractive interaction). The amplitude $\rho$ and phase $\varphi$ are both real functions of $x$ and $t$ to be determined later.

As a consequence, we find the amplitude $\rho(x,t)$, the similarity variable $X(x,t)$, the phase $\varphi(x,t)$, the external potential $V(x,t)$, the nonlinearity $g(x,t)$, and the gain or loss term $\Gamma(t)$ are required to satisfy the following system of nonlinear partial differential equations

[TABLE]

Since the solutions of the reduced equation (5) are known, thus if we can determine these variables $X(x,t),\,\rho(x,t),\,\varphi(x,t),\,g(x,t)$, and $V(x,t)$ from system (6), then we can obtain the solutions of Eq. (2) by using those of the reduced equation (5) and the similarity transformation (4).

It may be difficult to directly solve system (6). Here we only need its some nontrivial solutions. To solve system (6), we choose the similarity variable $X(x,t)$ in the form $X(x,t)=F(\xi)$, where $\xi(x,t)=\gamma(t)x+\delta(t)$ and $F(\cdot)$ is a differentiable function. Here $\gamma(t)$ is the inverse of the width of the localized solution, and $-\delta(t)/\gamma(t)$ is the position of its center of mass. After some algebra it follows from system (6) that

[TABLE]

where the prime denotes the derivative with respect to the variable $\xi$, $\alpha(t)$ is an arbitrary function of time, and the condition $\gamma/F^{\prime}(\xi)>0$ is required.

II.2 The choice of the parameters

Since the choice of the possible function $F(\xi)$ is very rich in the amplitude $\rho(x,t)$, the nonlinearity $g(x,t)$, and the external potential $V(x,t)$, considering the external potential $V(x,t)$ to be a combination of the harmonic and Gaussian traps and the possible application of our results to BECs j5 ; exp0 ; gau , we thus choose the Gaussian shaped interaction

[TABLE]

with $\xi=\gamma(t)x+\delta(t)$, which corresponds to $F^{\prime}(\xi)=\exp(-\xi^{2})$, i.e., $F(\xi)=\int^{\xi}_{-\infty}\exp(-s^{2})ds$. It follows from Eqs. (7d) and (8) that we have the external potential

[TABLE]

where we have introduced the following functions

[TABLE]

Taking $\delta(t)=0$ and $\alpha(t)=(1/2)\int\gamma(t)^{2}dt$, then we obtain the required potential, i.e., the combination of the harmonic and Gaussian potential

[TABLE]

and other variables

[TABLE]

Thus, if we choose the inverse of the width of the localized solution $\gamma(t)$ and the gain or loss term $\Gamma(t)$, then we can determine these variables $X(x,t)$, $\rho(x,t)$, and $\varphi(x,t)$ for which the solutions of Eq. (2) can be obtained from those of Eq. (5) using the similarity transformation (4).

III Stationary potential, nonlinearity and solutions

In this section we consider the simple case of the stationary nonlinearity $g(x)$ and potential $V(x)$ in Eq. (2), which means that $\gamma(t)\equiv\gamma={\rm const}$ and $\Gamma(t)\equiv 0$   (i.e., zero-gain or zero-loss term) based on Eq. (8), then Eq. (2) can be reduced to the GP equation with the external potential and nonlinearity depending only on the spatial coordinate bb07

[TABLE]

Thus, based on the results in Section II, Eq. (16) can be reduced to Eq. (5) with $X(x,t)\to X(x)$ via the following similarity transformation (c.f. Eq. (4))

[TABLE]

with the similarity variables being of the form

[TABLE]

with ${\rm erf}(\cdot)$ being the error function, and the external potential and nonlinearity in Eq. (16) satisfying the conditions

[TABLE]

from which the external potential $V(x)$ is a combination of both harmonic potential $(1/2)\omega^{2}x^{2}$ and the Gaussian potential $-(1/2)\varepsilon\omega e^{-2\omega x^{2}}$, where $\omega=\gamma^{2}$. The parameters $\omega^{2}/2$ and $|\varepsilon\omega/2|$ measure the relative strengths of the harmonic and Gaussian potentials, respectively.

Notice that when $\varepsilon=0$, the external potential $V(x)$ given by Eq. (19b) becomes the usual harmonic potential and the reduced equation (5) is changed into the ordinary differential equation

[TABLE]

whose solution has been considered for only the attractive interaction $G<0$ bb07 . In addition, the dynamical evolution of a Fermi super-fluid had been numerically studied based on the quasi-1D density-functional GP equation with a double-well potential similar to Eq. (19b) pu .

In the following, we consider exact solutions and their dynamics of Eq. (16) with the combination of the harmonic and Gaussian potential (19b) and Gaussian shaped interaction (19a) for the general case $\varepsilon\omega\not=0$.

III.1 The repulsive interaction

We have the periodic sn-wave solution of Eq. (5) with the repulsive interaction $G>0$ ($X(x,t)\to X(x)$) yan06

[TABLE]

where $k\in(0,1)$, $\nu$ is a constant, $\varepsilon=\nu^{2}(1+k^{2})>0$, $X(x)$ is given by Eq. (18a), the function ‘${\rm sn}$’ denotes the Jacobi elliptic function with $k$ being its modulus.

The Jacobi elliptic sn function is defined by using the inversion of the first kind of the incomplete elliptic integral, i.e., the first kind of the incomplete elliptic integral (also called Legendre normal form) sn

[TABLE]

generates the Jacobi elliptic sn function

[TABLE]

from which we have the limit properties: ${\rm sn}(x,k)=\tanh(x)$ for $k\to 1$, and ${\rm sn}(x,k)=\sin(x)$ for $k\to 0$. Moreover we further have another two basic Jacobi elliptic cn and dn functions: ${\rm cn}(u,k)=\cos\sigma$, and ${\rm dn}(u,k)=\sqrt{1-k^{2}\sin^{2}\sigma}>0$ from Eq. (23). These three basic Jacobi elliptic functions admit these two relations: ${\rm sn}^{2}(u,k)+{\rm cn}^{2}(u,k)=1$ and $k^{2}{\rm sn}^{2}(u,k)+{\rm dn}^{2}(u,k)=1$. In particular, when $\sigma=\pi/2$ the first kind of the ${\it incomplete}$ elliptic integral (22) becomes the first kind of the ${\it complete}$ elliptic integral

[TABLE]

These Jacobi elliptic functions have the special values ${\rm sn}(0,k)=0$,  ${\rm cn}(0,k)={\rm dn}(0,k)=1$, ${\rm sn}(K,k)=1$,  ${\rm cn}(K,k)=0,\,{\rm dn}(K,k)=\sqrt{1-k^{2}}$, and have the periodic properties ${\rm sn}(u+4K,k)=-{\rm sn}(u+2K,k)={\rm sn}(u,k)$,  ${\rm cn}(u+4K,k)=-{\rm cn}(u+2K,k)={\rm cn}(u,k)$, ${\rm dn}(u+2K,k)={\rm dn}(u,k)$. Based on the above-mentioned periodic properties of the sn function, we know that $2nK(k)\,(n=0,1,2,...)$ are its zero points sn ; sn2 , i.e.,

[TABLE]

Moreover, based on the reciprocals and quotients of three basic functions ${\rm sn}(u,k),\,$cn(u,k),  ${\rm dn}(u,k)$, there are other nine Jacobi elliptic functions: ${\rm ns}(u,k)=1/{\rm sn}(u,k$, ${\rm nc}(u,k)=1/{\rm cn}(u,k)$, ${\rm nd}(u,k)=1/{\rm dn}(u,k)$, ${\rm sc}(u,k)={\rm sn}(u,k)/{\rm cn}(u,k)$, ${\rm sd}(u,k)={\rm sn}(u,k)/{\rm dn}(u,k)$, ${\rm cd}(u,k)={\rm cn}(u,k)/{\rm dn}(u,k)$, ${\rm cs}(u,k)={\rm cn}(u,k)/{\rm sn}(u,k)$, ${\rm ds}(u,k)={\rm dn}(u,k)/{\rm sn}(u,k)$, ${\rm dc}(u,k)={\rm dn}(u,k)/{\rm cn}(u,k)$ sn ; sn2 .

It follows from Eq. (18a) that we have $X(x)=\sqrt{\pi}/2[1+{\rm erf}(\sqrt{\omega}x)]\in(0,\sqrt{\pi})$. The choice of $X(x)$ ensures that $\Phi_{sn}[X(x),k]\to 0$ for all parameters $\nu,\,k$ and $G$ as $x\to-\infty$. Imposing a zero boundary condition for $x\to+\infty$ we have the constraint for the parameter $\nu$ in the solution (21)

[TABLE]

based on the property (25), that is to say, $\Phi_{sn}[X(x),k]\to 0$ for all parameters $k,\,G$, and $\nu$ given by Eq. (26) as $x\to+\infty$.

Therefore, for the given Gaussian nonlinearity and the combination of the harmonic and Gaussian potentials given by Eqs. (19a) and (19b), we obtain a family of exact stationary solutions of Eq. (16)

[TABLE]

where $k\in(0,1),\,\theta_{n}(x)=nK(k)[1+{\rm erf}(\sqrt{\omega}x)],\,n=1,2,3,...$, and $\omega=\gamma^{2}$ is a positive real number. It can be shown that $|\psi_{sn,n}(x,t)|\to 0$ as $x\to\pm\infty$, which means that the family of exact stationary solutions of Eq. (16) is bounded.

Therefore, we have four free parameters $G,\,n,\,\omega,$ and $k$ for the obtained solutions (27). In order to investigate the dynamics of exact analytical solutions (27), we choose the parameters as $G=1,\,n=1,\,\omega=0.1\,(\omega=0.2)$, and $k=0.4$, then the profiles of the Gaussian nonlinearity $g(x)$ given by Eq. (19a) (solid line) and the combination of the harmonic and Gaussian potential $V(x)$ given by Eq. (19b) (dashed line) are plotted in Figs. 1a and  1d, among these the nonlinearity can be generated by controlling the Feschbach resonances optically or magnetically using a Gaussian beam j5 ; fr . In Figs. 1b and  1e we show the intensity distribution of the exact solutions (27) as $t=0$. Moreover, we find that when $\omega$ becomes large ($0.1\to 0.2$), the Gaussian deep for the nonlinearity becomes high ($0.16\to 0.22$), the width of the soliton becomes narrow ($12.4\to 9.0$), and the the height of the soliton becomes low ($0.34\to 0.47$). Also, we study the stability of the exact solutions in response to perturbation by initial stochastic noise through direct numerical simulations (see Figs. 1c and 1f).

In the same way, for the case $n=2$, we display two-peak solutions given by Eq. (27) in Fig. 2 for other parameters $G=1,\,\omega=0.1\,(\omega=0.3),$ and $k=0.3$. The profiles of the Gaussian nonlinearity $g(x)$ given by Eq. (19a) (solid line) and the combination of the harmonic and Gaussian potential $V(x)$ given by Eq. (19b) (dashed line) are plotted in Figs. 2a and 2d. Figs. 2b and 2e plot the intensity distribution of the exact solutions (27) as $t=0$.

It is observed from Figs. 1c, 1f, 2c, and 2f that the exact solutions (27) with an initial amplitude stochastic noise of level $0.005$ for $n=1$ and $n=2$ are dynamically stable for the chosen parameters. Moreover, for the given parameters, we find that the amplitude noise regimes for the stable solutions are $(0,0.011)$ and $(0,0.03)$ for the solutions (27) with $n=1$ and $n=2$, respectively. In addition, we also analyze the effect of the phase noise on the solutions such that we find that, for the given parameters, the phase noise regimes for the stable solutions are $(0,0.05)$ and $(0,0.06)$ for the solutions (27) with $n=1$ and $n=2$, respectively.

For the given parameters $n,\,G$, and $k$, we analyze the stability of the solutions for the parameter $\omega$, which are related to the external potential, nonlinearity, and solutions. We find that the parameter regimes for the stable solutions are $(0,0.33)$ and $(0,0.18)$ for the solutions (27) with $n=1$ and $n=2$, respectively.

III.2 The attractive interaction

For the attractive interaction $G<0$, we have the periodic sd-wave solution of Eq. (5) with $X(x,t)\to X(x)$ yan06 ; yan

[TABLE]

where $k\in(0,1)$, $\nu$ is a constant, the function sd denotes the Jacobi elliptic function, i.e., ${\rm sd}[\nu X(x),k]={\rm sn}[\nu X(x),k]/{\rm dn}[\nu X(x),k]$ (c.f. Eq. (23)) sn , and the condition is required

[TABLE]

Since ${\rm dn}[\nu X(x),k]>0$ for all $X$ and $k$, thus we know that ${\rm sd}(\nu X(x),k)={\rm sn}(\nu X(x),k)/{\rm dn}(\nu X(x),k)$ has the same zero points as ${\rm sn}[\nu X(x),k]$, i.e., ${\rm sd}[2nK(k),k]=0\,(n=0,1,2,...)$. Similarly, based on the constraint for $\nu$ given by Eq. (26) for the requirement of the zero boundary condition for $x\to\pm\infty$, for the combination of harmonic and Gaussian potentials given by Eq. (19b), we have another family of exact solutions of Eq. (16)

[TABLE]

where $k\in(0,1),\,\theta_{n}(x)=nK(k)[1+{\rm erf}(\sqrt{\omega}x)],\,n=1,2,3,...$, and $\omega$ is a positive real number. It can be shown that $|\psi_{sd,n}(x,t)|\to 0$ as $x\to\pm\infty$, which means that the family of exact stationary solutions of Eq. (16) is bounded.

The required condition (29) leads to the following three cases for different modulus $k\in(0,1)$:

[TABLE]

which make the external potential $V(x)$ given by Eq. (19b) admit three different types of potentials, i.e., the harmonic-like potential ($\varepsilon>0$, see Figs. 3 and 4), the harmonic potential ($\varepsilon=0$), and the double-well potential ($\varepsilon<0$, see Figs. 5 and 6). Exact solutions of Eq. (16) with the harmonic potential case ($\varepsilon=0$) have been studied bb07 . Here we consider solutions of Eq. (16) with other two cases ($\varepsilon>0$ and $\varepsilon<0$).

III.2.1 The case $\varepsilon>0$

In order to investigate the dynamics of exact analytical solutions (III.2), we choose the parameters as $G=-1,\,n=1,\,\omega=0.1$, and $k=0.5$ satisfying Case I given by Eq. (34), i.e., $0<k<1/\sqrt{2}$. Thus the profiles of the Gaussian nonlinearity $g(x)$ given by Eq. (19a) (solid line) and the combination of harmonic and Gaussian potentials $V(x)$ given by Eq. (19b) (dashed line) are plotted in Fig. 3a, the intensity distribution of the exact solutions (III.2) for $t=0$ is displayed in Fig. 3b, and the numerical result for the evolution of the exact solutions (III.2) is shown in Fig. 3c.

Similarly, for the case $n=2$, we display the two-peak solutions given by Eq. (III.2) in Fig. 4 for these parameters $\omega=0.1,\,k=0.2$, and $G=-1$.

It is observed from Figs. 3c and 4c that the exact solutions (III.2) for $n=1$ and $n=2$ are dynamically stable for the chosen parameters. Moreover, for the given parameters, we find that the amplitude noise regimes for the stable solutions are $(0,0.02)$ and $(0,0.018)$ for the solutions (III.2) with $n=1$ and $n=2$, respectively. In addition, we also analyze the effect of the phase noise on the solutions such that we find that, for the given parameters, the phase noise regimes for the stable solutions are $(0,0.07)$ and $(0,0.05)$ for the solutions (III.2) with $n=1$ and $n=2$, respectively.

For the given parameters $n,\,G$, and $k$, we analyze the stability of the solutions for the parameter $\omega$, which are related to the external potential, nonlinearity, and solutions. We find that the parameter regimes for the stable solutions are $(0,0.2)$ and $(0,0.46)$ for the solutions (III.2) with $n=1$ and $n=2$, respectively.

III.2.2 The case $\varepsilon<0$

In this case, we choose the parameters as $G=-1,\,n=1,\,\omega=0.12$, and $k=0.8$ satisfying Case III given by Eq. (34), i.e., $1/\sqrt{2}<k<1$, to investigate the dynamics of exact analytical solutions (III.2). Thus the profiles of the Gaussian nonlinearity $g(x)$ given by Eq. (19a) (solid line) and the double-well potential $V(x)$ given by Eq. (19b) (dashed line) are plotted in Fig. 5a, the intensity distribution of the exact solutions (III.2) for $t=0$ is displayed in Fig. 5b, and the numerical result for the evolution of the exact solutions (III.2) is shown in Fig. 5c.

Similarly, for the case $n=2$, we display the two-peak solutions given by Eq. (III.2) in Fig. 6 for these parameters $\omega=0.07,\,k=0.75$, and $G=-1$ also satisfying Case III given by Eq. (34), i.e., $1/\sqrt{2}<k<1$.

It is observed from Figs. 5c and 6c that the exact solutions (III.2) for $n=1$ and $n=2$ are dynamically stable for the chosen parameters. Moreover, for the given parameters, we find that the amplitude noise regimes for the stable solutions are $(0,0.04)$ and $(0,0.05)$ for the solutions (III.2) with $n=1$ and $n=2$, respectively. In addition, we also analyze the effect of the phase noise on the solutions such that we find that, for the given parameters, the phase noise regimes for the stable solutions are $(0,0.05)$ and $(0,0.04)$ for the solutions (III.2) with $n=1$ and $n=2$, respectively.

For the given parameters $n,\,G$, and $k$, we analyze the stability of the solutions for the parameter $\omega$, which are related to the external potential, nonlinearity, and solutions. We find that the parameter regimes for the stable solutions are $(0,0.13)$ and $(0,0.07)$ for the solutions (III.2) with $n=1$ and $n=2$, respectively.

IV Time-depended potential and nonlinearity

So far we have considered stationary solutions of the GP equation with the potential and nonlinearity depending only on the spatial coordinate. Now we focus on the case that the potential $V(x,t)$ and nonlinearity $g(x,t)$ depend on both time and spatial coordinates, and the gain or loss term $\Gamma(t)$ is a function of time (see Eq. (2)).

IV.1 The repulsive interaction

For the determined potential $V(x,t)$ and nonlinearity $g(x,t)$ and other conditions given by Eqs. (14) and (15), we have a family of exact solutions of Eq. (2) with the repulsive interaction $G>0$ by using the solution of Eq. (5) given by Eq. (21) and the determined similarity transformation (4)

[TABLE]

where the function sn denotes the Jacobi elliptic function with $k\in(0,1)$ being its modulus,  $K(k)=\int^{\pi/2}_{0}[1-k^{2}\sin^{2}(s)]^{-1/2}ds$, $\theta_{n}(x,t)=nK(k)[1+{\rm erf}(\gamma(t)x)],\,n=1,2,3,...$, and the phase $\varphi(x,t)$ is given by Eq. (15d).

To make sure that the nonlinearity $g(x,t)$ and the coefficients of the external potential $V(x,t)$ in space are bounded for realistic cases, we choose the inverse of the width of the localized solution $\gamma(t)$ and the gain or loss distribution $\Gamma(t)$ as the localized periodic functions in the form

[TABLE]

For the case $n=1$, we display the one-peak solutions given by Eq. (IV.1). In Fig. 7, we plot the potential $V(x,t)$ given by Eq. (14), the nonlinearity $g(x,t)$ given by Eq. (15c), the intensity distribution of exact solutions (IV.1), and the numerical result for the evolution of exact solutions (IV.1) for the chosen parameters.

Similarly, for the case $n=2$, we display the two-peak solutions given by Eq. (IV.1) in Fig. 8 for other chosen parameters $\gamma(t)=6.4[0.1+0.01\cos(t)],\,\Gamma(t)=-0.3\sin(3t),\,k=0.2$, and $G=1$.

It is observed that the amplitude of the potential $V(x,t)$ varies quasi-periodically with respect to time, so does the nonlinearity $g(x,t)$. It is also observed from Figs. 7d and 8d that exact solutions (IV.1) for $n=1$ and $n=2$ are dynamically stable for the chosen parameters. Moreover, for the given parameters, we find that the amplitude noise regimes for the stable solutions are $(0,0.011)$ and $(0,0.015)$ for the solutions (IV.1) with $n=1$ and $n=2$, respectively. In addition, we also analyze the effect of the phase noise on the solutions such that we find that, for the given parameters, the phase noise regimes for the stable solutions are $(0,0.05)$ and $(0,0.04)$ for the solutions (IV.1) with $n=1$ and $n=2$, respectively.

For the given parameters $n,\,G,\,k$, and $\gamma(t)$ given by Eq. (39), we analyze the stability of the solutions for the varying gain or loss term

[TABLE]

which are related to the external potential, nonlinearity, and solutions. We find that the parameter ($c)$ regimes for the stable solutions are $(-0.3,0)$ and $(-0.3,0)$ for the solutions (IV.1) with $n=1$ and $n=2$, respectively.

IV.2 The attractive interaction

We have a family of exact solutions of Eq. (2) with the attractive interaction $G<0$ by using the solution (28) of Eq. (5) and the similarity transformation (4)

[TABLE]

where the sd function ${\rm sd}(\cdot,k)\equiv{\rm sn}(\cdot,k)/{\rm dn}(\cdot,k)$ denotes the Jacobi elliptic function with $k\in(0,1)$ being its modulus, $K(k)=\int^{\pi/2}_{0}[1-k^{2}\sin^{2}s]^{-1/2}ds$, $\theta_{n}(x,t)=nK(k)[1+{\rm erf}(\gamma(t)x)],\,n=1,2,3,...$, and the phase $\varphi(x,t)$ is defined by Eq. (15d).

In Fig. 9, we plot the potential $V(x,t)$ given by Eq. (14), the attractive interaction $g(x,t)$ given by Eq. (15c), the intensity distribution of the exact solutions (IV.2) and the numerical result for the evolution of the exact solutions (IV.2) for the chosen parameters

[TABLE]

Similarly, for the case $n=2$, we show the two-peak solutions given by Eq. (IV.2) in Fig. 10 for another parameters $\gamma(t)=6.3[0.1+0.01\cos(t)],\,\Gamma(t)=-0.3\sin(3t),\,k=0.1$, and $G=-1$.

It is observed that the amplitude of the potential $V(x,t)$ varies quasi-periodically with respect to time, so does the nonlinearity $g(x,t)$. It is also observed from Figs. 9d and 10d that the exact solutions (IV.2) for $n=1$ and $n=2$ are dynamically stable for the chosen parameters. Moreover, for the given parameters, we find that the amplitude noise regimes for the stable solutions are $(0,0.015)$ and $(0,0.01)$ for the solutions (IV.2) with $n=1$ and $n=2$, respectively. In addition, we also analyze the effect of the phase noise on the solutions such that we find that, for the given parameters, the phase noise regimes for the stable solutions are $(0,0.05)$ and $(0,0.04)$ for the solutions (IV.2) with $n=1$ and $n=2$, respectively.

For the given parameters $n,\,G,\,k$, and $\gamma(t)$ given by Eq. (45), we analyze the stability of the solutions for the gain or loss term given by Eq. (40), which are related to the external potential, nonlinearity, and solutions. We find that the parameter ($c)$ regimes for the stable solutions are $(-0.3,0)$ and $(-0.3,0)$ for the solutions (IV.2) with $n=1$ and $n=2$, respectively.

V Conclusions

In conclusion, we have studied matter-wave solutions in quasi-one-dimensional Bose-Einstein condensates in the combination of the harmonic and Gaussian potentials. As a consequence, we have found several families of exact solutions of the GP equation with the combination of the harmonic and Gaussian potentials by using the similarity transformation. Moreover, the stability of the obtained exact solutions is investigated by using numerical simulations such that some stable solutions are found. In particular, we analyze the effects of both the amplitude noise and the phase noise. We also give the parameter regimes for the stable solutions. These results may raise the possibility of relative experiments and potential applications.

The used method can be also extended to investigate exact solutions of some other physical models in BECs such as the cubic-quintic, the two-component, the spinor-1, and/or higher-dimensional GP equations with the combination of the harmonic and Gaussian potentials. In fact, our idea can also be applied to the related NLS equation with varying coefficients and its extensions in nonlinear optics b1 ; b4 ; j8 .

Acknowledgements.

This work was supported by the NSFC under Grant Nos. 11071242 and 61178091.