Exact dynamics of quantum correlations of two qubits coupled to bosonic baths

# Exact dynamics of quantum correlations of two qubits coupled to bosonic baths

Chen Wang    Qing-Hu Chen Department of Physics, Zhejiang University, Hangzhou 310027, P. R. China
Center for Statistical and Theoretical Condensed Matter Physics, Zhejiang Normal University, Jinhua 321004, P. R. China
July 26, 2019
###### Abstract

The counter rotating-wave term (CRT) effects from the system-bath coherence on the dynamics of quantum correlation of two qubits in two independent baths and a common bath are systematically investigated. The hierarchy approach is extended to solve the relevant spin boson models with the Lorentzian spectrum, the exact dynamics for the quantum entanglement and quantum discord (QD) are evaluated, and the comparisons with previous ones under the rotating-wave approximation are performed. For the two independent baths, beyond the weak system-bath coupling, the CRT essentially changes the evolution of both entanglement and QD. With the increase of the coupling, the revival of the entanglement is suppressed dramatically and finally disappears, and the QD becomes smaller monotonically. For the common bath, the entanglement is also suppressed by the CRT, but the QD shows quite different behaviors, if initiated from the correlated Bell states. In the non-Markovian regime, the QD is almost not influenced by the CRT and generally finite in the long time evolution at any coupling, while in the Markovian regime, it is significantly enhanced with the strong coupling.

###### pacs:
03.67.Mn,03.65.Ud,03.65.Yz

## I introduction

Quantum correlation, originated from the coherent superposition in quantum mechanics einstein1 (); bell1 (), has been attracting persistent attention ranging from quantum information, quantum optics, to many-body physics nielsen1 (); taylor1 (); modi1 (). As the characteristic trait of the quantum correlation proposed by Schrödinger schrodinger1 () without the classical counterpart, quantum entanglement played the crucial role in the highly efficient quantum computation and quantum information processing vidal1 (); horodechi1 (). It has also been applied to identify the phase transition in quantum many-particle systems amico1 (). However, there exist exceptions where the entanglement seems unnecessary, such as in the Grover search algorithm grover1 (); ahn1 () and the deterministic quantum computation with one pure qubit (DQC)knill1 (). It is evident that other kind of the nonclassical quantum correlation rather entanglement dominates the computational implementation of DQC lanyon1 (). Zurek and Vedral described one particular measure to quantify all nonclassical correlation, termed as quantum discord (QD) ollivier1 (); zurek1 (); henderson1 (). Only when the QD disappears, the information can be safely obtained by locally measuring the distantly separate subsystem without disturbing the bipartite quantum system, where the corresponding state is fully classical. Recently, the QD has successfully been used to study the quantum phase transition of the critical systems werlang (); dillenschneider1 (); maziero1 (); allegra1 (); chen (); wang (), operational meanings for quantum processors gu1 (), and state preparation walther1 ().

It is known that the realistic system inevitably interacts with the environment, resulting in the decoherence, which is of fundamental importance in the quantum information processing and measurement breuer1 (); weiss1 (). Specially, Yu and Eberly  yu1 () found that the dynamical behavior of the entanglement is significantly different from the single qubit in a pair of qubits separately coupling with Markovian baths, which can be modeled from the well-known spin boson models leggett1 (). It shows finite time disentanglement called as entanglement sudden death (ESD) yu1 (). However, this was investigated under quantum master equation with the Born and Markov approximation, which assumes that system-bath interaction is weak and the relaxation of the bath is much faster than that of the system. Then further works have been done in the non-Markovian regime in the framework of the rotating-wave approximation (RWA) to reveal novel characters of the quantum correlations compared to the Markovian effect bellomo1 (); bellomo2 (); maniscalco1 (); werlang1 (); wang1 (); franco1 (), which is only valid under the condition of weak system-bath interaction. Besides, entanglement dynamics beyond the weak coupling spin boson models has also been exploited based on the optimal polaron transformation cao1 (), where the higher order terms are neglected. It is pointed out that the corresponding accuracy in the wide parameter zone should be carefully analyzed lee1 (). Therefore, It is highly desirable to solve the dynamics without performing any approximation.

Beyond both the Born-Markov approximation and RWA, Tanimura et al. tanimura1 () have developed an efficient hierarchy approach, which has been later extensively applied to some chemical and biophysical systems ishizaki1 (); xu1 (). Recently, this method has also been extended to study the dynamics of entanglement for two qubits coupled to a common bath arend1 (); ma1 (). However, the dynamics of QD, one recently mostly studied quantum correlations has not been investigated beyond the Born-Markov approximation and the RWA. Moreover, the exact study of the dynamics of both entanglement and QD for two qubits coupled to their own baths are not available in the literature, to the best of our knowledge. This case is designed for two remote qubits each interacting locally with its own environment. Such a physical condition is also relevant to the quantum information science and quantum computation.

In the present paper, to give a comprehensive picture of the two qubits coupled to their own bath and common bath, we will extend the hierarchy equation to these two kinds of baths to study the dynamics of the pairwise entanglement and the QD in both Markovian and non-Markovian regimes. The effect of the counter-rotating term on these dynamics in a wide coupling regime will be systematically studied. This paper is organized as follows. In Sec. II we briefly review the definition of the QD. In Sec. III we exactly solve the reduced qubits density matrix of two independent spin boson models, and investigate its quantum correlation compared with that in RWA. In Sec. IV the quantum correlation of two separate qubits interacting with a common bath is analyzed. Finally, a summary is given in Sec. V.

## Ii Quantum Discord

The QD is one main route to fully measure the nonclassical correlation, which can be extended from the classical information theory ollivier1 (); zurek1 (); henderson1 (); modi1 (). It is interpreted as the difference of two quantum mutual correlations of subsystems and , before and after local measurement operated on subsystem  ollivier1 (); zurek1 (). The total quantum correlation of two subsystems is determined by their joint density matrix as , where . The von Neumann entropy . The other quantum mutual information , which is equivalent to in classical information frame, is derived by locally measuring subsystem with a complete set of orthonormal projectors , where is one outcome state of . After this measurement, the joint density is reduced to conditional counterpart as , where the corresponding probability is . Then the mutual information based on specific performance is shown as . Hence, we should maximize it to capture all classical correlation by obtaining . Finally, the QD is defined as

 D(ρ)=I(ρ)−J(ρ). (1)

The quantum nature of the system can be explicitly observed from the QD. It is only zero in pure classical correlation condition, and has the same values as the entanglement in pure state. While it remains finite even in separable mixed states, where the corresponding entanglement may completely disappear ollivier1 ().

## Iii two independent spin-boson models

The Hamiltonian of two independent qubits coupled two independent bosonic baths and is given by , reads

 Hs = ∑a=A,Bωa2σaz, Hb = ∑a=A,B;kωakb†a,kba,k, Hsb = ∑a=A,B;kσax(gakba,k+ga∗kb†a,k). (2)

where , is the Pauli operator with Zeeman energy of spin , creates (annihilates) one boson with frequency in bath , and is the coupling strength between the system and the bath.

The spectra density of the bosonic bath is selected as the following Lorentz distribution

 Ja(ω)=12πλaγ2a(ω−ωa0)2+γ2a. (3)

which can describe an bosonic field inside an imperfect cavity mode with the system-bath coupling strength . The corresponding correlation function at zero temperature is shown as (31)

 Caa(ω)=λaγa2exp[−(γa+iωa0)t]. (4)

Then by applying the Feynman-Vernon influence functional procedure in Appendix, the hierarchy equation is derived by (A)

 ∂∂tρ→n,→m(t) = −(iH×S+→n⋅→μA+→m⋅→μB)ρ→n,→m(t) +∑k=1,2(−1)kQ×Aρ→n+→ek,→m(t) +λAγA4∑k=1,2nk[QoA+(−1)k+1Q×A]ρ→n−→ek,→m(t) +∑k=1,2(−1)kQ×Bρ→n,→m+→ek +λBγB4∑k=1,2mk[QoB+(−1)k+1Q×B]ρ→n,→m−→ek(t).

The commutation relations and . The unit vectors , , and , . in the system-bath interaction represents . At the initial state, only is not zero. For the simplest condition in the following calculations, the two subsystems are identical as , , .

For comparison, we thereby briefly review the main results derived in the RWA, where system-bath interaction becomes . When the initial system reduced density matrix in the basis is form

 ρAB(0)=⎛⎜ ⎜ ⎜ ⎜⎝ρ11(0)00ρ14(0)0ρ22(0)ρ23(0)00ρ32(0)ρ33(0)0ρ41(0)00ρ44(0)⎞⎟ ⎟ ⎟ ⎟⎠, (6)

then we have

 ρ11(t) = ρ11(0)P2t, ρ22(t)=ρ22(0)Pt+ρ11(0)Pt(1−Pt), ρ33(t) = ρ33(0)Pt+ρ11(0)Pt(1−Pt), ρ44(t) = 1−[ρ11(t)−ρ22(t)−ρ33(t)], ρ14(t) = ρ14(0)P(t), ρ23(t)=ρ23(0)P(t), (7)

where

 Pt=e−γt[cos(Rt)+γ2Rsin(Rt)]2 (8)

with . It is known that the bath relaxation time is , and the system-bath correlation time is . From the time dependence of the function , we immediately know that the Markovian regime corresponds to where decays exponentially, while the non-Markovian regime to where shows oscillating behavior describing the coherent process between the system and the bath. This non-Markovian condition has been realized in cavity quantum electrodynamics having Rydberg atoms inside the Fabry-Perot cavities with  kuhr1 (). We should point out here that the Markovian and non-Markovian regimes are only defined in the framework of the RWA, which can not be well defined in the exact treatment without the RWA. Throughout the present paper, we still use such notations. Doing so is only for comparisons and analysis.

As usual, the initial Bell states with anti-correlated spins and correlated spins will be considered in the present paper

 |Φ(α)⟩ = α|10⟩+√1−α2|01⟩, (9) |Ψ(α)⟩ = α|00⟩+√1−α2|11⟩, (10)

which both obey form in Eq. (6) during the evolution. So it is straightforward to derive the concurrence, a pairwise entanglementbellomo1 ()

 CΦ(t) = max{0,2α√1−α2Pt}, (11) CΨ(t) = max{0,2√1−α2Pt[α−√1−α2(1−Pt)]}.

The analytical expression for QD in this case is quite complicated and also only consists of and , which is implied in Refs. wang1 (); fanchini1 (); wang0 ().

It is shown analytically above that the dynamics of both entanglement and QD under the RWA can be distinguished in the non-Markovian and Markovian regimes, characterized by the ratio . In each regime, the essential feature only depends on the weights of the initial states, and independent of the system-bath coupling strength . In other words, in the framework of the RWA, all essential properties for the physical observable are independent of the system-bath coupling strength, which is not always true with the increasing coupling realized in many recent experiments niemczyk1 (). Hence, it is necessary to go beyond the RWA.

### iii.1 Exact dynamics and comparisons

The dynamics of the pairwise concurrence of the two qubits in their own baths without the RWA can be obtained exactly by the hierarchy approach outlined above. Initiated from the anti-correlated Bell state with , the time evolution of the concurrence from weak, intermediate, to strong coupling cases are presented in the upper panel of Fig. 1. The concurrence under the RWA has been investigated exactly in Ref. bellomo1 () both in the Markovian and non-Markovian regimes. The relevant results for the same parameters are also reproduced in the upper panel of Fig. 1 for comparison.

In the weak system-bath coupling (), the rotating-wave term dominates the entanglement evolution as indicated in Fig. 1(a), and the results including the CRT show negligible deviation from those in RWA. In the non-Markovian regime, the concurrence exhibits periodic oscillations with amplitude damping. The corresponding time of the zero points is . Whether it only vanishes at these discrete moments is very subtle, the periods of the time for zero entanglement is too short to be visible in this case . While in the Markovian regime, the concurrence decays exponentially and vanish asymptotically. The decay rate is larger than the non-Markovian one. So at such a weak coupling, we have not found any evident effect from the CRT, so the previous RWA description is really valid.

What happens if we increase the system-bath coupling to the intermediate regime, say . It is interesting to observe in Fig. 1(b) for that the dynamic behavior is evidently different from that in the RWA. From Eq. (11), one can see that if initiated from anti-correlated Bell state, in the RWA the ESD never happens, which was first found in Ref. bellomo1 () in the same system. But without the RWA, in both the Markovian () and the non-Markovian () regimes, the concurrence exhibits ESD obviously. It follows that the RWA can not describe the essential feature of the evolution of the entanglement beyond the weak coupling regime and is therefore broken down under strong coupling.

Importantly, this physical condition can be realized in the recent strong coupling experiments  niemczyk1 (). In the non-Markovian limit (), the spectral function of baths is reduced to the single mode case by . Hence, the system can be approximately simplified to two independent Rabi models. The corresponding effective system-bath coupling constant is given by , so for , which can be practically realized in the recent experiments of the superconducting qubits coupled to LC resonators, where the coupling constant already reached to ten percentages niemczyk1 ().

In the Markovian regime, the entanglement vanishes quickly and permanently beyond the weak coupling as shown in Fig. 1 (b) and 1(c), mainly due to the generation of extra phonons with increasing coupling, which disturb the pairwise entanglement. No revival of the entanglement appears in this case. Interestingly, the evolution behavior is almost unchanged as the increasing coupling strength for , which follows that is strong enough in this case. However, in the non-Markovian regime, the entanglement dynamics show essential different behavior with the system-bath coupling. The revival of the entanglement after the ESD occurs at , but this revival phenomena never happens in the strong coupling regime, as shown in Fig. 1(c) for . We propose that the more extra phonons activated by the CRT would suppress the feedback from the bath to the qubit chen1 () even in the non-Markovian regime, resulting in the permanent disentanglement. The effective system-bath coupling strength for in Fig. 1(c) can be estimated as , which is hopefully realized in the near future solano1 (); solano2 ().

Then we investigate the dynamics of the QD for two qubits in two independent baths to explore the influence of the CRT. The QD in the RWA has been discussed by Wang et al. wang1 () and Fanchini et al. fanchini1 (). In the lower panel of Fig. 1(d), we display the evolution of the QD for the same parameters as the upper panel by both the hierarchy approach and RWA. In the weak coupling regime, it behaves similar to the concurrence in Fig. 1(a) for both the Markovian and non-Markovian evolutions, where the CRT can be ignored. Beyond the weak coupling regime, as shown in Fig. 1(e) and 1(f), the CRT drives the QD to deviate from those in RWA obviously. Specifically for the non-Markovian case, the QD decays dramatically and shows weak revival signal, which mainly attributes to the excitation of the phonons in baths. For the Markovian case, the QD decreases faster than that in RWA and decay to zero asymptotically.

To see the comprehensive effect of the system-bath interaction, we plot the dynamics of the quantum correlation as a function of both and in Fig. 2 for initial anti-correlated Bell state. Both the entanglement and the discord are suppressed monotonically with the increasing interaction. The revival phenomena for entanglement is not present completely as the coupling strength exceeds a critical value, e.g. is about for the parameters in Fig. 2 (a). The region of ESD becomes wider with the increase of the system-bath coupling. But Fig. 2 (b) reveals that the QD never vanishes suddenly, consistent with the previous observation concluded from an arbitrary Markovian evolution  ferraro1 (). Therefore, such a limitation is removed in the present observation in the exact evolution which definitely includes the non- Markovian effect. No sudden death of the QD is universal in this sense, in sharp contrast with the entanglement. More over, in our exact solution for the full model without the RWA, we actually have not found the real zero QD at any time, despite arbitrary small. This nonzero nature for the QD maybe intrinsic at least for the spin-boson model. Previous observation in the RWA that the QD can vanish at some discreet times wang1 (); fanchini1 () may be an artificial result due to the absence of the CRT.

The essential feature for the dynamics of the quantum entanglement and discord under the RWA is not altered qualitatively with the system-bath coupling in the same non-Markovian or Markovian regime. However, it is not that case in the real exact solution without the RWA. For the weak system-bath interaction, we have shown above that the RWA can basically give the right description of the quantum correlation of two separate spin-boson models for both Markovian and non-Markovian regimes. While in the strong coupling regime, the CRT should be necessarily included for correctly depicting the phonons generation from the baths. And the quantum correlation shows significantly different behavior from that in RWA, especially in non-Markovian case. Hence, we focus on the dynamics of the quantum correlation in the strong coupling and non-Markovian regime in the following with various mostly used initial states.

### iii.2 Non-Markovian dynamics at strong coupling

First, we examine the dependence of the dynamics of the quantum correlation on the value of for the initial Bell state in Eq. (9), which characterize weights of the two superposed anti-correlated spin states. Since the quantum correlation is similar in the regime to that in , the first zone is selected to study the correlation evolution in Fig. 3(a)-3(c). In RWA, both the concurrence and the QD show damping oscillations with stronger amplitude as increases. The other evident feature is that they only vanish at some discrete moments. When the CRT is included, the concurrence exhibits ESD after finite time evolution in the whole regime. The starting time of ESD is much earlier than first vanishing moment in the RWA. The QD first decreases almost to zero and then revives with amplitude much too weaker than that in the RWA, mainly attributed to the suppress from the phonons.

Fig. 3(d)-3(f) presents the results for the another initial Bell state in Eq. (10). It is known that the discord dynamics for is similar to the concurrence for . Hence, we also select these initial states for comparison as in Ref. fanchini1 (). The exact dynamics shows that the ESD occurs very quickly and irrelevant with the value of . The QD first decreases almost to zero, and then revives weakly. The revival amplitude becomes stronger with the value of . For the same , the revival amplitude is much too weaker than that in the RWA, due to the CRT effect, as demonstrated in the inset.

Finally, we assume the extended Werner like state () as the initial state and see what happen to the evolution of these quantum correlations. The main results are list in Fig. 4. As one classical kind of the mixed states, it is important in quantum information processing werner1 (); popescu1 (); hagley1 (). As the initial state with , the concurrence in RWA shows complete death in , whereas the discord in the regime behaves similar with the concurrence in  wang1 (). Hence, we study the quantum correlation in the regime to explore the effect of the CRT in Fig. 4(a)-4(c). The entanglement vanishes suddenly and permanently for all values of , RWA ones show damping oscillations for large (e.g. , ). The QD is also strongly suppressed first, and then revives weakly. The revival amplitude increase slightly with . It is again shown that the QD is not so fragile as the entanglement. If is replaced by in the initial extended Werner like state, no essential different behavior can be observed except that the revival amplitude becomes much too weaker, as demonstrated in Fig. 4(d)-4(f).

For the three kinds of initial states studied above, it is generally find that the larger initial quantum correlations is, the stronger they evolves at the same moments.

## Iv two qubits coupled to one common bosonic bath

The Hamiltonian of two qubits coupled to the common bath reads

 Hs = ω02∑a=A,Bσaz, Hb = ∑kωkb†kbk, Hsb = V∑k(gkbk+g∗kb†k), (13)

where the notations are the same as those in Eq. (2). If the spectrum of the common bath is also Lorentz distribution , and the correlation function of the bath should be also . By choosing the proper auxiliary influence functional in Appendix, the hierarchy equation is given by (Ama1 ()

 ∂ρ→n(t)∂t = −(iH×S+→n⋅→ν)ρ→n(t)−i∑k=1,2V×ρ→n+→ek(t) −iγλ4∑k=1,2nk[V×+(−1)kVo]ρ→n−→ek(t),

where , the unit vector and . The reduced density matrix for qubits is , by which we can derive the quantum correlation numerically.

In RWA, the system-bath interaction becomes . If we represent the Hamiltonian in the basis , is found to be isolate from others. Hence we can rewrite the Hamiltonian in three level formgarraway4 ()

 H0 = 2ω0|2⟩⟨2|+ω0|+⟩⟨+|−ω0, (15) Hint = √2∑kgkak(|+⟩⟨0|+|2⟩⟨+|)+H.c.

For the single particle excitation or double excitation case, the motion equation of the system density matrix can be derived by the called pseudomode mater equation in the interacting picture garraway1 (); garraway2 (); garraway3 (); garraway4 (); mazzola1 (); fanchini1 ()

 ∂ρI(t)∂t = −i[V,ρI(t)] (16) −γ(a†aρI(t)+ρI(t)a†a−2aρI(t)a†),

where is the creator (annihilator) of the bosonic pseudomode,

 V=√λγ(a|+⟩⟨0|+a†|0⟩⟨+|+a|2⟩⟨+|+a†|+⟩⟨2|).

and . By integrating the bosonic part of the density matrix, the system reduced density is given by . For single excitation, an alternatively exact solution has also been done in Ref. maniscalco1 (). If the initial state is chosen by form in Eq. (6), the concurrence is derived by mazzola1 ()

 C(t) = 2max{0,|ρ23|−√ρ11ρ44,|ρ14|−√ρ22ρ33}.

For the common bath, the entanglement evolution initiated from the anti-correlated Bell state without the RWA has been studied by Ma et al. ma1 (). On the other hand, in the previous studies of the dynamics of QD under the RWA, the correlated Bell state is usually selected as the initial state wang1 (); fanchini1 (); garraway4 (); mazzola1 (). Therefore in the present exact study without the RWA, we will focus on the initial correlated Bell state.

First, we exhibit the dynamics of the pairwise concurrence of the two qubits coupled to a common bath with and without the RWA initiated from the Bell state in Eq.  with in the upper panel of Fig. 5 for various coupling strength in both Markovian and non-Markovian regimes. Similarly, for the weak coupling (), the CRTs take little effect for both dynamics in different evolution regimes, showing the same results as in Ref. fanchini1 (). In the intermediate coupling regime (), the ESD occurs more quickly than that in the RWA in the Markovian regime, and the evolution shows slight deviations from the RWA ones with damping oscillations in the non-Markovian regime. When the strong coupling regime is reached, the CRT takes remarkable effects in both evolution regimes, where concurrence vanishes completely and permanently without the RWA. Hence, the CRT also suppress the entanglement dramatically as the system-bath coupling for the common bath, similar to those for the independent baths.

Then we turn to the relevant QD, as shown in the lower panel of Fig. 5. It is interesting to note that in the non-Markovian regime, the CRTs take very limited effect with the increasing coupling. Even at the strong coupling, only slightly deviations are observed. More interestingly, the QD takes finite value at any time, and take a moderate finite value asymptotically in the long time limit, for arbitrary coupling strength. In the Markovian regime, the QD tends to a remarkable value in the intermediate and strong coupling regime. Actually, even at weak coupling like , the asymptotical value of QD takes a considerable non-zero value.

From the numerical iteration of the hierarchy equation in Eq. (IV) with Bell correlated initial state under Markov regime (e.g. in Fig. 5), we find that the form of the two qubits reduced density matrix under the same basis as in Eq. (6) always to be type as

 ρ(t)=⎛⎜ ⎜ ⎜ ⎜⎝a(t)00w(t)0b(t)b(t)00b(t)b(t)0w∗(t)00d(t)⎞⎟ ⎟ ⎟ ⎟⎠. (18)

In the weak system-bath coupling (say ), and decrease gradually to after the finite time evolution, whereas rises to reach . Hence both concurrence and the discord decay to exhibit monotonically as indicated in Fig. 5(a) and 5(d), similar to that in the RWA. When the coupling strength becomes stronger, e.g. , the phonons are excited considerably due to the strong coupling chen1 (). According to the CRT , the qubits are also excited accompanying the extra excitations of the phonons, which stabilize and in long time evolution. It is interesting to observe that they tend to and in the strong coupling limit at . The corresponding reduced density matrix is then described by . which is not altered by different . Under this separable mixed state, the disentanglement occurs naturally from Eq. (IV), and the QD is stabilized at .

To see the overall effect of the system-bath interaction on the dynamics of the quantum correlation beyond the weak coupling regime, we draw the corresponding D plots for in Fig. 6 in the non-Markovian regime. For the entanglement, one can find that there are two regimes. One is damping oscillation at . The other is the regime where the entanglement revives following the ESD, and the revival will be suppressed and finally disappear with the further increase of the coupling regime. Very interestingly, after weak oscillations, the QD is generally robust in the long time evolution in arbitrary coupling regime, which is very useful as the quantum information resource.

## V conclusion

In summary, without performing any approximation, we exactly investigate quantum correlations in two typical kinds of spin-boson models by applying the hierarchy approach. The corresponding RWA results are also reproduced for comparison. The Markovian and non-Markovian effects of two kinds baths in strong system-bath interaction on the dynamics of the quantum correlation are revealed in detail, which demonstrates that the CRT should be necessarily included beyond the weak coupling regime to capture the essential features of quantum correlations correctly.

For two independent baths, the CRT suppresses both the entanglement and the QD dramatically, due to the activations of phonons in baths as the system-bath coupling increases. If initiated from the anti-correlated Bell state, The ESD is driven by extra phonons activated with the increasing system-bath coupling, which is absent in the RWA study. It follows that previous picture under the RWA is essentially modified with the increasing coupling. The QD is found to be always higher than zero, despite extremely small, also different from that in the RWA. The dynamics with other initial conditions is also discussed at the strong coupling and non-Markovian regime. The general trend is that the larger initial quantum correlations, the stronger they evolve.

For the model of two qubits coupled to the common bath, the suppression of the entanglement by the CRT is also observed. It is interesting to find that the non-classical correlation reserved by the QD can be enhanced as the coupling strength, in contrast with the entanglement. It follows that while the entanglement becomes fragile with the increasing coupling to the surrounding environment, the QD is however robust in the long time evolution. In recently experiments, the strong system-bath coupling has been reached in circuit QED systems niemczyk1 (), where for multi modes and the RWA description is broken down. The present study for the evolution of the quantum correlations beyond the RWA may motivate the corresponding experimental studies based on these strong coupling systems.

## Vi Acknowledgement

This work was supported by National Natural Science Foundation of China under Grant No. 11174254 and the National Basic Research Program of China under Grant No. 2009CB929104.

## Appendix A basic concept and application of the Hierarchy equation

The Hamiltonian of the system-bath coupling system is modeled as  . Then the time dependent density matrix is given by . In the interacting picture, the corresponding density matrix reads , with . Under the specific representation ,

 ⟨xt,αt|ρI(t)|x′t,α′t⟩ (19) = ⟨xt,αt|UI(t)ρ(0)U†I(t)|x′t,α′t⟩ = ∫dx0dx′0∫dα0dα′0⟨xt,αt|UI(t)|x0,α0⟩ × ⟨x0,α0|ρ(0)|x′0,α′0⟩⟨x′0α′0|U†I(t)|x′t,α′t⟩.

It is known that . Hence by integrating the bath degree of freedom, the reduced system density is

 ⟨at|ρIS(t)|α′t⟩ (20) = ∫dxt⟨xt,αt|ρI(t)|xt,α′t⟩ = ∫dxtdx0dx′0∫dα0dα′0⟨xt,αt|UI(t)|x0,α0⟩ × ⟨x′0,α′0|U†I(t)|xt,α′t⟩⟨x0,α0|ρ(0)|x′0,α′0⟩.

Considering the propagator in functional form

and the initial factorized state  ( is the equilibrium state), the reduced density matrix is reexpressed by

 ⟨αt|ρIS(t)|α′t⟩ = ∫dα0dα′0⟨α0|ρS(t)|α′0⟩ × ∫αtα0D[α(τ)]∫α′tα′0D[α′(τ)]∫dxtdx0dx′0⟨x0|ρB|x′0⟩ × ∫xtx0D[x(τ)]∫xtx′0D∗[x′(τ)]eiS[x(τ),α(τ)]−iS[x′(τ),α′(τ)].

Then the functional field is introduced by

 F[α(τ),α′(τ)] (22) = ∫dxtdx0dx′0⟨x0|ρb|x′0⟩ × ∫xtx0D[x(τ)]∫xtx′0D∗[x′(τ)]eiS[x(τ),α(τ)]−iS[x′(τ),α′(τ)] = Trb{ρbU†b[α′]Ub[α]},

where . The functional field generally satisfies , where if there is no system-bath interaction.

Here, we include the harmonic bath with the Hamiltonian as . Without loss of generality, we firstly consider only one mode case as , and the system-bath coupling term is . In the interaction picture, , and . Through the tedious by standard Feynman-Vernon influence functional calculation, we derive the influence functional as , where the influence phase is

 Φ[α,α′] = ∫t0dt′∫t′0ds{Q[α(t′)]−Q[α′(t′)]} × {C(t′−s)Q[α(s)]−C∗(t′−s)Q[α′(s)]},

with corresponding bath correlation function

 C(τ)=1mΩ[coth(βΩ/2)cosΩτ−isinΩτ]. (24)

Here we should note that the variable is a real number. Then we generalize the results to the multi-mode case, where the system-bath interaction now is , and . The influence functional becomes

 F[α(τ),α′(τ)] = exp(−∑a∫t0dτ{Qa[α(τ)]−Qa[α′(τ)]} (25) ×{~Qa[α(τ)]−~Q∗a[α′(τ)]}),

where

 ~Qa[α(t)]=∑a′∫t0dτCaa′(t−τ)Qa′[α(τ)]. (26)

The correlation function in independent baths is

 Caa′(τ) = ⟨Fa(τ)Fa′(0)⟩ = ∑kδaa′cak22makΩak[coth(βaΩak2)cos(Ωakτ)−isin(Ωakτ)].

By using the spectral function , we finally have

 Caa′(t)=δaa′∫dωJa(ω)[coth(βaω2)cos(ωτ)−isin(ωτ)].

In the schrodinger picture, the motion equation of the system density is

 ∂ρS(t)∂t=−i[HS,ρS(t)]+e−iHSt∂ρIS(t)∂teiHSt. (28)

From Eq. (A) and (22), we have

 ∂ρIS(t)∂t = ∫dα0dα′0⟨α0|ρS(0)|α′0⟩ × ∫αtα0D[α(τ)]∫α′tα′0D∗[α′(τ)]∂F[α,α′]∂t.

The key step here is to connect time differentiation of the density matrix with that of the influence functional. The proper way is employing the hierarchy equation by introducing the auxiliary influence functional

 ρIS(n,t) = ∫dα0dα′0⟨α0|ρS(0)|α′0⟩ × ∫αtα0D[α(τ)]∫α′tα′0D∗[α′(τ)]Fn[α,α′].

While such auxiliary choice is closely dependent on the specific form of the bath correlation function. The hierarchy equation used in the present paper is described in the following.

For two independent baths with Lorentz type spectral function , the correlation function simply is

 Caa(ω)=λaγa2exp[−(γa+iωa0)t]. (31)

according to the following formula

 ∂∂t~Qa[α(t)]=λaγa2Qa[α(t)]−(γa+iωa0)~Qa[α(t)],

we introduce the auxiliary influence functional of each subsystem as

 Fa,→k[α,α′]=(~Qa[α(t)])k1(~Q∗a[α′(t)])k2F[α,α′], (32)

where and index . Then the whole influence functional becomes

 F→n,→m[α,α′]=FA,→n[α,α′]FB,→m[α,α′]. (33)

Finally, the evolution of the density matrix reads

 ∂∂tρ→n,→m(t) = −(iH×S+→n⋅→νA+→m⋅→νB)ρ→n,→m(t) +∑k=1,2(−1)kQ×Aρ→n+→ek,→m(t) +λAγA4∑k=1,2nk[QoA+(−1)k+1Q×A]ρ→n−→ek,→m(t) +∑k=1,2(−1)kQ×Bρ→n,→m+→ek(t) +λBγB4∑k=1,2m