Non-Markovian dynamics without using quantum trajectory
aa r X i v : . [ qu a n t - ph ] J a n Non-Markovian dynamics without using quantum trajectory
Chengjun Wu, Yang Li, Mingyi Zhu, and Hong Guo ∗ CREAM Group, State Key Laboratory of Advanced Optical Communication Systems and Networks (Peking University)School of Electronics Engineering and Computer Science, Peking University, Beijing 100871, China (Dated: August 28, 2018)Open quantum system interacting with structured environment is important and manifests non-Markovian behavior, which was conventionally studied using quantum trajectory stochastic method.In this paper, by dividing the effects of the environment into two parts, we propose a determin-istic method without using quantum trajectory. This method is more efficient and accurate thanstochastic method in most Markovian and non-Markovian cases. We also extend this method to thegeneralized Lindblad master equation.
PACS numbers: 03.65.Yz, 42.50.Lc
When an open quantum system interacts with environ-ment, it experiences decoherence and dissipation whichlead to loss of information. Such open quantum system isdepicted by a reduced density matrix which shows non-unitary evolution. On the other hand, the environment isclassified as Markovian with no memory effect, and non-Markovian with memory effect. In Markovian case, sincethere is no memory effect, the quantum trajectory basedMonte Carlo wave function (MCWF) method [1–3] andquantum state diffusion (QSD) method [4, 5] are applied.However, in non-Markovian case, due to memory effect,the information lost by the system during the interactionwith the environment will come back to the system in alater time and so shows much more complicated behav-iors than Markovian case.Non-Markovian systems are important for their ap-plications to many fields of physics, such as quan-tum information processing [6, 7], quantum optics [8],solid state physics [9], and chemical physics [10]. Re-cently, non-Markovian behaviors have also been studiedin biomolecules where the molecules are embedded ina solvent and/or in a protein environment [11]. Sincethere is no true pure state quantum trajectory due tothe memory effect [12], the quantum trajectory basedMarkovian methods do not work. Thus, doubled Hilbertspace (DHS) method [13], triple Hilbert space (THS)method [14], non-Markovian QSD method [15, 16], andnon-Markovian quantum jump (NMQJ) method [17, 18]are proposed to solve the non-Markovian dynamics of thesystem where the memory effect is taken into account.However, in order to obtain high accuracy, all these meth-ods, which are based on stochastic simulations, need tofulfill a large number of realizations and is very time-consuming. So, new methods which are more efficientand accurate are highly desired.In this paper, a deterministic method without us-ing quantum trajectory is proposed to solve the non-Markovian dynamics. The influence of the environmenton the system is divided into two parts, i.e., the non- ∗ Correspondence author: [email protected] unitary evolution of the states and the probability flowbetween these states. Moreover, we also extend thisapproach to the generalized Lindblad master equationwhich can deal with some strong coupling cases [19]. Thealgorithm and numerical efficiency are given, which showthat our method is more efficient and accurate than thosebased on stochastic simulation in most Markovian andnon-Markovian cases.The dynamics of the non-Markovian system is gov-erned by the following master equation [8]˙ ρ ( t ) = 1 i ~ [ H s , ρ ( t )] + X j γ j ( t ) C j ( t ) ρ ( t ) C † j ( t ) − X j γ j ( t ) { ρ ( t ) , C † j ( t ) C j ( t ) } , (1)where H s is the system Hamiltonian including the Lambshift, C j ( t ) are the jump operators which induce changes[e.g., jump from state ψ α ( t ) to ψ α ′ ( t ) i.e., | ψ α ′ ( t ) i = C j ( t ) | ψ α ( t ) i / || C j ( t ) | ψ α ( t ) i|| ] in the system, and γ j ( t )are the decay rates which may take negative values forsome time intervals. The reduced density matrix can bewritten as [17] ρ ( t ) = N eff X α =1 p α ( t ) | ψ α ( t ) ih ψ α ( t ) | , (2)where p α ( t ) is the probability of the system being in thestate | ψ α ( t ) i at time t . Further, it should be pointedout that the effective number of the states N eff is de-termined by C j ( t )’s [18], P N eff α =1 p α ( t ) = 1 and that thestate | ψ α ( t ) i is normalized.To solve the dynamics of the system, one should knowthe time evolution of | ψ α ( t ) i and its probability p α ( t ). Inour method, the time evolution of the state | ψ α ( t ) i is thesame as that in NMQJ [17]. In NMQJ, the probability p α ( t ) is calculated in a stochastic way by using quan-tum trajectory to N ensemble members. In our method,however, the evolution of probability p α ( t ) is given in adeterministic way:˙ p α ( t ) = − X j Γ jα ( t ) p α ( t ) + X ( α ′ ,j ) ′ Γ jα ′ ( t ) p α ′ ( t ) , (3)where Γ jα ( t ) = γ j ( t ) k C j ( t ) | ψ α ( t ) ik and P ( α ′ ,j ) ′ representsthe summation over all the pairs ( α ′ , j ) satisfying | ψ α ( t ) i = C j ( t ) | ψ α ′ ( t ) i / k C j ( t ) | ψ α ′ ( t ) ik . One finds thatthe probability of the state, p α ( t ), changes via the mech-anism of jumps for “out” ( α → α ′ ) and “in” ( α ′ → α ),respectively.The numerical simulation corresponding to Eq. (3) isstraightforward: p α ( t + δt ) = p α ( t ) − δt X j Γ jα ( t ) p α ( t )+ δt X ( α ′ ,j ) ′ Γ jα ′ ( t ) p α ′ ( t ) . (4)Note that there is no stochastic noise and no need toconsider the sign of the decay rate during the simula-tion. Additionally, the p α ( t )’s in our method do repre-sent the probability of the system actually being in thecorresponding pure state ensemble.Consider a particular transition: | ψ α ′ ( t ) i = C j ( t ) | ψ α ( t ) i / || C j ( t ) | ψ α ( t ) i|| , then the correspondingprobability change takes the form: p α ( t + δt ) = p α ( t ) − δtp α ( t )Γ jα ( t ) ,p α ′ ( t + δt ) = p α ′ ( t ) + δtp α ( t )Γ jα ( t ) . (5)When the decay rate γ j ( t ) is positive or negative, theprobability flow is from | ψ α ( t ) i to | ψ α ′ ( t ) i or reversed.This has been mentioned in Ref.[17]. However, it is moreexplicit in our method. From Eq. (5), it is clear that, inthe negative decay region, the amount of probability flowonly depends on the target state and the probability ofthe system being in the target state. This is similar tothe situation in NMQJ [17], where the jump probabilityin the negative decay region is proportional to the num-ber of particles in the target state. These indicate thatthe trajectory of a particle in NMQJ can not be inter-preted as true trajectory since the jump process dependson the status of other particles in the system. Becausetrue pure state quantum trajectories do not exist in thenon-Markovian dynamics [12], it is not necessary to cal-culate p α ( t ) in a stochastic way.Next, we extend our method to the recently proposedgeneralized Lindblad master equation which can solve thedynamics of some highly non-Markovian systems[19], ddt ρ i = − i [ H i , ρ i ] + X jλ (cid:18) R ijλ ρ j R ijλ † − { R jiλ † R jiλ , ρ i } (cid:19) , (6)where i, j = 1 , , · · · , n , H i are any Hermitian operators,and R ijλ are any system operators. It should be indicatedthat ρ ( t ) = n P i =1 ρ i ( t ).The i th density matrix is decomposed as: ρ i ( t ) = N ieff X α =1 p αi ( t ) | ψ αi ( t ) i h ψ αi ( t ) | , (7) where N ieff is determined in the same way in Eq. (2)by taking all the jump operators R ijν , s and all the states ψ αj ( t ) , s in each ρ j ( t ) into consideration.The evolution of state | ψ αi ( t ) i is governed by the non-linear differential equation[21] i ddt | ψ αi ( t ) i = ˆ G ( ψ αi )( t ) | ψ αi ( t ) i , (8)where ˆ G ( ψ αi )( t ) = H i − i P jν R jiν † R jiν + i P jν (cid:13)(cid:13) R jiν | ψ αi ( t ) i (cid:13)(cid:13) . By combining Eqs. (6), (7),(8) and noting that (cid:12)(cid:12) ψ i ( t ) (cid:11) (cid:10) ψ i ( t ) (cid:12)(cid:12) , (cid:12)(cid:12) ψ i ( t ) (cid:11) (cid:10) ψ i ( t ) (cid:12)(cid:12) , · · · , (cid:12)(cid:12)(cid:12)(cid:12) ψ N ieff i ( t ) (cid:29) (cid:28) ψ N ieff i ( t ) (cid:12)(cid:12)(cid:12)(cid:12) are linearly independent, theevolution of p αi ( t ) is given by˙ p αi ( t ) = − X jν Γ jiνα p αi ( t ) + X ( j,ν,α ′ ) ′ Γ ijνα ′ p α ′ j ( t ) , (9)where Γ ijνα = (cid:13)(cid:13) R ijν (cid:12)(cid:12) ψ αj ( t ) (cid:11)(cid:13)(cid:13) and P ( j,ν,α ′ ) ′ represents thesummation over all the pairs ( j, ν, α ′ ) satisfying | ψ αi ( t ) i = R ijν (cid:12)(cid:12)(cid:12) ψ α ′ j ( t ) E.(cid:13)(cid:13)(cid:13) R ijν (cid:12)(cid:12)(cid:12) ψ α ′ j ( t ) E(cid:13)(cid:13)(cid:13) . It can be easily seen thatby setting n = 1 and taking the decay rates γ ( t ) into theequation, Eq. (9) degenerates to Eq. (3). Example : Detuned Jaynes-Cummings model. –Consider a system with a two-level atom in a detuneddamped cavity, which is governed by the time convolu-tionless master equation [8]˙ ρ ( t ) = − i S ( t ) { σ + σ − , ρ ( t ) } + γ ( t ) { σ − ρ ( t ) σ + − σ + σ − ρ ( t ) − ρ ( t ) σ + σ − } . (10)The spectral density of the cavity is supposed to be ofLorentzian profile, i.e., J ( ω ) = γ λ π [( ω − ∆ − ω ) + λ ] , where∆ = ω − ω c is the detuning between the cavity modeand the atom. To second order approximation, theLamb shift and the decay rate take the form [8] S ( t ) = γ λ ∆ λ +∆ { − e − λt [cos(∆ t ) + λ ∆ sin(∆ t )] } , γ ( t ) = γ λ λ +∆ { − e − λt [cos(∆ t ) − ∆ λ sin(∆ t )] } . In this model, there is onlyone jump operator C = σ − = | g i h e | , which is a lower-ing operator. We assume that ρ (0) = | ψ (0) i h ψ (0) | andchoose | ψ (0) i = (4 | e i +3 | g i ) /
5. Acting the jump opera-tor on the state | ψ (0) i , we get | ψ (0) i = | g i . Accordingto Eq. (4), at time t + δt , the probabilities become p ( t + δt ) = p ( t ) − δtp ( t )Γ ( t ) ,p ( t + δt ) = p ( t ) + δtp ( t )Γ ( t ) , (11)where Γ ( t ) = γ ( t ) |h e | ψ ( t ) i| . In this example, ρ ee ( t ) is proportional to the energy ofthe system and p ( t ) represents the probability for onephoton being in the environment. Although p ( t ) and p ( t ) can be solved analytically, in order to illustrate ourmethod, we use Eq. (11) to do the simulation. The pa-rameters are chosen as ∆ = 12 λ, γ λ = 4 , λδt = 0 . FIG. 1: (color online) Dynamics of detuned Jaynes-Cummings model. The initial state is | ψ (0) i = (4 | e i +3 | g i ) / λ, γλ = 4 , λδt = 0 . | ψ ( t ) i and | ψ ( t ) i . (b) The population of the excited state ρ ee (ini-tially higher line) and the absolute value of the coherence ρ eg (initially lower line) with three methods: analytic (red solidcurve), our method (blue long-dashed curve) and NMQJ (with N = 10 particles in the system, green dash-dot curve). bility flow. We can see from Fig. 1 (a) and (b) thatwhen the probability flow gets reversed, the energy andcoherence of the atom increase. These show explicitly thememory effect that the reduced system restores the in-formation lost earlier. In Fig. 1 (b), the result of NMQJ(with N = 10 particles in the system) is also given,which shows that our method is more accurate. Example 2: Application to generalized Lindblad mas-ter equation. – To illustrate our method for this kind ofequation, we consider a two-state system coupled to anenvironment consisting of two energy bands, each witha finite number of evenly spaced levels. This may beviewed as a spin coupled to a single molecule or a singleparticle quantum dot [20]. By using time-convolutionlessprojection operator technique, to the second order, thegeneralized Lindblad master equation takes the form [21] ddt ρ = Z t dt h ( t − t )[2 γ σ + ρ σ − − γ { σ + σ − , ρ } ] ,ddt ρ = Z t dt h ( t − t )[2 γ σ − ρ σ + − γ { σ − σ + , ρ } ] , (12)where γ i h ( t − t ) , ( i = 1 , , is the environment corre-lation function with h ( t ) = δε sin ( δεt/ π ( δεt/ where δε is thewidth of the upper and lower energy bands. The reduceddensity matrix for the system is given by ρ = ρ + ρ . We assume that ρ (0) = | e i h e | and ρ (0) = 0. The pa-rameters are chosen as δǫ = 0 .
31 and γ = γ = 1. InFig. 2 we compare the results of our method, analyticalsolution and Monte Carlo simulation which is based onthe unraveling of the master equation (with N = 10 tra-jectories) [21]. Apparently, our method is more accurate FIG. 2: (color online) A two-state system coupled to an envi-ronment consisting of two energy bands. Comparison of ourmethod (blue long-dashed curve) and Monte Carlo simulation(with N = 10 trajectories, green dash-dot curve) to analyt-ical result (red solid curve). The parameters are δǫ = 0 . ,γ = γ = 1 and time step δt = 0 . than Monte Carlo simulation method.According to Eqs. (4), we only need to calculate N eff states and change the probabilities deterministically. Thetime cost is almost determined by the calculation of N eff states. However, the evolution of N eff states is indepen-dent with each other, so we can calculate them parallelly.In addition, if the jump operators can be represented bysparse matrixes, we only need to calculate the evolutionof the states appearing in the decomposition of ρ (0) anduse the jump operators to obtain other states. Moreover,since the sign of the decay rate makes no difference dur-ing the simulation, in non-Markovian case, our methodis as efficient as it behaves in Markovian case.Similar to our method, the NMQJ method [17, 18]needs to calculate N eff states. However, in addition tothat, NMQJ has to consider the sign of the decay ratesand generate N random numbers ( N ≫ N eff ) to decidethe jump process at each time step δt . Apparently, ourmethod is more efficient than NMQJ in any case.In Markovian case, the MCWF [1] and QSD [4] methodneed to realize a large number of trajectories for everystate appearing in the decomposition of ρ (0). Whenthe number of these trajectories is larger than N eff ,which is always the case, our method is more efficientthan them. In non-Markovian case, the DHS method[13], THS method [14] and non-Markovian QSD method[15, 16] all introduce additional cost for computationalefficiency compared to MCWF or QSD. However, in non-Markovian case, our method is as efficient as it behaves inMarkovian case. Thus, when the number of these trajec-tories is larger than N eff , our method is obviously moreefficient than them, too.As for the accuracy, since there is no statistical noisein our method and the error caused by finite time step δt is the same, compared with all the methods basedon stochastic simulation, our method is more accurate.Actually, our method is the limit case when the numberof realizations in the stochastic based methods tends toinfinite.In conclusion, by dividing the influence of the environ-ment on the system into two parts, i.e., the non-unitaryevolution of these states and the probability flow betweenthem, we propose a deterministic method to solve thenon-Makovian dynamics. Compared with the methodbased on stochastic simulation, our method has advan- tages in efficiency and accuracy. Additionally, we ex-tended this approach to the generalized Lindblad masterequation , which is useful to solve the dynamics of somehighly non-Markovian systems.This work is supported by the Key Project of the Na-tional Natural Science Foundation of China (Grant No.60837004). [1] J. Dalibard, Y. Castin, and K. Mølmer, Phys. Rev. Lett. , 580 (1992);[2] H. Carmichael, An Open System Approach to Quan-tum Optics , Lecture Notes in Physics (Springer-Verlag,Berlin, 1993), Vol. m18.[3] M. B. Plenio and P. L. Knight, Rev. Mod. Phys. , 101(1998).[4] N. Gisin and I. C. Percival, J. Phys. A , 5677 (1992); , 2233 (1993); , 2245 (1993).[5] I. Percival, Quantum State Diffusion (Cambridge Univer-sity Press, Cambridge, England, 2002).[6] M. A. Nielsen and I. L. Chuang,
Quantum Computationand Quantum Information (Cambridge University Press,Cambridge, England, 2000)[7] Y. Li, J. Zhou, and H. Guo, Phys. Rev. A , 012309(2009).[8] H.-P. Breuer and F. Petruccione, The Theory of OpenQuantum Systems (Oxford University Press, Oxford,2002).[9] See, e.g., C. W. Lai, P. Maletinsky, A. Badolato, andA. Imamoglu, Phys. Rev. Lett. , 167403 (2006), andreferences therein.[10] J. Shao, J. Chem. Phys. , 5053 (2004); A. Pomyalov and D. J. Tannor, J. Chem. Phys. , 204111 (2005)and references therein.[11] P. Rebentrost, R. Chakraborty, and A. Aspuru-Guzik, J.Chem. Phys. , 184102 (2009).[12] H. M. Wiseman and J. M. Gambetta, Phys. Rev. Lett. , 140401 (2008).[13] H.-P. Breuer, B. Kappler, and F. Petruccione, Phys. Rev.A , 1633 (1999).[14] H.-P. Breuer, Phys. Rev. A , 012106 (2004).[15] W. T. Strunz, L. Di`osi, and N. Gisin, Phys. Rev. Lett. , 1801 (1999).[16] J. T. Stockburger and H. Grabert, Phys. Rev. Lett. ,170407 (2002).[17] J. Piilo, S. Maniscalco, K. H¨ark¨onen, and K.-A. Suomi-nen, Phys. Rev. Lett. , 180402 (2008).[18] J. Piilo, S. Maniscalco, K. H¨ark¨onen, and K.-A. Suomi-nen, Phys. Rev. A , 062112 (2009).[19] H.-P. Breuer, Phys. Rev. A , 022103 (2007).[20] J. Gemmer and M. Michel, Europhys. Lett. , 1 (2006).[21] M. Moodley and F. Petruccione, Phys. Rev. A79