Entanglement sudden death and revival in quantum dark-soliton qubits
Muzzamal I. Shaukat Eduardo. V. Castro Hugo Terças
quant-ph cond-mat.quant-gas





CeFEMA, Instituto Superior Técnico, Universidade de Lisboa, Lisboa, Portugal
University of Engineering and Technology, Lahore (RCET Campus), Pakistan
muzzamalshaukat@gmail.com

CeFEMA, Instituto Superior Técnico, Universidade de Lisboa, Lisboa, Portugal

Instituto de Plasmas e Fusão Nuclear, Lisboa, Portugal
hugo.tercas@tecnico.ulisboa.pt




67.85.Hj 42.50.Lc 42.50.-p 42.50.Md 



We study the finite time entanglement dynamics between two dark soliton qubits due to quantum fluctuations in a quasi one dimensional Bose-Einstein condensates. Recently, dark solitons are proved to be an appealing platform for qubits due to their appreciably long life time. We explore the entanglement decay for an entangled state of two phonon coherences and the qubits to be in the diagonal basis of so called Dicke states. We observe the collapse and revival of the entanglement, depending critically on the collective damping term but independent of the qubit-qubit interaction for both the initial states. The collective behavior of the dark soliton qubits demonstrate the dependence of entanglement evolution on the interatomic distance.
 

 
 


Entanglement sudden death and revival in quantum dark-soliton qubits
    Hugo Terças
    December 30, 2023
====================================================================





§ INTRODUCTION

Two basic features to distinguish the classical world from the quantum world are the superposition and entanglement, being the cornerstone of the rapidly developing field of quantum information and computation <cit.>. In recent years, the interaction between qubits and the environment has attracted a great deal of attention motivated by the possibility of controlling and exploiting the entanglement dynamics <cit.>. The occurence of spontantaneous emission due to the coupling of qubits with the reservoir, leading to the irreversible loss of quantum information therein encoded, has been regarded as the main obstacle in practical usage of entanglement. 

Much work has been done to understand the decoherence dynamics for a pair of qubits interacting with different reservoirs <cit.>. As it is currently known, two spatially separated qubits, initially prepared in a product of two pure states, get entangled as time evolves, producing to the creation of a so-called transient entanglement in the system <cit.>. Conversely, Yu and Eberly discovered a finite-time disentanglement of two initially-entangled qubits in contact with pure dissipative environments <cit.>. This effect is currently known under the name of entanglement sudden death, and has been confirmed with experiments performed both with photonic <cit.> and atomic systems <cit.>. 

In this context, Bose-Einstein condensates (BECs) have attracted a great deal of interest during the last decades, since the macroscopic character of the wavefunction allows BECs to display pure-state entanglement at macroscopic scales <cit.>. A scheme to generate entangled states in a bimodal BEC has been firstly announced in Refs. <cit.>; the investigation of the macroscopic superposition based on matter waves has been achieved with BEC Josephson junctions <cit.>. Moreover, light scattering with BECs have been used to enhance their nonlinear properties in superradiance experiments <cit.>, and to show the possibility of matter wave amplification <cit.> and nonlinear wave mixing <cit.>. Entanglement dynamics for coupled BECs have been investigated in <cit.>. 

Another important manifestation of the macroscopic nature of BEC is the dark soliton (DS), a structure resulting from the detailed balance between the dispersive and nonlinear effects, appearing also in other physical systems <cit.>. The dynamics and stability of DSs in BECs have been a subject of intense research over the last decade <cit.>.  The dynamical evolution of DS entanglement and how its stability is affected by quantum fluctuations has been studied in Ref. <cit.>. Collision-induced entanglement between fast moving matter-wave solitons using the Born approximation has been studied in BECs displaying attractive interactions <cit.>. Moreover, the study of collective aspects of soliton gases <cit.> bring DSs towards applications in many body physics <cit.>. In a recent publication, we have shown that DSs trapping an impurity can behave as qubits in quasi-1D BECs <cit.>, being excellent candidates to store information as consequence of their appreciably long lifetimes (∼ 0.01-1s). DSs qubits thus offer an appealing alternative to quantum optical of solid-state platforms, as information processing involves only phononic degrees of freedom: the quantum excitations on top of the BEC state.

Ghasemian et al. <cit.> described the collapse-revival phenomenon, after studying the dynamics of BEC atoms interacting with a single-mode laser field. The possibility of generating entangled Schrödinger cat states by using a BEC trapped in a double-well potential has been investigated <cit.>. In these works, it is reported that the revival of the initial state can be used as an unambiguous signature of the coherent macroscopic superposition, as opposed to an incoherent mixture. Also, two-impurity qubits surrounded by a BEC reveal the influence of quantum reservoirs on effects such as sudden death, revival, trapping and generation of entanglement <cit.>.


In this paper, we investigate the dynamics of the entanglement produced between two dark-soliton qubits in a quasi-one dimensional Bose-Einstein condensate. As described in Ref <cit.>, the qubits are produced by trapping impurities inside the potential created by the dark solitons. Moreover, the phonons (quantum fluctuations on top of the background density) play the role of a quantum reservoir. We show the occurrence of entanglement sudden death by computing the time evolution of the Wootters' concurrence and showing that it vanishes for a finite time. We further demonstrate that the concurrence dynamics critically depends on the distance between the DS qubits. 		

The paper is organized as follows: In sec. II, we model the mean-field dynamics of the BEC and the impurities by using the Gross-Pitaevskii (GP) and Schrödinger equations, respectively. We compute the coupling between the DS-qubits and the phonons. In sec. III,  we describe the Markovian master equation and extract the density matrix elements for the collective DS qubit states. The use of concurrence as a measure of the entanglement is discussed in sec. IV. Finally, a summary and discussion of the investigation is presented in Sec. V.







§ THEORETICAL MODEL AND QUANTUM FLUCTUATIONS


At the mean field level, the system is governed by Gross-Pitaevskii (GP) equations,

    iħ∂Ψ/∂ t   =   -ħ ^2/2m_1∂^2Ψ/∂ x^2+g|Ψ| ^2Ψ
    +χ|Φ| ^2Ψ, 
    
    iħ∂Φ/∂ t   =   -ħ ^2/2m_2∂^2Φ/∂ x^2+χ|Ψ| ^2Φ,

where, χ is the BEC-impurity coupling constant, g is the BEC particle-particle interaction strength, and m_1 and m_2 denote the BEC particle and impurity masses, respectively. Here, the discussion is restricted to repulsive interacions (g>0) where the dark solitons are assumed to be not disturbed by the presence of impurities. To achieve this, the impurity gas is considered to be much less massive than the BEC particles and sufficiently dilute, i.e. |Ψ|^2 ≫|Φ|^2 (experimental realization can be found in <cit.>). Therefore, the soliton behave like a potential for the impurities (considered to be free particles) i.e.,


    iħ∂Φ/∂ t=-ħ ^2/2m_2∂^2Φ/∂ x^2+χ|ψ _ sol| ^2Φ,

where a singular nonlinear solution corresponds to the ith (i=1,2) soliton profile is given by <cit.> 


    ψ^(i)_ sol(x)=√(n_0)tanh[ (x -x_i)/ξ],

where x_i=± d/2 determines the position of the solitons, ξ =ħ /√(mn_0g_11) is the healing length and n_0∼ 10^8-10^9 m^-1 is the typical linear density. For well separated solitons, d≫ξ, the internal level structure of the two qubits is assumed to be equal. Each qubit is characterized by its ground g_i and excited levels e_i (i=1,2), separated by a gap frequency ω _0=ħ(2ν -1)/(2mξ ^2) (see Fig. (<ref>)). Here, 


    ν=-1+√(1+4χ/g)


is a parameter controlling the number of bound states created by the DS, which operates as a qubit (two-level system) in the range 0.33≤ν <0.80 <cit.>. Typical experimentally accessible conditions provide longitudinal and transverse trap frequencies of ω _z/2π∼ (15-730) Hz ≪ω_r/2π= (1-5) kHz, and the corresponding length l_z=(0.6-3.9) μm <cit.>. More recent experiments lead to much less pronounced trap inhomogeneities by creating much larger traps, l_z∼ 100 μm <cit.>. 



 §.§ Quantum fluctuations
 

The total BEC quantum field includes the DSs wave functions and quantum fluctuations, Ψ_i(x)=ψ^(i)_ sol(x)+δψ_i(x), where 

 
    δψ_i(x)=∑_k (u^(i)_k(x) b_k +v^*(i)_k(x)b^†_k ).


Here, b_k are the bosonic operators verifying the commutation relation [b_k,b^†_q]=δ_k,q. The amplitudes u_k(x) and v_k(x) satisfy the normalization condition | u_k(x)| ^2 -| v_k(x)| ^2=1 and are explicitly given by <cit.> 

    . u^(i)_k(x)= e^ik(x-x_i)√(1/4πξ)μ/ϵ_k
    . ×
          . [((kξ)^2+2ϵ_k/μ)(kξ/2+itanh(x-x_i/ξ))+kξ/
    cosh^2(x-x_i/ξ)] .,


    . v^(i)_k(x)= e^-ik(x-x_i)√(1/4πξ)μ/ϵ_k
    . ×
          . [((kξ)^2-2ϵ_k/μ)(kξ/2+itanh(x-x_i/ξ))+kξ/
    cosh^2(x-x_i/ξ)] ..

The total Hamiltonian then reads

    H=H_ q+H_ p+H_ int.

The term H_ q describes the dark-solitons (qubits) Hamiltonian, which is given by


    H_ q=∑^2_i=1ħω _0σ _z ^(i)


with σ_z^(i)=a_1^(i)† a_1^(i)- a_0^(i)† a_0^(i) being the effective spin operator of the respective qubit. The phonon (reservoir) Hamiltonian is given by

    H_ p =∑_k ϵ _kb_k^† b _k,


with the Bogoliubov spectrum ϵ _k=μξ√(k^2(ξ^2k^2+2)) and the chemical potential μ=gn_0.  Finally, the interaction Hamiltonian H_ int can be expliclity written as 

    H_ int=∑_i,jχ∫ dxΦ _j^†Ψ _i^†Ψ _iΦ _j,

where Φ_j(x) describes the impurity wave function in the presence of DS potential and spannable interms of bosonic operators a_l   

    Φ_j(x)=∑_n=0^1 ϕ^(j)_n(x) a^(j)_n,

with the ground state ϕ_0(x)= sech[ (x -x_i)/ξ]/√(2ξ) and the excited state ϕ_1(x)=i √(3)tanh[ (x -x_i)/ξ]ϕ_0(x).
Therefore, Eq. (<ref>) can be decomposed as

    H_ int=H_ int^(0)+H_ int^(1)+H_ int^(2),


containing zeroth, first and second order terms in the bosonic operators b_k and b_k^†. The higher-order term ∼𝒪(b_k^2) is ignored, in consistency with the Bogoliubov approximation performed in Eq. (<ref>) owing to the small depletion of the condensate. This approximation is well justified in the case of two-level systems, as the inexistent higher excited states cannot be populated via two-phonon processes. The first part of Eq. (<ref>) corresponds to 


    H_ int^(0)=n_0χ∑_i=1^2∑_n,n^'=0^1a_n^† (i) a_n^'^(i)f_n,n^'^(i) .

with f^(i)_n,n^'=∫ dxϕ _n^† (i) (x) ϕ_n^'^(i)(x)  tanh^2[(x -x_i)/ξ], which can be omitted by renormalizing the qubit frequency ω̃_0≈ω_0+n_0χ. The first order term is given by


    H_ int^(1)   =   ∑_k,i=1^2∑_n,n^'=0^1a_n^(i)†a^(i)_n^'[ b_kg_n,n^'^(i)(k).	
       +   . b_k^†g_n,n^'^(i)∗(k)]


where


    g_n,n^'^i,j(k) = √(n_0)χ∫ dxϕ^(j) † _n(x) ϕ^(j)_n^'(x) tanh(x -x_i/ξ) u^(i)_k.
 

Eq. (<ref>) contains both interband (n ≠ n') and intraband (n = n') terms. However, within the RWA, the qubit transition can only be driven by near-resonant phonons, for which the intraband terms | g_00^(i)(k) | and | g_11^(i)(k)| are much smaller than the interband term  | g_01^(i)(k)|=| g_10^(i)(k)^*|. As such, we obtain 


    H_ int^(1)   =   ∑_k,i=1^2( g^(i)(k)σ^(i) _+b_k+ g^(i) ∗(k) σ^(i) _-b_k^†)+ h.c. .


Here, g^(i)(k)=g_n,n^'^i,j(k), σ^(i) _+=a_1^(i)†a^(i)_0 and vice versa. The counter-rotating terms proportional to b_k σ_- and b_k^†σ_+ that do not conserve the total number of excitations are dropped by invoking the RWA. The accuracy of such an approximation can be verified in Ref. <cit.>, where it is shown that the emission rate γ is much smaller than the qubit transition frequency ω_0 for DS qubits.
 


§ MASTER EQUATION


We derive the master equation (see appendix) to describe the dynamics of the DS qubit density matrix ρ_q after taking trace over the phonon's degrees of freedom <cit.>

    ∂ρ_q(t)/∂ t   =    -i/ħ[H_q,ρ_q(t)] -∑^2_i≠ jη_ij[σ_+^iσ_-^j,ρ_q(t)]
       +   ∑^2_ij=1Γ_ij[σ_-^jρ_q(t)σ_+^i-1/2{σ_+^iσ_-^j,ρ_q(t) }],

where

    Γ_ij   =    2L∫_0^∞dk g^(i)_kg_k^(j)*δ(ω_k-ω_0) 
    η_ij   =   L/2π℘∫_0^∞dk g^(i)_kg_k^(j)*1/(ω _k-ω _0),

with ℘ standing for the principal value of the integral.

For i=j, Γ_ii≡γ denotes the spontaneous emission rate of the qubits due to the Bogoliubov excitations (phonons); for i≠ j, Γ_ij≡Γ is the collective damping resulting from the mutual exchange of phonons. The coherent term η_ij≡η represents the interaction between DS qubits. Both the coherent and incoherent terms are dependent on the distance d between DS qubits. Fig. (<ref>) depicts the dependence of the collective damping  Γ and the coherent interaction η as a function of the DS qubit distance d. For large separations, i.e. d ≫ξ, both parameters are very small (i.e. Γ = η≈ 0).

The main concern of the present work is the study of the time evolution of the entanglement. To study this, the most adequate sates are the collective, or the so-called Dicke, states <cit.>


    | g⟩   =   | g_1,g_2⟩ , 
    |±⟩   =   (| e_1,g_2⟩±| g_1,e_2⟩)/√(2),  
    | e⟩   =   | e_1,e_2⟩ ,


as schematically represented in Fig. (<ref>). The density matrix of Eq. (<ref>) becomes


    ρ =(                                                   
    [ ρ _ee ρ _e+ ρ _e- ρ _eg; ρ _+e ρ _++ ρ _+- ρ _+g; ρ _-e ρ _-+  ρ _– ρ _-g; ρ _ge ρ _g+ ρ _g- ρ _gg ]
    ) ,

where ρ_ij =|ψ_i⟩⟨ψ_j|. The elements of the above density matrix can be determined by using Eq. (<ref>),

    ρ_ee(t)   =    e^-2γ tρ_ee(0)
    ρ_++(t)   =   e^-(γ+Γ) tρ_++(0) 
       +   (γ+Γ)/(γ-Γ)(e^-(γ+Γ) t-e^-2γ t)ρ_ee(0)
    ρ_–(t)   =   e^-(γ-Γ) tρ_–(0)
       +   (γ-Γ)/(γ+Γ)(e^-(γ-Γ) t-e^-2γ t)ρ_ee(0)
    ρ_eg(t)   =    e^-γ tρ_eg(0)
    ρ_+-(t)   =    e^-(γ-2iη) tρ_+-(0) 
    ρ_e+(t)   =    e^-1/2(3γ +Γ-2iη) tρ_e+(0)
    ρ_e-(t)   =    e^1/2(3γ +Γ-2iη) tρ_e-(0)
    ρ_g+(t)   =    e^-1/2(γ +Γ-2iη) tρ_g+(0)
       +   (γ+Γ)/(γ+2iη)2e^-1/2(2γ +Γ) tsinh(t/2(γ+2iη))ρ_+e(0) 
    ρ_g-(t)   =    e^-1/2(γ -Γ+2iη) tρ_g-(0)
       -   (γ-Γ)/(γ-2iη)2e^-1/2(2γ -Γ) tsinh(t/2(γ-2iη))ρ_-e(0)

with ρ_jk=ρ^*_kj and ρ_gg=1-ρ_ee-ρ_++-ρ_–. Eq. (<ref>) depicts that all transition rates to and from
the state ρ_++ are equal to γ+Γ while from state ρ_– are γ-Γ. Therefore, the state ρ_++ decays with an enhanced
(superradiant) rate and ρ_– with a reduced (subradiant) rate.



§ MEASUREMENT OF ENTANGLEMENT

The amount of entaglement can be determined by using the Wootter's concurrence <cit.>


    C(t)=max(0, √(ε)_1-√(ε)_2-√(ε)_3-√(ε)_4),


where the ε_i's are the eigenvalues in decreasing order of magnitude of the matrix

    ζ=ρσ_y⊗σ_yρ^∗σ_y⊗σ_y.

Here, ρ^∗ represents the complex conjugate of ρ and σ_y is the Pauli matrix. Depending on the initial state, concurrence can reach a value equal to zero assymptotically or at some finite time. It is interesting to observe that locally equivalent initial states with the same concurrence can disentangle at different times, depending on the parameters Γ and η. In what follows, we investigate this aspect. 



 §.§ Entangled State


Let assume that, initially, both or neither of the DS qubits are excited, i.e., the qubits are choosed to be prepared initially in an entangled state


    |Ψ⟩=√(1-α)| g ⟩ + √(α)| e⟩ ,


with 0≤α≤ 1. Therefore, the density matrix Eq. (<ref>) becomes,

    ρ =(                                                   
    [ ρ _ee     0     0 ρ _eg;     0 ρ _++     0     0;     0     0  ρ _–     0; ρ _ge     0     0 ρ _gg ]
    ) .

The eigenvalues of the respective matrix ζ are thus given by


    √(ε_1,2)   =   √(ρ_ee(t)ρ_gg(t))±|ρ_ge(t) |
    √(ε_3,4)   =   1/2[ ρ_++(t)+ρ_–(t)] ±1/2[ρ_++(t)-ρ_–(t)].


It is easy to verify that, depending on the largest eigenvalue (either √(ε_1) or √(ε_3)), the concurrence C(t) can be defined in two alternative ways, i.e.,


    C_1(t)   =   2|ρ_ge(t) | - [ρ_++(t)+ρ_–(t)] 
    
    C_2(t)   =   |ρ_++(t)-ρ_–(t)| - 2√(ρ_ee(t)ρ_gg(t)),


where C_1(t)  measures the entanglement produced by the state of Eq. (<ref>) with the necessary condition ρ_eg(t) ≠ 0, while C_2(t) provides the entanglement formed by the states |±⟩. Notice that the positiveness of C_2(t) is guaranteed if the latter states are not equally populated. 
At t=0, the system is entangled by the amount C_1(0)=2√(α(1-α)). By inspecting Eqs. (<ref>), (<ref>) and (<ref>), it is possible to observe that DS qubits radiating independently (Γ=0) can not be entangled by following the criterion C_2(t) because ρ _++(t)=ρ _–(t). Thus, the time for the entanglement death due to spontaneous emission can be found via the condition C_1(t)=0, which provides


    t_ death=1/γ ln(α/α - √(α(1-α))).


It is also pertinent to mention here that the condition for a finite-time disentanglement for independent DS qubits is α>1/2 (see Fig. (<ref>)). 
The situation changes when we allow the qubits to interact. In this case, the entanglement death is followed by its revival at a larger time, t_ revival. Fig. (<ref>) depicts the time evolution of the concurrence for the collective interactive system. It is shown that the entanglement dies as consequence of the spontaneous emission, but revives after a time t_ revival≃ 8/γ, for α≃ 1/4 and qubits placed at a distance d = 6ξ/5. After a careful inspection of Eq. (<ref>), it is observed that the concurrence C_1(t)<0 at long times. Therefore, finite time (t_ revival) entanglement is determined by following C_2(t) which yields. 


    t_ revival=2/3Γ ln(4γ/√(α)(γ-Γ)).


Moreover, It can be analyzed from Fig. (<ref>) that entanglement vanishes around the time at which ρ_–(t) is maximally populated, and that it is not undergo any revival. The latter is due to the impartiality between the term  ρ_eg(t) and ρ_–(t). In other words, ρ_eg(t) and ρ_–(t) go to almost zero at long times while the population of ρ_++(t) accumulates on the time scale t=1/(γ+Γ) which is sufficiently large at Γ≈ -0.5γ. The collective damping Γ≈ 0 at farthest distance d ≃ 5ξ for which DS qubits act like an independent qubits (fig. (<ref>)).









 §.§ Mixed State


We now consider a two-qubit system to be initially prepared in a diagonal basis of the collective states. Therefore, the initial density matrix has the form


    ρ(0) =1/3(                                                   
    [   α   0   0   0;   0   2   0   0;   0   0   0   0;   0   0   0 1-α ]
    ) ,

The initial concurrence is determined by C_2(0)=2(1-√(α(1-α)))/3, and the sudden-death time for independent DS qubits can be described by using criterion C_2(t)=0, which provides


    t_ death=1/γ ln(α/√(3α^2+5α))-(1+α)).


It is obvious from Eq. (<ref>) that the ESD is possible only for α≳ 1/3 (see Fig. (<ref>)). The time evolution of the concurrence for interacting qubits is depicted in Fig. (<ref>). It is observed that the entanglement first decay and then revives for α≳1/2 at d≈ 4ξ. ESD happens at t_ death∼ 0.75/γ, while entanglement revival is obtained at t_ revival∼ 1.4/γ. 









 §.§ Experimental estimates
 

We consider a quasi 1D BEC of ^85Rb with a chemical potential (μ) of a few kHz. This yields a qubit gap frequency ω _0/2π∼ 0.5 kHz, the spontaneous decay rate γ/2π∼ 29 Hz and a collective decay Γ/2π∼ 6 Hz, at d ∼ 6ξ/5∼ 1 μm. These rates validate a posteriori the RWA and Markovian approximations. Therefore, the sudden death-time for the maximally entangled state is t_ death∼ 19 ms and the revival time is t_ revival∼ 35 ms, where the period Δ t ∼ 16 ms is the dark period, i.e. the time interval during which C(t)=0. For large separations, ESD occurs at t_ death∼ 2 ms due to the balanced population of ρ_eg and ρ_– whereas, for the mixed initial state, the dark period occurs for ∼ 3.6 ms. Prolonging entanglement is essential for practical realization of quantum information and computation protocols based on entanglement. Therefore, the dark period of entanglement can be delayed or averted by carrying out local unitary operation on qubits <cit.>.  




§ SUMMARY AND DISCUSSION


To summarize, we investigate the finite-time disentanglement, or the entanglement sudden death, between two dark-soliton qubits produced in a quasi one-dimensional BEC. We derive the master equation and extracted the time evolution of the relevant density matrix elements. The Wooter's concurrence is used as a measure of entanglement and we show the collapse and revival behavior, depending on the collective damping and on the initial state. For initial an entangled state, the concurrence can not be revived at large distances in the range of 2-5 μm due to the impartial behavior of the populated states, while it revives for a mixed state. Therefore, it can be concluded that the collective behavior of the dark soliton qubits reveals the dependence of entanglement evolution on the interatomic distance and it becomes quite different from that of independent dark soliton qubits.






§ DERIVATION OF THE MASTER EQUATION


We begin by writing the total Hamiltonian by spanning the Hilbert space as


    H(t)= H_q⊗ I + I ⊗ H_p + H_ int(t),


The key ingredient for the application of the Born-Markov approximation is the assumption that H_q⊗ H_p is small compared to the remaining terms, so that a perturbative treatment of the interaction is possible. To make this more explicit, we take H_0=H_q + H_p and move into the interaction picture

    H_I(t)=e^i H_0 tH_Ie^-i H_0 t.

We make the Born approximation and assume that the density operator factorizes at all times as


    ρ(t)≅ρ_q(t) ⊗ρ_p,


where the reservoir density operator is assumed to be time independent i.e., ρ_p=ρ_p(0). Therefore, within the interaction picture, the density operator evolves according to


    dρ (t)/dt=-[H_I(t),ρ(t)],

and the respective formal integral solution is given by

    ρ(t)=ρ(0)-∫_0^t[H_I(s),ρ(s)]ds.


By plugging Eq. (<ref>) into Eq. (<ref>), and by taking the partial trace over reservoir degrees of freedom (phonons), we obtain


    dρ_q(t)/dt=-∫_0^tds tr_P[H_I(t),[H_I(s),ρ_q(s)ρ_p]],


This equation is called the Redfield equation where the term Tr_P[H_I(t),ρ(0)] is disregarded. Further, we make use of the Markov approximation to put Eq. (<ref>) into a more amenable form. This assures that the behaviour of ρ_q(t) is local in time. This master equation is still depends on the choice of the initial state. However, making the substitution s → t-τ  and letting the upper integration limit to go to infinity, we obtain


    dρ_q(t)/dt=-∫_0^∞dτ tr_P[H_I(t),[H_I(t-τ),ρ_q(t) ρ_p]]


The latter is the Born-Markov master equation. To gain a little more insight into the structure of Eq. (<ref>), it is useful to be more specific about the form of the interaction picture Hamiltonian 

    H_I(t)= S^†B+ SB^†,

where S^†(t)=∑_i=1^2σ^(i) _+ = e^ω _0tS^†, B(t)=∑_k g(k)e^-ω _ktb and vice versa. Here, we use the identity


    e^αASe^-αA= S+α[A,S]+α^2/2! [A,[A,S]]+⋯

Moreover, we invoke the cyclic property of trace to write the revervoir correlation function as Tr_P(b_kb^†_qρ_P)=δ_k,q. As such, the Born-Markov master Eq. (<ref>) can be finally rewritten as 


    dρ_q(t)/dt   =    -γ[ ρ_q(t) S^†S-Sρ_q(t)S^†]
       -   γ[ρ_q(t)SS^†-Sρ_q(t)S^†] +  h.c.,

where 

    γ  = ∑_k g(k)g(k)^*∫_0^tdτ e^-(ω _k-ω _0)(t-τ).

The sum over the phonon k- modes can be computed by taking the continuum limit

    ∑_k→∫_0^∞ D(k)dk,

where D(k)= L/2 π is the density of states, L is the size of the system, and 

    ∫_0^tdτ e^-i(ω _k-ω _0)(t-τ)=(πδ(ω_k-ω_0)-i℘/(ω _k-ω _0)).

Transforming Eq. (<ref>) back in the Schrödinger picture, we finally obtain

    dρ_q(t)/dt   =    -i/ħ[H_q,ρ_q(t)]-∑^2_i≠ jη_ij[σ_+^iσ_-^j,ρ_q(t)] 
       +   ∑^2_ij=1Γ_ij[σ_-^jρ_q(t)σ_+^i. 
       -   .1/2{σ_+^iσ_-^j,ρ_q(t) }]

Eq. (<ref>) is the final form of the master equation used in this work.











§ ACKNOWLEDGEMENTS


The authors acknowledge the support from DP-PMI programme and Fundação para a Ciência e a Tecnologia (Portugal), namely through the scholarship number SFRH/PD/BD/113650/2015 and the grant number SFRH/BPD/110059/2015 and No. IF/00433/2015. E.V.C. acknowledges partial support from FCT-Portugal through Grant
No. UID/CTM/04540/2013.

99


Nielsen2000  M. A. Nielsen and I.L. Chuang,  Quantum Computation and Quantum  Information,  Cambrige  University  Press: Cambrige, 2000.
Malinovsky2004 V. S. Malinovsky and I. R. Sola, Phys. Rev. Lett. 93, 190502 (2004).
Plenio1999M. B. Plenio, S. F. Huelga, A. Beige, and P. L. Knight, Phys. Rev. A 59, 2468 (1999).
Serafini2006 A. Serafini, S. Mancini, and S. Bose, Phys. Rev. Lett. 96, 010503 (2006).
Ficek2002 Z. Ficek and R. Tana ́s, Phys. Rep. 372, 369 (2002).


Ficek2003Z. Ficek and R. Tana ́s, Journal of Modern Optics, 50 (2003).
Ficek2004Z. Ficek and R. Tana ́s, Journal of Optics B, Quantum semiclass. Opt, 6 (2004).
Ficek2003aZ. Ficek and R. Tana ́s, FortSchr. Phys. 51 (2003).
Eberly2004 T. Yu and J. Eberly, Phys. Rev. Lett. 93, 140404 (2004).
Eberly2006 T. Yu and J. Eberly, Phys. Rev. Lett. 97, 140403 (2006).
Almeida2007 M. P. Almeida et al., Science 316, 579 (2007).
Laurat2007 J. Laurat et al., Phys. Rev. Lett. 99, 180504 (2007).

Davis1995 K. B. Davis et al , Phys. Rev. Lett. 75, 3969 (1995).
Anderson1995M. H. Anderson, J.R. Ensher, M.R. Matthews, C.E. Wieman, and E.A. Cornell, Science 269, 198 (1995).

Parkins1998 A. S. Parkins and D.F. Walls, Phys. Rep. 303, 1 (1998).

Hines2003 A. P.  Hines, R. H. McKenzie and G. J. Milburn,  Phys. Rev. A 67, 013609 (2003).
Micheli2003 A. Micheli, D. Jaksch, J. I. Cirac, and P. Zoller, Phys. Rev. A 67, 013607 (2003).

Gordon1999D. Gordon and M. C. Savage,  Phys. Rev. A 59, 4623 (1999).
Cirac1998J. I. Cirac, M. Lewenstein, K. Mølmer and P. Zoller,  Phys. Rev. A 57, 1208 (1998).
Ruostekoskil1998J. Ruostekoskil, M. J. Collett, R. Graham and D. F. Walls,  Phys. Rev. A 57, 511 (1998).


Inouye1999 S. Inouye et al., Science 285, 571 (1999).
Kozuma1999 M. Kozuma et al., Science 286, 2309 (1999).
Deng1999L. Deng et al.,  Appl. Phys. Lett. 68, 218 (1999).
Sanz2003 L. Sanz,M. R. Angela and K. Furaya, J. Phys. A: Math. Gen. 36, 9737 (2003).
Ng2009H. T. Ng and S. Bose, New Jour. Phys. 11, 043009  (2009). 
kivshar Y. S. Kivshar and G. P. Agrawal, Optical Solitons: From Fibers to Photonic Crystals (Academic Press, San Diego, USA, 2003).


Denschlag J. Denschlag, J. E. Simsarian, D. L. Feder, C. W. Clark, L. A. Collins, J. Cubizolles, L. Deng, E. W. Hagley, K. Helmerson, W. P. Reinhart, S. L. Rolston, B. I. Schneider, and W. D. Phillips, Science 287, 97 (2000).
burger S. Burger, S. Dettmer, W. Ertmer, K. Sengstock, A. Sanpera, G. V. Shlyapnikov, and M. Lewenstein, Phys. Rev. Lett. 83, 5198 (1999).



Dziarmaga2003J. Dziarmaga, Z. P. Karkuszewski, and K. Sacha, J. Phys. B: At. Mol. Opt. Phys. 36, 1217  (2003). 
Jackson2007B. Jackson, N. P. Proukakis, and C. F. Barenghi, Phys. Rev. A 75, 051601  (2007). 



Mishmash2009 R. V. Mishmash and L. D. Carr, Phys. Rev. Lett. 103, 140403 (2009).
Lewenstein2009M. Lewenstein and B. Malomed, New J. Phys. 11, 113014 (2009).

gael2005 G. A. El and A. M. Kamchatnov, Phys. Rev. Lett. 95, 204101 (2005).

tercas H. Terças, D. D. Solnyshkov and G. Malpuech, Phys. Rev. Lett. 110, 035302 (2013); ibid 113, 036403 (2014).
muzzamal2017 M. I. Shaukat, E. V. Castro and H. Terças, Phys. Rev. A 95, 053618 (2017).
Ghasemian2016E. Ghasemian, M.K. Tavassoly, Phys. Lett. A, 380 3262  (2016).
Fisher1989M. P. A. Fisher, B. P. Weichman, G. Grinstein, and D. S. Fischer, Phys. Rev. B 40, 546 (1989).
Huang2006 Y. P. Huang and M. G. Moore,  Phys. Rev. A 73, 023606 (2006).
McEndoo2013 S. McEndoo, P. Haikka, G. De Chiara, G.M.Palma and S. Maniscalco,  Euro Phys. Lett. 101, 60005 (2013).

zakharov72 V. E. Zakharov and A. B. Shabat, Sov. Phys. JETP 34, 62 (1972); ibid 37, 823 (1973).

huang G. Huang, J. Szeftel, and S. Zhu, Phys. Rev. A 65, 053605 (2002). 

parker2004 N. Parker, Numerical Studies of Vortices and Dark Solitons in atomic Bose Einstein Condensates, Ph.D Thesis (2004).
schmiedmayer2010 P. Krüger, S. Hofferberth, I. E. Mazets, I. Lesanovsky, and J. Schmiedmayer, Phys. Rev. Lett. 105, 265302 (2010). 

Martinez2011  J. C. Martinez, Euro. Phys. Lett. 96, 14007 (2011). 


Agarwal1974 G. S . Agarwal, Quantum Statistical Theories of Spontaneous Emission and their Relation to other Approaches (Springer Tracts in Modern Physics vol 70) ed G Hohler (Berlin: Springer)  (1974).
Lehmberg1970R. H. Lehmberg, Phys. Rev. A 2, 883  (1970); 889  (1970).



Dicke1954 R. H. Dicke, Phys Rev. 93, 99 (1954).
Wootters1998W. K. Wootters, Phys. Rev. Lett. 80, 2245 (1998)
Ficek2008Z. Ficek, R. Tanas, Phys. Rev. A 77, 054301  (2008).

Rau2008A. R. P. Rau, M. Ali and G. Alber, Euro Phys. Lett. 82, 40002  (2008).
Singh2017 A. Singh, S. Pradyumna, A.R.P. Rau and U. Sinha,  J. Opt. Soc. Am. B 34, 681 (2017).








