Quantum breathers associated with modulational instability in 1D ultracold boson in optical lattices involving next-nearest neighbor interactions

Optik 164 (2018) 575–589

Contents lists available at ScienceDirect

Optik journal homepage: www.elsevier.de/ijleo

Short note

Quantum breathers associated with modulational instability in 1D ultracold boson in optical lattices involving next-nearest neighbor interactions Z.I. Djoufack a,d,∗ , E. Tala-Tebue a , F. Fotsa-Ngaffo b,c , A.B. Djimeli Tsajio a , F. Kapche-Tagne a a University Institute of Technology Fotso Victor, Department of Telecommunication and Network Engineering, University of Dschang, P.O. Box 134 Bandjoun, Cameroon b Intitute of Wood Technologies, University of Yaounde I, P.O. Box 306 Mbalmayo, Cameroon c Department of Physics, Faculty of Science, University of Buea, Cameroon d African Institute for Mathematical Sciences, 6 Melrose, Muizenberg, Cape Town 7945, South Africa

a r t i c l e

i n f o

Article history: Received 17 January 2018 Received in revised form 17 March 2018 Accepted 18 March 2018 Keywords: Quantum breathers Ultracold boson in optical lattices Modulational instability Next-nearest neighbor interactions

a b s t r a c t The dynamics and modulation instability of an ultracold gas of bosonic atoms in an optical lattice can be portrayed by a Bose–Hubbard model and the system parameters are mastered by laser light. Based on the time dependent Hartree approximation combined with the semi-discrete multiple-scale method, the equation of motion for single-boson wave function is found analytically, the existence conditions of appearance of bright stationary localized solitons solutions of this quantum Bose–Hubbard model are discussed. We find that the introduction of the next nearest neighbor interactions (NNNI) may change the stability property of the plane waves and may predict the formation of modulational instability in the wave number k = kmax and k = keBZ in the system. With the help of stationary localized single-boson wave functions obtained, the quantized energy level and the quantum breather state are determined. The performance of the analytical results are checked by numerical calculations. Furthermore, we have shown that the presence of the NNNI affect significatively the shape of the region of modulational instability and it is responsible of the appearance of new region of modulational instability that occurs for the k = kmax carrier wave. The formation conditions of the modulational instability region predicted by the analytical analysis of stationary localized solutions in the wave number k = kmax and k = keBZ are in good agreement with the forecast respectively. © 2018 Elsevier GmbH. All rights reserved.

1. Introduction During the last few years, ultracold atomic gases have attracted tremendous attention. This is due to the experimental achievement of quantum degeneracy for bosonic [1] and fermionic [2] gases. Atoms in optical lattice are considered as

∗ Corresponding author at: Fotso Victor University Institute of Technology, Department of Telecommunication and Network Engineering, University of Dschang, P.O. Box 134 Bandjoun, Cameroon. E-mail address: [email protected] (Z.I. Djoufack). https://doi.org/10.1016/j.ijleo.2018.03.059 0030-4026/© 2018 Elsevier GmbH. All rights reserved.

576

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

one of more promising systems for which special emphasis has opened a new window to simulate many-body quantum systems. The Bose–Hubbard (bosonic gase) and Hubbard (fermion gase) are two most physical models introduced that can be simulated in optical lattice. Recently, the experimental developments in the manipulation of ultracold atoms have opened a new fascinating research field, for instance, the analysis of strongly-correlated atomic gases. In this way, experimental results on Bose–Einstein condensates in optical lattices [3] and experiments on ultracold gases in optical lattices, such as, the observation of the Mott-insulator to superfluid transition have already been reported [4]. The Bose–Hubbard model which can be used to describe the dynamics of an ultracold gas of bosonic atoms [7], has attracted large interest and can be considered today as a promising candidate for their practical implications in quantum computing because quantum computation is very sensitive to even smallest disturbances from environment [6,5]. In the classical limit, it is equivalent to a large number of bosons and can be well approximated by the Discrete Nonlinear Schrodinger (DNS) which is one of the most studied nonlinear lattice systems [8]. On experimental development in optical lattices, ultracold atoms offers the unprecedented potential to study the nonlinear properties of many-body interacting [9]. In classical nonlinear lattices, discrete breathers which are localized nonlinear excitation, generic time periodic and spatially localized solutions of the nonlinear systems [10,11], can play important role in both energy storage and transport [12] and have been observed experimentally and theoretically in many physical systems [13–18]. MacKay and Aubry were the first to establish the existence of breathers excitations in 1D nonlinear lattices [19]. However, their quantum counterparts called quantum breather, according to Fleurov [20], can be defined as superpositions of nearly degenerate many quanta bound states characterized with long times to tunnel from one lattice site to another. Quantum breather is nowadays, one of the most interesting domain that attracts the attention of searchers in many physical systems. In the Bose–Hubbard model, few studies have been devoted to prove the existence of quantum breathers. For instance, Tang [21] investigated the existence of quantum two-breathers formed by ultracold bosonic atoms in optical lattices for a large number of bosons. More recently, we probed the quantum signature of breathers in 1D ultracold bosons in optical lattices involving next-nearest neighbor interactions for a few number of bosons [22] and the importance of nearest and next-nearest-neighbor off-site interactions has also been emphasized experimentally in the extended Bose–Hubbard model [23]. Quantum breathers states which are the solitons bright solution are directly linked to the modulational instability of plane waves as demonstrated by Lai and Siever [24] that, the modulational instability is the basic mechanism for the formation of nonlinear localized structures in spin systems. The modulational instability is a basic process in the theory of nonlinear wave [25–29]. Modulational instability is characterized by an exponential growth of amplitude perturbations occurring under the interplay between nonlinearity and dispersion during propagation. The best way to predict the formation of intrinsic localized modes in nonlinear lattice is Modulational instability [30–32]. The generation of localized states in nonlinear lattices by the modulational instability was established firstly by Kivshar and Peyrard [33]. Kivshar showed that, the formation of bright intrinsic localized mode is instantly linked to the modulational instability of carrier waves [30]. However, in the case of the Bose–Hubbard model, few studies have been devoted just to analyze the non-equilibrium quantum dynamics of a Bose-Einstein condensate and to study the density modulations associated with the dynamical instability [34,35]. From the foregoing, it is clear that the study of the modulation instability in 1D ultracold boson in optical lattices is not completely achieved and needs more attention especially when the NNNI are included in the system. We will in this paper study the dynamics of ultracold bosonic atoms loaded in an optical lattice involving the NNNI and will look for the existence condition of quantum breathers associated with modulational instability in the system. The organization of this paper is as follows: In the next section, we rewrite the Hamiltonian of an ultracold gas of bosonic atoms in an optical lattice including the NNNI describes by the Bose–Hubbard model. In section III, we use the time dependent Hartree approximation to derive the equation of motion describing the dynamics of ultracold gas of bosonic atoms in an optical lattice. From the dispersion curves, the existence regions of bright and dark solitons stationary localized solutions are discussed in the half of the Brillouin zone. In section VI, a numerical simulation is used to check our analytical prediction. In section V, the modulational instability conditions are studied through the linear stability analysis and numerical verifications of analytical forecast are performed. The last section is the conclusion.

2. The model Hamiltonian involving the NNNI Since Jaksch et al. [7] have shown that the dynamics of the bosonic atoms on the optical lattices realize a Bose-Hubbard model, describing the hopping of bosonic atoms between the lowest vibrational states of the optical lattice sites for which the system can be controlled by the laser parameters and configurations [36]. It is important to use the Bose–Hubbard model to describe the dynamics of an ultracold gas of bosonic atoms because this model provides good description of utracold bosonic atoms in the limit of sufficiently deep optical lattice. An optical lattice is formed by pairs of counter-propagating laser beams, which creates effective potential that traps ultracold atoms. An ultracold atom loaded into a three-dimensional (3D) optical

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

577

lattice provide a realization of a quantum lattice gas [37,38]. Taking into account the next-nearest neighbor hopping, the Bose–Hubbard Hamiltonian [39,22] adopted in this work can be described as

H=

−J1

N  

†

†

aj aj+1 + aj aj+1



j=1

− J2

N  

†

j=1 N

+ V2



†

aj aj+2 + aj aj+2



+ V1

N 

nj nj+1

(1)

j=1

U nj (nj − 1), 2 N

nj nj+2 +

j=1

j=1

†

where j is atom site in the lattice, aj and aj are the bosonic creation and annihilation operators satisfying the commutation † [ai , aj ]

† † ıij , [ai , aj ]

= = [ai , aj ] = 0. J1 stands for the nearest-neighbor hopping strength, J2 denotes the hopping relations strength to the next-nearest-neighbor sites, N represents the number of sites, U is known as the interaction energy of two bosons at the same site. If U > 0, the system is repulsive and attractive in the contrary case, V1 and V2 are the interaction † energy of two bosons at two nearest-neighboring sites and at two next-nearest-neighboring sites. nj = aj aj represents the boson number.

3. Stationary localized solutions The Schrödinger picture is used in this work to study the quantum nonlinear ultracold bosons in deep optical lattices. As we can see in (1), the Hamiltonian operator of the system is not time dependent. But the temporal evolution of the lattice system state vector |(t) which is time dependent follows the following equation

i

d |(t) = H|(t) dt

(2)

Due to the fact that the chain is tanslationally invariant, the Hamiltonian of this quantum system commutes with the number † operator N = aj aj , whose eigenvalue is n. Therefore, the total boson number is conserved. In this case we can use as a basic representation of the lattice, the Fock state representation to expand any quantum state of the system. We can write a general n-boson system state vector [40] in the Fock space as

1   † † |n  = √ ··· ˇn (j1 , · · ·, jn , t)aj · · ·aj |0, n 1 n! N

N

j1 =1

jn =1

(3)

where |0 ≡ |01 |02 · · · |0N is the vacuum state and ˇn are fn time dependent coefficients of corresponding number states. ˇ(j1 , j2 , · · ·, jn ) is the n-boson wave function, which satisfies the following normalized condition

N N N   

···

j1 =1 j2 =1

jn =1

ˇn |(j1 , j2 , · · ·, jn , t)|2 = 1.

(4)

578

Z.I. Djoufack et al. / Optik 164 (2018) 575–589 †

Using the Bose commutation relations [ai , aj ] = ıij and inserting (3) into (2), we get the Schrodinger equation for the n-boson wave function given as follows

i

d ˇn (j1 , j2 , · · ·, jn , t) dt

= −J1

n  

ˇn (j1 , j2 , · · ·, jk−1 , jk + 1, jk+1 , · · ·, jn , t)



k=1

+ˇn (j1 , j2 , · · ·, jk−1 , jk − 1, jk+1 , · · ·, jn , t)

− J2

n  

ˇn (j1 , j2 , · · ·, jk−1 , jk + 2, jk+1 , · · ·, jn , t)



k=1

+ˇn (j1 , j2 , · · ·, jk−1 , jk − 2, jk+1 , · · ·, jn , t)

+ V1

n n  

(5)

ıjl ,jk+1 ˇn (j1 , j2 , · · ·, jl , · · ·, jk , · · ·, jn , t)

k=1 l = / k n n

+ V2



ıjl ,jk+2 ˇn (j1 , j2 , · · ·, jl , · · ·, jk , · · ·, jn , t)

k=1 l = / k n n

+

U  ıjl ,jk ˇn (j1 , j2 , · · ·, jl , · · ·, jk , · · ·, jn , t). 2 k=1 l = / k

We note here that the interaction between pairs of bosons is Kronecker delta-function. (5) is very difficult to solve, to reach our objective, we use the time-dependent Hartree approximation. This approximation is well known in quantum field theory [41] and it is applicable when the number of quanta is large. In the Hartree approximation, n-boson wave function ˇn (j1 , j2 , · · ·, jn , t) can be written as a product which takes the following form [42,43]

(H ) ˇn (j1 , j2 , · · ·, jn , t) =

n 

n,jk (t),

(6)

k=1

where n,jk (t) is the single-boson wave function, jk = 1, 2, . . ., f and k = 1, 2, . . ., n represent the boson. These single-boson wave functions are independent of k and they can be written as n,j (t), where j = 1, 2, . . ., f. With the aid of the Hartree approximation wave function in (6), the n-boson state vector in (3) can be reduced to

⎛

|n (t)(H)

⎞n

 1 † = √ ⎝ n,j (t)aj ⎠ |0, n! N

(7)

j=1

and from (3), the normalization condition is also reduced to

 j

|n,j |2 = 1.

(8)

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

579

According to [44], the equation of motion for n,j (t) can be determined by extremizing the following functional



∞

S=

dt

N N  

···

−∞

j1 =1 n

+ J1

ˇn∗

d i

jn =1



dt

ˇn (j1 , j2 , · · ·, jn , t)

ˇn (j1 , j2 , · · ·, jk−1 , jk + 1, jk+1 , · · ·, jn , t)



k=1

+ˇn (j1 , j2 , · · ·, jk−1 , jk − 1, jk+1 , · · ·, jn , t)

+ J2

n  

ˇn (j1 , j2 , · · ·, jk−1 , jk + 2, jk+1 , · · ·, jn , t)



k=1

(9)

+ˇn (j1 , j2 , · · ·, jk−1 , jk − 2, jk+1 , · · ·, jn , t)

− V1

n n  

ıjl ,jk+1 ˇn (j1 , j2 , · · ·, jl , · · ·, jk , · · ·, jn , t)

k=1 l = / k n n

− V2



ıjl ,jk+2 ˇn (j1 , j2 , · · ·, jl , · · ·, jk , · · ·, jn , t)

⎫ ⎬

k=1 l = / k

U  ıjl ,jk ˇn (j1 , j2 , · · ·, jl , · · ·, jk , · · ·, jn , t) . − 2 ⎭ n

n

k=1 l = / k

Substituting (6) into (9), one can get the following expression

S H

=

∞

n

dt −∞

N 

 ∗n,j

i

j1 =1

dn,j dt

+ J1 (∗n,j n,j−1

+ ∗n,j n,j−1 ) + J2 (∗n,j n,j−2 + ∗n,j n,j−2 ) − V1 (n − 1)∗n,j ∗n,j−1 n,j n,j−1

(10)

− V2 (n − 1)∗n,j ∗n,j−2 n,j n,j−2 −

U (n − 1)|n,j |4 2



.

Requiring ıS H /ı∗n,j = 0 for the optional Hartree solution, we get for the single-boson wave function, the following equation by neglecting the subscript n and adopting the transformation n,j → j : i

dj dt

+ J1 (j+1 + j−1 ) + J2 (j+2 + j−2 )

− (n − 1)V1 (|j+1 |2 + |j−1 |2 )j

(11)

− (n − 1)V2 (|j+2 |2 + |j−2 |2 )j − (n − 1)U|j |2 j = 0. (11) is a complicated discrete nonlinear equation to solve analytically and get exact solution. To get round, we adopt the semidiscrete approximation combined with multiple-scale method [45], we can solve approximate solutions of Eq. (11). Taking into account the normalization condition, when N is large enough, |j | is a small quantity. In this case, we must first make a scale transformation as j = εj , where ε is small but finite parameter. With the aid of this transformation, (11) becomes i

dj dt

+ J1 (j+1 + j−1 ) + J2 (j+2 + j−2 )

− (n − 1)ε2 V1 (|j+1 |2 + |j−1 |2 )j − (n − 1)ε2 V2 (|j+2 |2 + |j−2 |2 )j − (n − 1)ε2 U|j |2 j = 0.

(12)

580

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

Next, one can set a solution of (12) in the form j = j (1 , 2 , j ) exp(iˇj ),

(13)

= ε2 t,

where  1 = εt,  2  j = εja are multiple-scale variables, a is the lattice constant, ˇj = kja − ωt represents the phase of the carrier wave, k is the wave number and ω is the frequency of the carrier wave. d/dt and d/dx are derivatives operators which can be expended as follows d ∂ ∂ ∂ = + ε2 +ε dt ∂t ∂1 ∂2

(14)

d ∂ =ε . dx ∂j Substituting (14) and (13) into (12), we have [ω + 2J1 cos(ka) + 2J2 cos(2ka)]



× j exp(iˇj ) + iε +2J2 a sin(2ka))]

∂j + 2(J1 sin(ka) ∂1

∂j exp(iˇj ) + ε2 ∂j



i

∂j ∂2

(15)

2

∂ j + 2a2 [J1 cos(ka) + 4J2 cos(2ka)] exp(iˇj ) ∂j2



−(n − 1)(2V1 + 2V2 + U)|j |2 j exp(iˇj ) + o(ε2 ) = 0. Comparing the powers of ε in (15), then we get the following equations: ω + 2J1 cos(ka) + 2J2 cos(2ka) = 0,

(16)

∂j ∂j + 2[J1 sin(ka) + 2J2 a sin(2ka)] = 0, ∂1 ∂j

(17)

2

i

∂j ∂ j + 2a2 [J1 cos(ka) + 4J2 cos(2ka)] ∂2 ∂j2 − (n − 1)(2V1 + 2V2 + U)|j |2 j

(18)

= 0.

(16) is the dispersion relation of linear waves in which we can derive the group velocity Vg =

∂ω = 2a[J1 sin(ka) + 2J2 sin(2ka)]. ∂k

Let us set the value of J1 = 1 and on varying the value of J2 from 0.0 to 0.8 as used in Refs. [39,22] to study the effect of NNN hopping on ultracold gas of bosonic atoms in optical lattices. It is clearly seen in Fig. (1) that the presence of NNNI significatively changes the curvature of the linear dispersion curve in the half Brillouin zone. As we can see in Fig. (1), we will notice that, when the NNNI coefficient is small i.e. J2 /J1 < 1/4, it appears on the spectrum curves only a maximum and a minimum at the edge (wave number k = keBZ = /a) and at the center (wave number k = kcBZ = 0) of Brillouin zone respectively. However when the NNNI coefficient is strong enough, for i.e. J2 /J1 > 1/4, a second maximum appears on the spectrum at k = kmax = /a − arccos(J1 /4J2 )/a especially at the middle of the Brillouin zone. The inflection point (k = k2 ) which occurs in the dispersion spectrum for the critical value J2 /J1 = 1/4, represents the zero-dispersion point of the dispersion spectrum. The corresponding frequency obtained at kmax is given as follows ωmax (kmax ) = (J12 + 8J22 )/4J2 .

(19)

An inflection point located between kmax and keBZ appears on the dispersion spectrum at k = k2 , where



 1 k2 = − arccos a a

J1 +



J12 + 128J22

16J2



.

Eq. (19) presented above, displays the relationship between the temporal and spatial derivatives of wave deduced from (19), which can take the following form (1 , 2 , j ) = j (2 , z1 ).

(20) .

is a traveling (21)

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

581

Fig. 1. Dispersion curves for different values of J2 /J1 representing the coupling strength between nearest and NNNI.

Fig. 2. Dispersion curves of

d2 ω dk2

representing the appearance of the region for different values of J2 /J1 .

We introduce the new scale as follows z1 =  j − Vg  1 . With the help of new scale variable, Eq.(19) can be rewritten as 2

i

∂j ∂ j +P + Q |j |2 j = 0, ∂2 ∂z12

where P =

d2 ω dk2

(22)

= 2a2 [J1 cos(ka) + 4J2 cos(2ka)] and Q =− (n − 1)(2V1 + 2V2 + U). We can modify the form of Eq. (22) by creating

new variables that are z1 = εz, j = A/ε and  2 = ε2 t, then Eq. (22) becomes 2

i

∂A ∂ A + Q |A|2 A = 0. +P ∂t ∂z 2

(23)

Eq. (23) stands for the NLS equation. The solution of NLS (23) depends on the sign of the product QP when it is positive or negative. In this work, we have to make sure that the single-boson wave function is bounded thus we consider only the case 2 d2 ω where QP > 0. Since Q < 0, the regions of k are just determined by the sign of P = d ω 2 . The dispersion curves of 2 displaying dk

dk

regions are shown in Fig. 2 in the half of Brillouin zone for J2 /J1 = 0.24 and J2 /J1 = 0.27, respectively. In presence of NNNI, the Table 1 exhibits different regions of possible nonlinear localized excitations when J2 /J1 > 1/4. However, we realize that Table 1 Different regions of possible nonlinear localized excitations in the ultracold boson and in presence of NNNI when J2 /J1 = 0.27. Regions of k

d2 ω dk2

Stability

Type of soliton

(0, k1 ) (k1 , k2 ) (k2 , /a)

>0 <0 >0

stable unstable stable

dark bright dark

582

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

when J2 /J1 < 1/4, the stable region (k2 , /a) disappears and k2 takes the value /a. The half of the Brillouin zone in this case is divided into two regions, where stable and unstable regions occur in (0, k1 ) and (k1 , /a), respectively. Therefore, we can conclude here by saying that, the introduction of NNNI may change the stability property of the plane waves in the ultracold boson system. In the case of negative dispersion region (k1 , k2 ), the solution of Eq. (23) corresponding to the bright soliton is given by [46]

  A(t, z) =

2 sech 2 × exp



 Q (z − z0 ) 2P

(24)



i(2 2 Qt

+ ˇ0 ) ,

where stands for height or width, z0 is the initial position and ˇ0 is the initial phase. From Eq. (24), we define the singleparticle wave function on the jth site as follows: n,j (t)

= A(t, z) exp[i(kja − ωt)]

 

= 2 sech 2



Q (j − j0 )a − 2 2P



 Q Vg t 2P

(25)



× exp[ikja − i(ω − 2 2 Q )t + iˇ0 . If we chose a particular wave vector kmax = /a − arccos(J1 /4J2 )/a belonging to the region (k1 , k2 ), ωmax (kmax ) = (J12 + 8J22 )/4J2 . PQ = −a2 (n − 1)(2V1 + 2V2 + U)( Smax = Q/(2P) = −J2

max =

J12 − 16J22 2J2

) (26)

(n − 1)(2V1 + 2V2 + U) a2 (J12 − 16J22 )

2 ωmax − 2 max Q.

In this case where J2 /J1 > 1/4, the sign of the product PQ > 0 and the solution (25) can be rewritten as follows: n,j (t)

= 2 max sech[2 max



Smax (j − j0 )a]

(27)

× exp(ikmax ja − i max t + iˇ0 ). Eq. (27) is the bright soliton solution whose eigenfrequency max > ωmax lies above the top of the linear wave spectrum. From (4), n,j (t) should satisfy the normalization condition. The discrete sum can be changed into continue integral, because the envelope part of n,j (t) is looked as continuum variable in the approximate treatment, then we have



|n,j (t)|

2

=

j

=

where  = 2 max max



2

max

a

4 max

a

2

sech ()d

Smax



+∞

Smax

−∞

(28)

= 1,



Smax (z − z0 ). From Eq. (28), the value of max is given by:



a 1 = Smax = 4 4

−J2

(n − 1)(2V1 + 2V2 + U) J12 − 16J22

(29)

.

Substituting (27) into (7) and using (29), we can construct the Hartree product eigenstates as

|n (t)(H) =

1 exp[−i( max t + ˇ0 )n] √ n!

⎛



×⎝

j=1



  a

n Smax

2



⎞n

a2 Smax † sech (j − j0 ) aj ⎠ |0. 2

(30)

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

583

From (30), the average number of bosons on the jth freedom are the stationary states which can be evaluated by the following formula †

nj (H)

=(H) n (t)|aj aj |n (t)(H) a2 nSmax = sech2 4





(31)

a2 Smax (j − j0 ) . 2

The Hartree n-boson states given by (30) can be considered as quantum breathers states and its quantized energy (30) found is



= n max = n

Emax,n

J12 + 8J22



4J2

− J2

n(n − 1) (2V1 + 2V2 + U)2

(32)

2

8(J12 − 16J22 )

.

For k = keBZ = /a, in this particular case where J2 /J1 < 1/4, PQ =− a2 (n − 1)(2V1 + 2V2 + U)(−J1 + 4J2 ) > 0 and we have ωeBZ = ω(k = /a) = 2(J1 − J2 ) SeBZ = Q/(2P) =

(33)

−(n − 1)(2V1 + 2V2 + U) 2a2 (−J1 + 4J2 )

2 Q.

eBZ = ωeBZ − 2 eBZ

The stationary localized solution in this case is given by = 2 eBZ sech[2 eBZ

n,j (t)



SeBZ (j − j0 )a]

(34)

× exp(ikeBZ ja − i eBZ t + iˇ0 ), where eBZ =

1 4



−(n − 1)(2V1 + 2V2 + U) . 2(−J1 + 4J2 )

(35)

(34) stands for the quantum breathers solution whose eigenfrequency eBZ > ωeBZ also lies above the top of the linear wave spectrum. Replacing (34) into (7) and using (35), we construct the Hartree product eigenstates as (H)

|n (t)

1 = exp[−in( eBZ t + ˇ0 )] √ n!

⎛



×⎝



SeBZ (−1) sech 2 j

a2

  a

n

SeBZ

2

⎞n

 (j − j0 )

† aj

(36)

⎠ |0.

j=1

From (36), the average number of bosons on the jth freedom are the stationary states can be given by: (H)

nj 

SeBZ a2 n sech2 = 4





SeBZ a2 (j − j0 ) . 2

(37)

The corresponding Hartree energy of quantum breather states (36) in this case becomes EeBZ,n = 2n(J1 − J2 ) − n

(n − 1)2 (2V1 + 2V2 + U)2 . 8(−J1 + 4J2 )

(38)

Eq. (38) indicates that, the energy of this quantum breathers states is quantized. 4. Numerical simulations In order to check numerically our analytical results, we use the split-step (Fourier) to test the behavior of quantum breathers solutions obtained in Section 3. All analytical solutions found in Section 3 are numerically checked and we noted that these localized quantum solutions are spatially localized excitation with a long lifetime which is full analogous to their classical counterparts [47]. Let us have a look at the region k1 < k < k2 , according to our analytical outcomes mentioned in section III, this region displays a bright soliton excitation. For this case, taking Eq. (25) as an initial solution of Eq. (23), we show in Fig. 3(b) the time evolution of the average number of bosons on the jth freedom with k = kmax = /a − arccos(J1 /4J2 )/a. The shape of the

584

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

Fig. 3. Time evolution of the average number of bosons on the jth freedom, where a = 1, n = 15, V1 = 5 ×10−3 , V2 = 10−3 , U = 9 ×10−3 , f = 195 and j0 = 0 .: (a) J1 /J2 = 0.2 and k = keBZ = /a; (b) J1 /J2 = 0.27 and k = kmax = /a − arccos(J1 /4J2 )/a.

time evolution of the average number of bosons confirms that the solution Eq. (30) is stable and localized. It is clearly seen in Fig. 3 (b) that, many bosons are localized in sites of the vicinity at the central position for j = j0 . We also realize that the site j = j0 depends on the initial exciting condition, this indicates that the system is translational symmetry and that may be excited at any lattice site in the system. Next, for the case J2 /J1 < 1/4, the unstable region (k1 , /a) attracts our attention. According to our analytical result, we also have a bright soliton solution. This solution is tested by our numerical simulation for k = keBZ . As we can see in Fig. 3(a), the solution Eq. (36) is also stable and localized. If we compare the amplitude seize obtained in this last case to the one presented in Fig. 3(b) when the NNNI are involved in the system, we note that the amplitude seize is high. Thus the intensity of localization of quantum breathers can be reduced by including the NNNI in the ultracold boson in optical lattices. From these latter figures, we conclude that the analytical results obtained in section III have been tested by numerical simulation and good agreement is found. Lai and Siever [24] demonstrated that the modulational instability is the basic mechanism for the formation of nonlinear localized structures in spin systems. This implies that the formation of nonlinear localized bright solutions in ultracold bosons is immediately linked to the modulational instability of carrier waves. Then, we expect to have the modulational instability at the wave number k = kmax and at k = keBZ . In section V, these expectations we will be probed.

5. Modulational instability In this section, we adopt the standard linear stability analysis to study the modulational instability of (11). We consider the propagation of wave signal along the ultracold bosons in optical lattice and we find that (11) admits a constant amplitude solution in the form j = 0 exp[i(kaj − ωt)],

(39)

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

585

where a is the lattice constant, 0 is the initial amplitude, k and ω are respectively the wave number and the angular frequency of the plane wave. Substituting (39) into, (11), one can get the nonlinear dispersion relation given as ω

= −2J1 cos(ka) − 2J2 cos(2ka)

(40)

+ (n − 1)(2V1 + 2V2 + U)|0 |2 .

In order to study the linear stability of the initial plane waves in the system, we introduce a small perturbation ıj (t) and set a solution as j (t) = [0 + ıj (t)] exp[i(kaj − ωt)].

(41)

Inserting (41) into (11) and keeping in mind only the linear terms, we obtain the linear differential equation, which is written as i

dıj dt



= J1 2 cos(ka)ıj + ıj−1 exp[−ika]





+ıj+1 exp[ika] + J2 2 cos(2ka)ıj +ıj−2 exp[−2ika] + ıj+2 exp[2ika]

 (42)

− (n − 1)V1 (ı∗j+1 + ı∗j+1 + ı∗j−1 + ı∗j−1 )|0 |2 − (n − 1)V2 (ı∗j+2 + ı∗j+2 + ı∗j−2 + ı∗j−2 )|0 |2 − (n − 1)U(ı∗j + ı∗j )|0 |2 . The general solution of (42) has the form ıj (t) = f+ cos[Qj − t] + if − sin[Qj − t],

(43)

where k and are the wave number and the frequency of the perturbation, respectively. f+ and f− stand for constant perturbations amplitudes. Inserting (43) into (42), we have the linear equations for f+ and f−



M11 −

M12

M21

M22 −



f+ f−



  =

0 0

,

(44)

with M11 = M22 = 2J1 sin(ka) sin Q + 2J2 sin(2ka) sin(2Q ) M12 = 2J1 cos(ka) + 2J1 cos(ka) cos Q 2J2 cos(2ka) + 2J2 cos(2ka) cos(2Q ) M21 = 2J1 cos(ka) + 2J1 cos(ka) cos Q + 2J2 cos(2ka)

(45)

+2J2 cos(2ka) cos(2Q ) − 4(n − 1)V1 20 cos Q −4(n − 1)V2 20 cos(2Q ) − (n − 1)U20 . The condition for the existence of nontrivial solutions of (44) can be found by the following quadratic equation for the frequency 2

(M11 − ) − M12 M21 = 0.

(46)

The modulational instability can be measured by the power gain which at any wave number is defined by gain(G) = Im( ) and the growth rate is defined by 1

gain(G) = Im(M12 M21 ) 2 ,

(47)

where Im denotes the imaginary part of the frequency . The existence of localized modes can be possible if the constant amplitude of the plane wave solution is unstable. Eq. (47) determines the stability/instability of a plane wave as function of the wave number Q in discrete ultracold bosons chain. The importance of the linear stability analysis is first to determine the instability domain in space parameter and secondly to predict the evolution of the amplitude of a modulation at the onset of the instability.

586

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

Fig. 4. Gain spectrum of modulational instability in the (k, Q) plane for different values of the parameter J2 /J1 controlling the strength of the NNNI in the lattice, where a = 1, 0 = 0.1, V1 = 5 ×10−3 , V2 = 10−3 , U = 9 ×10−3 : (a) J2 /J1 = 0.0; (b) J2 /J1 = 0.2; (c) J2 /J1 = 0.4 and (d) J2 /J1 = 0.8.

It is not difficult to express the initial amplitude 0 with respect to a threshold amplitude T0 . We find that the plane wave introduced in the chain is modulationally unstable if the initial amplitude |0 | exceeds the threshold amplitude |T0 |, which is defined as follows |0 |2 ≥ |T 0 |2 = −

2J1 cos(ka)(1 + cos Q ) + 2J2 cos(2ka)(1 + cos(2Q )) (n − 1)(4V1 cos(ka) + 4V2 cos(2ka) + U)

(48)

From (48), we mention that the parameter J1 representing the nearest-neighbor hopping strength and J2 denoting the hopping strength to the next-nearest-neighbor may influence the magnitude of the threshold amplitude for a given wave number k. Fig. 4 portrays the stability/instability regions in the (k, Q) plane and also presents the corresponding effect of the parameter J2 /J1 . As we will notice in Fig. 4, the dark bluish area stands for a region where the nonlinear plane waves are stable whereas the region with bright yellowish orange area corresponds to the region of modulational instability. It is clearly shown in Fig. 4 for the wave number k = 0 that, the plane wave is stable of any perturbation Q. According to [30], this implies that a dark-type localized mode may exist in the ultracold bosons in optical lattices system. This case corresponds to our analytical results when PQ < 0. Consequently, we conclude that the formation of the stationary bright type intrinsic localized mode can not appear at the center of the Brillouin zone of this quantum system. However in the absence of the parameter J2 representing the hopping strength to the next-nearest-neighbor, we realized that for k = /a, the modulational instability can occur (see Fig. 4(a)). This result agrees with the formation of the bright soliton solution predicted in our analytical investigation presented in Section 3. In the presence of NNNI, i.e when the parameter J2 is not null, as we can see in Fig. 4(b), the plane wave begins to be unstable for a wave number k = kmax as predicted analytically in Section 3. When the value of the parameter J2 /J1 becomes sufficiently large as we can see in Fig. 4(c) and (d), the shape of the region of modulational instability changes significatively as a very interesting impact of the introduction of NNNI in the chain. Therefore, the inclusion of the NNNI in chain is responsible to the appearance of new region of modulation instability for the k = kmax carrier wave. Fig. 5 depicts the plot 3D of the modulational instability gain in the (k, Q) plane for different values of the parameter J2 /J1 controlling the strength of the NNNI in the lattice, verifying the effect of modulational instabilities regions described in Fig. 4 where we notice that, the power gain of modulational instabilities regions increases when the value J2 /J1 becomes strong. Furthermore, we analyze the instability growth rates for different values of the parameter J2 /J1 in order to understand well the dependence of the modulation instability for k = kmax and k = keBZ the plane wave, respectively. In Fig. 6, we portrays the growth rates curves of the modulation instability for the k = kmax and k = keBZ plane wave. We find that, the instability growth rate may be significantly affected by the intensity of the parameter J2 /J1 controlling the strength of the NNNI in the lattice. It is clearly shown in Fig. 6(a) corresponding to the case of k = kmax that, in the absence of the NNNI (J2 /J1 = 0), the instability growth rates is nil. Consequently for this case, there is not formation of modulational instability region as we can

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

587

Fig. 5. Gain spectrum of modulational instability in the (k, Q) plane for different values of the parameter J2 /J1 controlling the strength of the NNNI in the lattice, where a = 1, 0 = 0.1, n = 15, V1 = 5 ×10−3 , V2 = 10−3 , U = 9 ×10−3 : (a) J2 /J1 = 0.0; (b) J2 /J1 = 0.2; (c) J2 /J1 = 0.4 and (d) J2 /J1 = 0.8.

Fig. 6. Growth rate versus the wave number of the perturbation Q for different values of J2 /J1 , where a = 1, 0 = 0.1, n = 15, V1 = 5 ×10−3 , V2 = 10−3 , U = 9 ×10−3 : (a) k = kmax and (b) k = keBZ .

588

Z.I. Djoufack et al. / Optik 164 (2018) 575–589

confirm by the shape of the Fig. 4(a). However, in the presence of the NNNI (J2 /J1 = / 0) and when this value is not large, the modulational instability region appears in two positions along the Q wave number axis. The modulational instability region appears in only one position along the Q axis and the amplitude of the instability growth rate slightly increases with an increasing value of the NNNI parameter, while the maximum value of the gain appears when the value of J2 /J1 becomes very large. In Fig. 6(b), we realize that the instability growth rate is not nil in the absence of the NNNI in the chain for the case of k = keBZ as we showed in Fig. 4(a), where one can notice that while the value of J2 /J1 increases, the amplitude curves of the instability growth rate is shifted along the Q axis. This scenario not longer occur when J2 /J1 ≥ 2 and the formation of the maximum value of the gain appears for J2 /J1 ≥ 0.24. 6. Conclusion In this work, we have investigated quantum breathers associated with modulational instability in 1D ultracold boson in optical lattices involving NNNI. By means of the time dependent Hartree approximation combined with the semi-discrete multiple-scale method and the equation of motion for the single-boson wave function was found. It was found that the dispersion spectrum in the half of the Brillouin zone displays different regions of possible nonlinear localized excitations for J2 /J1 > 1/4 when the NNNI are included in the system, for J2 /J1 < 1/4, the stable region (k2 , /a) disappears and k2 takes the value /a. For J2 /J1 > 1/4, in the wave number k = kmax , the bright type intrinsic localized mode occurs above the top of the linear wave spectrum and for J2 /J1 < 1/4 at the edge of the Brillouin zone in the wave number k = keBZ , the bright type intrinsic localized mode also appears above the top of the linear wave spectrum. With the help of those bright stationary localized single-boson wave functions, we established quantum breather states and their corresponding quantized energy in each case. The performance of our analytical analysis is accredited by numerical simulations. As Kivshar pointed out that the existence of the stationary localized solutions in nonlinear lattices is directly linked to the discrete modulational instability of plane waves, therefore, in order to verify the stability analytical analysis in our model was carried out. We have shown that the presence of the NNNI affect significatively the shape of the region of modulational instability and it is responsible of the appearance of new region of modulational instability that occurs for the k = kmax carrier wave. The formation conditions of the modulational instability region predicted by the analytical analysis of stationary localized solutions in the wave number k = kmax and k = keBZ are in good agreement with the forecast respectively. Finally, we have found that, the instability growth rate may be significantly affected by the intensity of the parameter J2 /J1 controlling the strength of the NNNI in the lattice. We noticed that for the case of k = kmax that, in the absence of the NNNI (J2 /J1 = 0), the instability growth rates is nil and consequently for this case, there is not formation of modulational instability region. However, in the presence of the NNNI (J2 /J1 = / 0) and when this value is not large, the modulational instability region appears in two positions along the Q wave number axis. The modulational instability region appears in only one position along the Q axis and the amplitude of the instability growth rate slightly increases with an increasing value of the NNNI parameter, while the maximum value of the gain appears when the value of J2 /J1 becomes very large. We also realize that the instability growth rate is not nil in the absence of the NNNI in the chain for the case of k = keBZ and while the value of J2 /J1 increases, the amplitude curves of the instability growth rate are shifted along the Q axis. This scenario not longer occur when J2 /J1 ≥ 2 and the formation of the maximum value of the gain appears for J2 /J1 ≥ 0.24. We point out that results obtained in this work may improve our understanding of the generation of stationary localized states in nonlinear optical lattices by the modulational instability. Another interesting feature of this is that the long range interaction can be very important in many-body quantum systems for their practical implication in quantum computing. We expect that such a theoretical work on ultracold gas of bosonic atoms in an optical lattice needs experimental investigation. Acknowledgements Zacharie Isidore Djoufack would like to thank the African Institute for Mathematical Sciences (AIMS) South Africa for hospitality, where this work was initialized. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13]

M.H. Anderson, et al., Science 269 (1995) 198. B. DeMarco, D.S. Jin, Science 285 (1999) 1703. I. Bloch, J. Dalibard, W. Zwerger, Many-body physics with ultracold gases, Rev. Mod. Phys. 80 (2008) 885. M. Greiner, et al., Nature (London) 415 (2002) 39. O. Morsch, M. Oberthaler, Dynamics of Bose–Einstein condensates in optical lattices, Rev. Mod. Phys. 78 (2006) 179. G.K. Brennen, C.M. Caves, P.S. Jessen, I.H. Deutsch, Phys. Rev. Lett. 82 (1999) 1060. D. Jaksch, C. Bruder, J.I. Cirac, C.W. Gardiner, P. Zoller, Cold bosonic atoms in optical lattices, Phys. Rev. Lett. 81 (1998) 3108. A.C. Scott, J.C. Eilbeck, H. Gilhøj, Quantum lattice solitons, Physica D 78 (1994) 194. I. Bloch, Nat. Phys. 1 (2005) 23. S. Flach, C.R. Willis, Discretes breathers, Phys. Rep. 295 (1998) 181–264. S. Flach, A.V. Gorbach, Discrete breathers-advances in theory and applications, Phys. Rep. 467 (2008) 1. R. Sarkar, B. Dey, Energy localization and transport in two-dimensional Fermi-Pasta-Ulam lattice, Phys. Rev. E 76 (2007) 016605. K. Yoshimura, Existence and stability of discrete breathers in diatomic Fermi Pasta Ulam type lattices, Nonlinearity 24 (2011) 293.

Z.I. Djoufack et al. / Optik 164 (2018) 575–589 [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47]

589

S. Flach, Existence of localized excitations in nonlinear Hamiltonian lattices, Phys. Rev. E 51 (1995) 1503. A.J. Sievers, S. Takeno, Intrinsic localized modes in anharmonic crystals, Phys. Rev. Lett. 61 (1988) 970. E. Trías, J.J. Mazo, T.P. Orlando, Discrete breathers in nonlinear lattices: experimental detection in a Josephson array, Phys. Rev. Lett. 84 (2000) 741. M. Sato, A.J. Sievers, Direct observation of the discrete character of intrinsic localized modes in an antiferromagnet, Nature 432 (2004) 486. M. Sato, B.E. Hubbard, A.J. Sievers, B. Ilic, D.A. Czaplewski, H.G. Craighead, Observation of locked intrinsic localized vibrational modes in a micromechanical oscillator array, Phys. Rev. Lett. 90 (2003) 044102. R.S. MacKay, S. Aubry, Proof of existence of breathers for time-reversible or Hamiltonian networks of weakly coupled oscillators, Nonlinearity 7 (1994) 1623. V. Fleurov, Discrete quantum breathers: what do we know about them? Chaos 13 (2003) 676. B. Tang, Quantum two-breathers formed by ultracold bosonic atoms in optical lattices, Int. J. Theor. Phys. (2016), http://dx.doi.org/10.1007/s10773-015-2903-9. Z.I. Djoufack, A. Kenfack-Jiotsa, J.-P. Nguenang, Quantum signature of breathers in 1D ultracold bosons in optical lattices involving next-nearest neighbor interactions, Int. J. Mod. Phys. B 31 (2017) 1750140, 21 p. C. Verdozzi, M. Cini, Phys. Rev. B 51 (1995) 7412. R. Lai, A.J. Sievers, Phys. Rev. B 57 (1998) 3433. F. Baronio, S. Chen, P. Grelu, S. Wabnitz, M. Conforti, Baseband modulation instability as the origin of rogue waves, Phys. Rev. A 91 (2015) 033804. G. Gor, T. Macrí, A. Trombettoni, Modulational instabilities in lattices with powerlaw hoppings and interactions, Phys. Rev. E 87 (2013) 032905. L. Wang, Y.-J. Zhu, F.-H. Qi, M. Li, R. Guo, Modulational instability, higher-order localized wave structures, and nonlinear wave interactions for a nonautonomous Lenells–Fokas equation in inhomogeneous fibers, CHAOS 25 (2015) 063111. L. Wang, J.-H. Zhang, Z.-Q. Wang, C. Liu, M. Li, F.-H. Qi, R. Guo, Breather-tosoliton transitions, nonlinear wave interactions, and modulational instability in a higher-order generalized nonlinear Schrödinger equation, Phys. Rev. E 93 (2016) 012214. L. Wang, J.-H. Zhang, C. Liu, M. Li, F.-H. Qi, Breather transition dynamics, Peregrine combs and walls, and modulation instability in a variable-coefficient nonlinear Schrödinger equation with higher-order effects, Phys. Rev. E 93 (2016) 062217. Y.S. Kivshar, Localized modes in a chain with nonlinear on-site potential, Phys. Lett. A 173 (1993) 172. J. Stockhofe, P. Schmelcher, Modulational instability and localized breather modes in the discrete nonlinear Schrödinger equation with helicoidal hopping, Physica D 328 (2016) 9. C.B. Tabi, A. Mohamadou, T.C. Kofane, Modulational instability in the anharmonic Peyrard–Bishop model of DNA, Eur. Phys. J. B 74 (2010) 151. Y.S. Kivshar, M. Peyrard, Modulational instabilities in discrete lattices, Phys. Rev. A 46 (1992) 3198. R. Asaoka, H. Tsuchiura, M. Yamashita, Y. Toga, Density modulations associated with the dynamical instability in the Bose–Hubbard model, J. Phys. Soc. Jpn. 83 (2014) 124001. R. Asaoka, H. Tsuchiura, M. Yamashita, Y. Toga, Dynamical instability in the Bose–Hubbard model, Phys. Rev. A 93 (2016) 013628. D. Giuliano, D. Rossini, P. Sodano, A. Trombettoni, XXZ spin-1/2 representation of a finite-U Bose–Hubbard chain at half-integer filling, Phys. Rev. B 87 (2013) 035104. D. Jaksch, P. Zoller, The cold atom Hubbard toolbox, Ann. Phys. 315 (2005). I. Bloch, Ultracold quantum gases in optical lattices, Nat. Phys. 1 (2005). Y. Gao, F. Han, Effet of the next-nearest-neighbor hopping in optical lattices, Mod. Phys. Lett. B 22 (2008) 33. E. Wright, J.C. Eilbeck, M.H. Hays, P.D. Miller, A.C. Scott, The quantum discrete self-trapping equation in the Hartree approximation, Physica D 69 (1993) 18. J.W. Negele, The mean-field theory of nuclear structure and dynamics, Rev. Mod. Phys. 54 (1982) 913. Y. Lai, H.A. Haus, Time-dependent Hartree approximation, Phys. Rev. A 40 (1989) 844. P. Bonche, S. Koonin, J.W. Negele, One-dimensional nuclear dynamics in the time-dependent Hartree–Fock approximation, Phys. Rev. C 13 (1976) 1226. E.M. Wright, Quantum theory of soliton propagation in an optical fiber using the Hartree approximation, Phys. Rev. A 43 (1991) 3836. M. Remoissenet, Low-amplitude breather and envelope solitons in quasi-one-dimensional physical models, Phys. Rev. B 33 (1986) 2386. M. Remoissenet, Waves called solitons, in: Concepts and Experiments, 2nd ed., Springer, Berlin, Heidelberg/New York, 1996. R.A. Pinto, S. Flach, Phys. Rev. B 77 (2008) 024308.