Bruno Zambon

Roberto Franzosi
email: franzosi@fi.infn.it

C N I S M Udr

Dynamics of atoms in a time orbiting potential trap: consequences of the classical description

Keywords: numbers: 03.75.Kk, 05.45.-a, 03.65.Sq

The classical model that describes the motion of an atom in a magnetic trap is solved in order to investigate the relationship between the failure of the usual adiabatic approximation assumption and the physical parameters of the trap. This allows to evaluate the effect that reversing of the bias field rotation produces on the vertical position of the atomic orbit, a displacement that is closely related to the adiabatic character of the trap motion. The present investigation has been motivated by a similar experimental test previously carried out in the actual magnetic time orbiting potential trap. We find that the non-adiabatic effects provided by the classical model are extremely small. Thus, we conclude that the theoretical explanation of the experimental measures, requires a quantum description of the dynamics in magnetic traps.

I. INTRODUCTION

In the realization of Bose-Einstein condensates (BECs) in alkali atoms, the gas is confined and cooled in magnetic traps. Among the different kinds of setups, in many experiments have been used time-orbiting-potential (TOP) traps [1], which are realized with two magnetic components: a static inhomogeneous magnetic-field and an uniform rotating one. The static inhomogeneous component originates a potential well for the atoms, around a local minimum of the magnetic field. The smaller this minimum is, the tighter the atomic confinement is. Nevertheless, at zero magnetic field, atoms can leak from the trapping potential as a consequence of Majorana transitions. In order to avoid these atomic losses, in the TOP traps is added the second magnetic component, that eliminates the zero of the magnetic field at the trap center. The tight and stable confinement of the atomic clouds, achieved in the TOP traps, makes these systems well-suited for the Bose-Einstein condensation as well as for quantum statistical/dynamical studies. With respect to the theoretical investigation of the dynamics of atoms in TOP traps, the adiabatic approximation, in which the atomic magnetic moment is assumed to be constantly anti-aligned with respect to the magnetic field, is usually assumed. A further possible mathematical simplification, comes from the timeaverage approximation, in which one assumes as potential for the atoms the time average of the true one. Under these two approximations, the potential seen by an atom in a TOP trap results harmonic and then the system Hamiltonian has trivial quantum eigenstates.

However, such systems display a residual micromotion due to the fact that the trapping potential is, indeed, time-dependent. This phenomenon has been theoretically predicted in Ref. [2] on the base of a classical model that prescinds from the adiabatic approximation, and it has been experimentally observed in Refs. [3,4]. These latter papers, have recently met the interest of several authors [5][6][7][8][9][10]. Other theoretical investigations on the quantum dynamics of atoms (condensates) in TOP traps with the adiabatic approximation [11,12] ( [13]), or without this approximation [14], agree with the experimental observations about this phenomenon.

A second surprising phenomenon has been observed in Ref. [3], in this experiment, and for the first time, this team has been able to reveal an anomalous shift in the vertical position of the atomic cloud, as the rotating uniform magnetic field reverses its rotation. In the same paper it has been hypothesized that the failure of the adiabatic approximation could be at the origin of this phenomenon. To our knowledge, this anomalous shift has never been explained from a theoretical point of view. Thus, the main purpose of the present paper is to establish if, within a classical description of the dynamics of an atom in a TOP trap, such phenomenon can be explained. The definite outcome of the present paper is that the non-adiabatic effects in the classical description of this system are extremely small. Thus, such anomalous shift, cannot be explained within the classical scenario commonly used to give a simple explanation of the mechanism behind cloud stability in the magnetic traps.

Consequently, we conclude that, the theoretical explanation of this phenomenon, requires a quantum description of the dynamics of atoms in a TOP trap. This work is organized as follows. In section II and III we review some results and general ideas concerning the atomic motion in a magnetic trap. In particular in section II A are introduced the classical equations of motion describing the dynamics of an atom in a magnetic trap, whereas section II B is devoted to a general discussion about the adiabatic approximation that simplifies drastically the resolution of the equations of motion for this kind of systems. In section III, we derive several periodic orbits that take place in a magnetic trap that has a cylindric symmetry around the vertical axis. In section IV we discuss the periodic motion in a triaxial TOP trap similar to that used in the experimental investigation of Ref. [3] and we put a special emphasis on the reversing bias field rotation test. Therefore, in Sec. IV A we derive the periodic orbit that takes place in a triaxial trap for a very small asymmetry. In Sec. IV B we perform a perturbative calculation that, by using the result of App. A, allows us to evaluate the effect that the asymmetry produces on the dipole orientation in space and the consequent reposition of the atomic trajectory along the vertical axis. Finally, in Sec IV B, we evaluate the consequent shift in the vertical position of the atomic trajectory as the bias field reverses its rotation. In Section V we review some fundamental point concerning the above experimental results, and the connected numerical simulation. Thus we analyze the perspective emerging from the fact that the results provided by the classical model do not agree with what is observed experimentally.

II. ATOMIC MOTION IN A MAGNETIC TRAP

The confinement in space of a neutral particle carrying a magnetic moment, can be obtained by means of the gradient magnetic-field forces experienced by a magnetic dipole, in the presence of an inhomogeneous magnetic field. This mechanism is at the basis of the atomic traps used by many experimental groups working on the subject of Bose-Einstein condensation [1]. In these traps, the atoms are first collected by means of laser forces in a limited region of the space where such inhomogeneous magnetic field plays its active role, and, afterwards, cooled down by an evaporative mechanism so as to make the transition to a Bose-Einstein condensate possible. So far the physics of the confinement by magnetic field has been explained mainly by means of a set of equations in which the particle is seen as a classical point-like magnetic dipole that obeys to the classical equations of motion [2]. Some quantum mechanical calculations exists for the more simple configuration without the rotating bias field [15]. Nevertheless, this rotating field plays an important role in the TOP traps and makes the physics of these systems more intriguing and complicated. In fact, an important matter inherent this motion is the validity of the adiabatic approximation.

In the case where this approximation holds, the magnetic dipole is constantly anti-aligned with respect to the magnetic field and, thus, the magnetic energy of the atoms assumes the role of a space confining potential. Existing quantum mechanical approaches [11][12][13] of atoms (or Bose-Einstein condensates) within a magnetic trap, assumes from the start the fulfillment of the adiabatic condition which is, actually, classical in nature. Indeed the question concerning its range of validity has been tested experimentally in [3]. From these investigations it turns out that such an approximation could become less stringent in the proximity of the marginal stability configuration of a TOP trap. The critical test performed in [3] concerns the observation of an anomalous shift in the vertical position of the atomic cloud as the bias radio-frequency field reverses its rotation. The issue of adiabatic motion has also been theoretically investigated for a symmetric trap, from which the triaxial trap used in [3] is obtained by setting the symmetry axis in the rotating bias field horizontal plane.

The dynamics in the cylindric trap has been studied by means of a variational technique approach to the full quantum problem in [14]. Whereas a semi-classical equations including quantum atomic correlation has been used to investigate the more general problem of the motion in a inhomogeneous magnetic field in [16].

However, as said above, genuine quantum approaches to the study of the dynamics of atoms in the TOP trap of Ref. [3] does not exists. These would have the merit of shedding light on the condition of validity of the adiabatic approximation in the presence of the spatial motion that a magnetic dipole confined in a magnetic trap undergoes. These questions can also be related to the existence of a Berry's phases phenomenology [17]. While acknowledging the need of a more detailed quantum investigation we will discuss here some features of the magnetic confinement in a magnetic trap entirely from the classical point of view and in particular we will focus on the test of the reversal of the bias field rotation. We believe that this has to be carried out considering the importance and the intuitive support that this model provides in explaining the trapping mechanism in the BEC experiments. Thus before of proceding with a quantum computation it is necessary to assess clearly the scope of validity of the already existing classical-model. This will be done here by checking the validity of this model against the experimental results in [3]. In particular when we will refer generically to the experimental work or to the experimental measurements we will mean experimental work and experimental measurements as done in Ref. [3].

A. Trap classical equations

As a basic ingredient of a magnetic trap there is a magnetic field configuration having a minimum of its intensity at a point in the space generally named trap center. This is obtained with a static quadrupole magnetic field whose components along the Cartesian (

Moreover, in order to prevent the spin flip transition that would occur when the atom approaches the trap center, a rotating radio-frequency uniform magnetic-field, usually known as bias field, is superimposed to the static inhomogeneous one. In such a way the atom trajectory is kept quite far away from this point. Let B b be the bias field intensity and ω its angular velocity in the x -y plane. A different configuration for the inhomogeneous field can be obtained by setting the symmetry axis to coincide with the x axis and by inverting the currents, this leads to

b x = 2 b x , b y = -b y , b z = -b z , (2)
which also provide a stable confinement. We denote by B the total magnetic field experienced by an atom which is given by the sum of the inhomogeneous magnetic field and the bias rotating field. The configuration where the inhomogeneous magnetic field is given by Eqs.

(1), is known as symmetric or cylindric TOP trap, while the one corresponding to the configuration of Eqs. ( 2) is named triaxial trap and is the one used in the experimental setup of Ref. [3].

Within the classical realm, the dynamics for an atom in a magnetic trap is given by the equations of motion

m ¨ r = F = ( µ • ∇) B + m g , ˙ µ = γ µ × B , (3)
where γ connects the magnetic moment to the angular momentum of the particle, i.e. µ = g e 2me L = γ L := µ n, here m e is the electron mass, g is the gyro-magnetic factor and n is the unit vector parallel to the magnetic moment. In the discussion that follows we will use normalized units for the quantities involved in the trap dynamics. First of all we will assign a unit value to the period of the rotating bias field and we adopt a length unit which coincide with the characteristic length of the trap. Furthermore, in order to lift the constant in the second of Eqs. (3), we will include the γ constant in the magnetic field B, that will be measured in angular frequency units. Consistently we will absorb the factor 1/γ in constant µ/m. For a typical TOP trap, we have a bias field frequency of the order of 10kHz and a characteristic length of the order of 1.µm. Thus the time unit will correspond to T 0 = 10 -4 s, the length unit to L 0 = 10 -6 m, and the magnetic field intensity in the new unit, will be obtained by multiplying the old value for the factor γT 0 . E.g., a magnetic field of 10.Gauss (order of magnitude of the bias field used in [3]) transforms to a field of about 10 -3 1.6×10 - 19 9.1×10 -31 10 -4 = 1.8 × 10 4 , a field gradient b of 100.Gauss cm -1 transform to

10 -2 ×1.8×10 11 ×10 -4 10 -2 ×10 6
= 18. The parameter µ/m of the first of Eqs. ( 4), in the new unit is obtained by multiplying the original value for the factor T 0 /(L 2 0 γ). Thus, by doing so for a trap, as the one of Ref. [3], operating with 87 Rb atoms we have the following numerical values µ m = 0.036, which is independent from the unit chosen, g = 0.0981, ω = 2π, B b ≈ 10 4 and a value of the field gradient b ranging from 0.1 up to a value of 10. We will refer to these values in the following as the the standard trap parameters. From Eqs. (3) we finally obtain the normalized evolution equations

¨ r = µ m ( n • ∇) B + g , ˙ n = n × B , (4)
where n is the unit vector corresponding to the atomic magnetic dipole.

B. Adiabatic approximation

In the present section we summarize the conditions of validity of the adiabatic approximation and its consequences about the dynamics of atoms in magnetic traps. The adiabatic approximation, that greatly simplifies the analysis of the atomic motion, comes, pretty naturally, from an analysis of equations ( 4) and, in particular, from the second one. By indicating with n the projection of the vector n on the magnetic field B, i.e. n = n • B B , and with n ⊥ the component of n perpendicular to this same field, that is

n ⊥ = n -( n • B B ) B B
, it is easy to show, by using the second of Eqs. ( 4), that the following equality holds true

ṅ = n ⊥ • ˙ B B . (5)
From this equality it results clear that the validity of the adiabatic regime is a consequence of the fast rotation of n ⊥ around B, that takes place at an angular velocity of the order of

B. In fact, if ˙ B
B has the same direction of B, i.e. B changes only in magnitude, the quantity n is exactly constant. If this is not the case, but ˙ B is small, ṅ is also a small quantity with zero average because of the fast rotation of n ⊥ around B. Thus, under the hypothesis of fast rotation of n ⊥ (that is B >> 1) and small ˙ B, we can conclude that n is a quasi-constant of motion and the magnetic energy of a particle in the trap results to be

U = -µ • B = -n µ | B| . (6)
Now the magnetic energy only depends on the projection of the magnetic moment on the field B and not from the orientation of the dipole in space, since the other degrees of freedom have been frozen as a result of the adiabatic approximation. Thus, a particle anti-aligned with the field, or more in general a particle for which n is negative, experiences a confining magnetic potential in the proximity of a minimum of the magnetic field intensity.

In the case of the TOP traps, where the total magnetic field B is time-dependent, it is more convenient to examine the torque equation (second of ( 4)) in a reference frame rotating around the z axis together the bias magnetic field. In this frame, the effective magnetic field experienced by the particle is given by B 0 = B + ω, where ω is directed along the z axis. In the case of a cylindric trap the latter will be a time-independent quantity, whereas for the triaxial trap B 0 remains time-dependent. Thus, in the case of the TOP traps, the adiabatic approximation has to be applied in the reference frame rotating with the bias field and, therefore, it is the field B 0 that plays the role of the field B of Eqs. ( 5) and (6). Some questions related to the fact that a proper choice of the reference frame really improves the analysis of a dynamical system when the adiabatic approximation is to be used have been discussed in [18]. As a matter of fact it is out of doubt that in the TOP traps the frame co-rotating with the bias field is the one that has such property.

In the following we will summarize briefly the peculiarities of the motion in the symmetric and triaxial traps. In this last case, we will calculate the shift of the vertical position of the atomic cloud, as the bias field reverses its rotation. In the case of the symmetric trap this shift is shown to be exactly ∆z = -ω b . In fact, a -2ω change in the magnetic field B 0 , is exactly compensated by the z component of the inhomogeneous field in Eq. ( 1), as the orbit displaces itself of ∆z in the vertical position. However, in the case of the triaxial trap, this argument does no more hold true because its inhomogeneous field configuration does not enjoy the rotational symmetry around the z axis. This produces a different dynamical evolution of the atomic magnetic moment, that influences also the trap vertical position beyond the simple argument given above. The calculation of this contribution, will show how this shift is connected to the failure of the adiabatic approximation, and to which extent this mechanism can be understood within the limit of this classical model.

III. CYLINDRIC TRAP MOTION

A. The periodic orbits

Stationary solutions of the classical equations (4) for the case of a cylindrical trap, have been given in Ref. [2]. In order to derive these solutions, and to find the stable orbits, Eqs.

(4) are transformed, and then solved in the reference frame rotating around the z vertical axis which is opposite to the gravity vector g. This frame is the most convenient one because here the bias field results to be time-independent as well as the fictitious magnetic field B 0 = B + ω. In the latter expression ω is directed along the z axis, and B includes both the bias field and the inhomogeneous field. Remarkably, this latter has the cylindrical symmetry around the z axis. For the first of Eqs. (4), also the apparent forces must be taken into account. There exists a stable solution in which the magnetic dipole is anti-parallel to field B 0 , this complying with the torque equation, and the particle is at rest in this frame. Thus the total force on the particle must be zero. The projection of the equilibrium equation along the z axis yield F z = -2 µ b n z -m g = 0, from which we derive

n z = -cos ϑ = - m g 2 b µ . (7)
Let us consider now the equilibrium in the x-y plane. To fix the ideas we chose the rotating x axis coincident with the bias field (which has intensity B b ), thus B 0 = (B b +bx, by, -2bz +ω)

and the gradient force has components in the x-y plane given by F x = bµn x and F y = bµn y .

The balance with centrifugal force leads to n x = -mω 2 bµ x and n y = -mω 2 bµ y. By combining this latter equation with the anti-alignment condition, that is

B 0 × n = 0, from ( B 0 × n) z = 0
we obtain y = 0. This means that the particle rotates in phase with the bias field. The above arguments lead, therefore, to

n x = -1 -n 2 z = -sin ϑ = - m ω 2 b µ r , (8)
where n = (-sin ϑ, 0,cos ϑ) and r being the distance of the equilibrium point from the origin, i.e. the radius of the particle orbit in the original frame at rest. This latter quantity can be derived by the condition of normalization of n, and it results

r = µb mω 2 1 - mg 2µb 2 . (9)
The magnetic moment configuration for the dynamical regime of this trap, and its relative position with respect to the bias field and total field are shown in Fig. 1. Finally the anti-alignment condition between n and B 0 leads to

B 0z B 0x = -2 b z + ω B b + b r = 1 tan ϑ . ( 10
)
By combining this latter equation and the one derived from the Eqs. ( 7) and ( 8), we can obtain the z equilibrium position. From Eqs. ( 9) and ( 10) we can observe that if ω changes its sign, by keeping fixed all the other parameters of the trap, ϑ does not undergo any change, as can be gather from Eqs. ( 7) and ( 8), and a displacement in the z position of amplitude

ω b will therefore occur.

B. The pseudo-periodic orbits

We derive in the present subsection a particular kind of orbits that are important for the discussion we will make in the following about the results contained in Ref. [3]. In the case in which the bias field B b is much greater of its rotation frequency ω, a class of orbits which are almost periodic do exist, even if the magnetic moment is not exactly anti-aligned with total field B 0 . As we have recalled in the discussion about the adiabatic approximation, whatever the initial angle between the dipole and the field is, it will remain unchanged since the dipole performes a very fast processional motion around the magnetic field. As far as the orbital motion is concerned, the fast oscillation of the dipole, that happens on shorter time scales, will be averaged out, and, consequently, the force on the dipole will depend only on its time averaged value. When a motion similar to that described above will take place, the average vector n replaces the original vector n, in the expressions ( 7)- (10) for the periodic motion given in the previous section. Thus, n will be anti-aligned with respect to B 0 , and will have intensity | n | = | cos ψ| where ψ is the initial angle between the dipole and the field. Thus, from Eqs. ( 7), ( 8) and ( 10) we deduce

n z = -| cos ψ| cos ϑ = - mg 2bµ , n x = -| cos ψ| sin ϑ = - mω 2 bµ r , (11)
that combined with the equation

-2bz + ω B b + br = 1 tan ϑ ( 12
)
determine the orbits height. It is worth to pointing out that, in this case, the orbit z coordinate depends also on the dipole-field angle ψ as also ϑ does, in virtue of Eqs. (11). In other words this orbit z-coordinate depends on the dipole initial conditions.

One may ask if, such an orbit, has any relation with the trapping mechanism operating in the BEC experiments. In this respect it should be remarked that the only classical orbits which correspond to a stationary quantum state are periodic orbits. For these trajectories, the wave associated to the orbital motion undergoes constructive interference as it returns in the same position in space, after that a period has elapsed. This condition gives account for the old Sommerfeld's quantization rule. On the contrary, in the quasi-periodic orbit, the magnetic dipole evolution, is quite independent of the orbital evolution and the particle returns in the same spatial position with its dipole state differing from that of the previous passage through that point. This would lead to destructive interference of the wave associated to the particle. Therefore, the pseudo-periodic orbits can hardly be associated to any stationary quantum state of the system.

C. Oscillations in the vertical position

In the presence of a strong adiabatic field, we can easily describe the oscillations which take place when the initial z-position of the trapped particle is different from the z position given by Eq. ( 12), in the case of a pseudo-periodic trajectory with initial angle ψ. Indeed, the potential energy in the case of strong adiabatic motion is given by

U = mgz + µ| cos ψ| (B b + br) 2 + (-2bz + ω) 2 (13)
to which it corresponds the force

F z = -mg -2bµ n z = -mg + 2bµ| cos ψ| (-2bz + ω) (B b + br) 2 + (-2bz + ω) 2 . ( 14
)
Thus, when r is small such as br can be neglected with respect to B b , this force depends only on z, and the resulting motion is an oscillatory one around the z equilibrium position derived by Eq. ( 12). This oscillating regime is important for the discussion that will follow. In fact the z equilibrium position of a BEC cloud is experimentally evaluated just by measuring the average value of these oscillations.

IV. TRIAXIAL TRAP MOTION

In the case of a triaxial trap it is no more possible to obtain an exact solution for the periodic orbit as done above. Indeed, the rotational symmetry about the vertical axis is lost and, consequently, the magnetic field B 0 in the frame rotating with the bias field, remains a time-dependent quantity. In order to overcome this problem, we will perform a perturbative calculation in which we will take as a basis the periodic orbit that takes place in this trap, for a very small b parameter. Therefore, in Sec. IV A, we will derive such an orbit, then, by using the result of calculation of App. A, we will evaluate the effect that the small inhomogeneous static magnetic field, seen as a varying one by the dipole during its motion, produces on the dipole orientation in space. We shall show that in general, the dipole orientation, after that the particle returns in its initial position in space, will be changed.

Thus, an additional rotation that brings the magnetic dipole back to its initial state, can be obtained by a static magnetic field. This magnetic field compensation will be obtained by a vertical shift of the orbit in space. Of course, such an adjustment is not needed for the symmetric trap, because the inhomogeneous field is seen as constant one during the orbital motion. It is noteworthy that, such a repositioning of the orbit, will turn out to be different for a bias field rotating in the positive or in the negative sense. This latter point will be shown in Sec. IV B

A. Zero order solution

In this section, we derive the periodic orbit that takes place for a very small b parameter.

In order to accomplish this, we neglect the inhomogeneous field components in the x -y plane, only in the torque equation. Thus we write

˙ n = n × (B b cos(ωt), B b sin(ωt), -bz) .
We can assume the adiabatic approximation i.e that n is anti-aligned with respect to the approximate magnetic field, which is given by the bias field and by the z component of the inhomogeneous field. This condition allows us to write, in the laboratory frame, n x = -n ⊥ cos ω t and n y = -n ⊥ sin ω t, where n ⊥ indicates the intensity of the component of the dipole unit vector in the plane x -y. Let n z be the intensity of the z component of the dipole unit vector. Thus, the anti-alignment condition gives the equation

n ⊥ n z = B b (bz -ω) . (15)
The forces on the particle are

F x = -2 µ b n ⊥ cos(ωt) , F y = µ b n ⊥ sin(ωt) , F z = -µ b n z -m g , (16)
by solving the Newton's equations in the x -y plain, we get the time law for the position coordinates x and y in the laboratory frame

x(t) = 2 µ b m ω 2 n ⊥ cos(ωt) , and
y(t) = - µ b m ω 2 n ⊥ sin(ωt) , (17)
where n ⊥ is still an unknown quantity. Thus, for small b the orbit takes place entirely in a plane parallel to the x -y one on an elliptic trajectory, and is counter-rotating with respect to the bias field. The force equation along the z axis provides the angle by which n must be tilted with respect to this axis, and the value of n ⊥ . Indeed, the equation F z = 0 gives, from the third of Eqs. ( 16),

-n z = m g µ b := cos ϑ , and n ⊥ = sin ϑ .

By combining the latter equations with the one in (15), we get

-b z + ω B b = 1 tan ϑ , (19)
which with (18) makes it possible to determines the orbit's z level.

In order to include in the torque equation the contribution given by inhomogeneous field in the plane x -y, we have to compute the magnetic field b (t), in the x -y plane, that the dipole experiences along its trajectory for effect of the inhomogeneous field. In the laboratory frame this is given by b (t) = (2bx(t), -by(t)), where x(t) and y(t) are given in Eq. (17). It is convenient to refer this field respect to the system frame rotating with the bias field, it is not difficult to find

b x = 3 2 µ b 2 m ω 2 n ⊥ cos 2 ω t + 5 2 µ b 2 m ω 2 n ⊥ , b y = - 3 2 µ b 2 m ω 2 n ⊥ sin 2ω t . (20)
Thus, the inhomogeneous field component on the dipole results to be the sum of a constant field and a small time varying component. Apart form the constant term which resets the value of the bias field, as occurred also in the previous case concerning the symmetric trap, the sinusoidal components produces deviations from the anti-alignment condition. In order to evaluate these deviations, and how they can be compensated, we will solve the torque equation under the presence of a constant field complemented with two small oscillating components.

B. Difference of height between ω and -ω

In the present section we will apply the results found in appendix IV A, making explicit Eq. (A6), in the case of the motion in the triaxial trap. Before making these equations explicit for the case of the triaxial trap, we evaluate the order of magnitude of the deviation from the adiabatic approximation, a fair indication of this being the angle between n 0 and B 0 . From Eq. ( A6) is clear that this quantity increases as B 0 gets closer to ω. However, in the usual experimental setups, this condition is carefully avoided because it is known to lead to instabilities in the trap. Indeed, reasonable experimental values for the quantities that appear in (A6) are B 0 ≈ B b of the order of 10 4 and ω = 2π in our units. Thus, if we assume the condition in which ω B 0 , the second term of Eq. (A6) dominates over the first one. Its intensity is of the order of , and then, it results to be smaller. However, its dependence on ω is important and it may perhaps become more relevant for trap where the adiabatic condition is not so strongly enforced, as in the case of the experimental setup being considered here.

We will derive explicitly the shift entailed from this term in the following of the section. It must be noticed also that, given the relative size of the two terms in the above equation, the adiabatic approximation is broken earlier than the effect of the reversed rotation can be observed.

In order to use the result of appendix IV A, we adopt a frame rotating with the bias field and we set

b 0 = 3 2 µ b 2 m ω 2 n ⊥ . (21)
The only components of the fields b 1 and b 2 of Eq. (A1), different from zero in this frame are b 1,x = b 0 and b 2,y = -b 0 , as they can be read from Eq. ( 20). Moreover we have B 0 = (B b + 5 3 b 0 , 0, -b z + ω) and n 0 = (-sin ϑ, 0.,cos ϑ), where the angle ϑ is that one corresponding to the unperturbed motion and given by Eq. ( 18). Since initially n 0 and B 0 are opposite to each other, Eq. (A6) prescribes how one of them should be modified. However since n 0 , and consequently the angle ϑ is already determined by the equilibrium equation along the z axis, only B 0 can be modified. This modification amounts to a small time-independent contribution obtained, in the specific case of the triaxial trap, by a shift of the whole orbit in the z direction. In order to calculate the angle between n 0 and B 0 we note that the only component of n 0 × B 0 different form zero is given by

( n 0 × B 0 ) y = 1 2 b 0 ∆ 2 (-2 ω sin ϑ + B 0 sin ϑ cos ϑ) , (22)
where one should remind that the angular frequency of the oscillating field component is 2ω instead of ω (Cf. Eqs. ( 20) and (A1)). The above equation makes evident that n 0 still lies in the x -z plane, and that the tilting angle of B 0 with respect to n 0 , whose sign is determined by the y axis, is now given by ∆ϑ = -

1 2 b 0 ∆ 2 - 2 ω B 0 sin ϑ + sin ϑ cos ϑ , (23)
where B 0 is the intensity of the total magnetic fields and ∆ 2 = B 2 0 -4 ω 2 . The new component B 0,z must now satisfy the condition

B 0 ,z B b + 5 2 b 0 = 1 tan (ϑ+∆ϑ)
. By using the approximation

1 tan (ϑ+∆ϑ) = 1 tan ϑ -∆ϑ sin ϑ 2 we obtain z = z 0 + 1 2 b 0 ∆ 2 z 0 + ω b , (24)
where z 0 is the z position of the orbit in the zero approximation. From this expression it results that the shift in z, caused by the change ω → -ω, for which ∆z 0 = -2ω b and ∆ω = -2ω, is given by ∆z

= - 2 ω b 1 + b 0 ∆ 2 . ( 25
)
This expression depends on b but, considered the small size of the perturbing field b 0 =

3 2 µ b 2
m ω 2 sin ϑ, it remains very small anyway even if B 0 becomes close to the bias field frequency ω where, on the other hand, the perturbative calculation presented here would no longer apply. It is important to point out that the correction found above depends, quite naturally, on the ratio between the strength of the average inhomogeneous field experienced by the particle during its motion and the total static field B 0 . Then, it is approximatively given by the ration between b 0 , and B b . The predicted shift is shown in Fig. 2 two is really a very small one even for the parameter used in a real trap of Fig. 2.

C. Non-adiabatic regime and bi-stability of the trap

In order to compare the above theoretical derivation with the results obtained by numerically integrating the equations of motion (4), we have determined, by numeric computations, the periodic orbits for this system. The method, which is customarily used in investigations of this kind, assumes that a periodic, and then closed, orbit, named reference orbit, is known for certain values of the system parameters. A small change of one of these parameters opens in space, by a slight amount, the reference orbit when evolved over a period of the timedependent force. Thus, a small variation of the initial condition could recover the periodicity, and close the orbit in space again. A linear map connecting initial and final values, for all orbits that lay close to the reference one, can thus be determined by performing evolutions of the system with initial conditions very close to those of the reference orbit. This map does not appreciably change as the parameters undergo small changes. Thus the new periodic orbit can be found as the fixed point of this map by inverting numerically the related matrix.

In our case the initial reference orbit is that corresponding to a small value of the parameter In order to make more evident this comparison the zero order shift ∆z 0 = 2ω b has been subtracted to both of these two quantities. At small values of the 1 b the regularity of this dependence is somewhat broken as shown by the little terminating wiggle. This corresponds also to orbits which are on the point of becoming unstable. The parameters used for this figure are B r = 50., µ m =0.036, g = 0.1.

Fig. 3, we observe that the agreement is quite good even for values of b in a more wide range of those commonly encountered in experimental setups. The calculated shift is extremely small such as being, at the present, out of the reach of experimental verification, and certainly far away from the values observed in the experiment of Ref. [3].

By using the numerical routines that we have designed for the analysis of the periodic orbits, we have performed a coarse examination of the unstable dynamical regimes of the trap. It is found that signatures of instability appears, quite naturally, as the intensity of the inhomogeneous field experienced by the particle becomes of the same order than that of the bias field. For these values the peculiar scenario leading to unstable behavior with orbit bifurcations, orbit bi-stability and eventually chaotic behavior will emerge. However the instability is always connected, at least for the range of the parameters observed here, to the loss of the system adiabatic condition. We have numerically determined the angle between the effective field B 0 and the dipole vector at time t = 0, i.e. when the bias field crosses the positive x axis. This is shown in Fig. 4 for both directions of bias field rotation. As it can be seen here the adiabatic approximation is quite resilient even if the trap parameters are rather loose, in terms of the adiabatic condition enforcement, with respect to those normally used in a TOP trap. For the same range of values of the b parameter where an appreciable dis-alignment occurs we have observed a bi-stability regime. The two orbits involved in this regime are shown in Fig. 5. The analysis of the unstable regimes of this system would require a systematic study that is out of the purpose of the present work.

V. COMPARISON WITH THE EXPERIMENTAL RESULTS

There is an evident discrepancy between the results obtained here and the experimental ones reported in Ref. [3]. One of the hypothesis that we can follow to deal with this situation is to assume that a classical model is not an adequate one to describe the dynamics of an atom in the triaxial TOP trap and that a quantum calculation is in order. Before to assume the above hypothesis as the final one, let us to try an alternative analysis that take in account of a procedure for the bias-field reversal, similar to that one actually performed in the experiment of [3]. In the reasoning which follow we will refer to a symmetric trap considered that the perturbation to the magnetic dipole motion produced by the inhomogeneous field of the triaxial trap are really negligible for the experimental parameters. Thus the effect that we are going to describe will occur indifferently for the symmetric and triaxial trap.

Let us suppose that the system moves actually on a quasi-periodic orbit with an angle ψ between dipole and field. We remind that this angle does determine the z position of the orbit according to Eq. ( 11) and Eq. ( 12). Now, if the sense of rotation of the bias field is abruptly changed the system will experience the oscillations along the z axis discussed in Section III C, because the z coordinate corresponding to the new pseudo-periodic orbit will be different from the one in which the system is actually in. In the experiment, the averaging of these oscillations is assumed as the z coordinate of the BEC stationary state corresponding to the new bias field rotation sense. On the other hands, the classical equations will lead to a shift that coincide with a zero order compensation shift only if the angle ψ, and by consequence the equilibrium ϑ, does not change. On the contrary, if this does not occur an anomalous contribution should be added to the standard compensation shift. Indeed, this would be the case if we assume that during the switching between the two rotation senses the dipole orientation in space does not change. If this is the case, the angle ψ will change uniquely as a result of the change of the direction of B 0 in space, as a consequence of the change of its z projection by an amount of ∆ B 0,z = ∆ω = -2ω. In appendix B, we have calculated the shift obtained under these conditions, assuming that both the magnetic dipole and the field lies in the plane x -y, and that ψ is the initial angle between the two. This latter quantity can, of course, assume positive or negative values. After a bit of algebra, reported in appendix B we have found

∆z = - 2ω b 1 + mg 2bµ 1 sin ϑ sin ψ cos 2 ψ . (26)
This expression evidences how by choosing appropriately the angle ψ, the trend of the experimental data, according to which the deviation increases as b approaches its minimum value, could be reproduced. We notice that in this case, also ϑ approaches the zero value, thus amplifying even further the expected effect. The above expression indicates that the anomalous shift effect disappear in absence of gravity, a particularity which coincides also with the conclusions reported in the experimental work. Here we are presuming that the switching between positive and negative ω does occur so fast that the dipole has no time to adapt to the new field B 0 , present after the switch has taken place. Although this is rather unlikely to occur, an alternative mechanism wouldn't make sense. In fact, it would be possible to imagine that a loss of the adiabaticity could take place during the switching process, if at certain point the intensity of the magnetic field would get close to zero, in such a way that the angle between the field and the dipole emerges out the switching process changed. This, of course, would depend on how technically this switching is carried out.

If this would be the case, the new oscillation point would be displaced by a quantity that contains a further anomalous contribution. However, if it would be so, the anomalous ∆z a would be only a fictitious one being dependent just on how the switching process is carried out. Moreover the question of why the system selects a particular value of ψ = 0 and, in addition, of how a pseudo-periodic orbit can be put in correspondence with a quantum state,

would still remain open ones. In this respect we believe that although our arguments, based on periodic orbits, are particular sound it might still be possible that this is not jet sufficient to properly describe quantum stationary states of the BEC cloud and quantum calculation should be in order. In any case further evidences in support of some of the possibilities expressed here should also be searched for.

VI. CONCLUSION

In this work, we have used the classical equations of motion to investigate some properties of the magnetic confinement of a neutral atom in the TOP trap. Our initial purpose was that of understanding the origin of an anomalous shift in the vertical position of a BEC cloud that has been observed experimentally in Ref. [3]. We have calculated analytically this shift and found that this is not in agreement with these measurements. In spite of this, our calculations show that, in principle, there is no reason for this shift being limited only to the standard compensation shift 2ω b , as it is the case in the symmetrical trap. We have shown that an anomalous shift exists in any case in the triaxial trap. In fact, the anomalous shift calculated here is produced by the interplay between the particle space dynamics, and that of its magnetic moment driven by the torque equation. However its magnitude results to be extremely small for the parameters used in the normal trap. On the contrary, if such a shift were to exist in line with the calculations made in [3], of which account has been given here, it would be independent of the kind of trap used, i.e. symmetrical or triaxial one, and would depend basically by the peculiarity of the switching precess. In addition we have observed a bi-stability regime in the range of values of the b parameter where the adiabatic regime is lost. The systematic study of this bi-stability regime and of the consequent unstable-dynamics, deserves a further investigation.

In conclusion, the general feeling is that, the approach based on the classical model, leaves so many unanswered questions, that we have the only possible alternative to explain the experimental results by means a quantum mechanical calculation. In fact, in the quantum mechanical approach, the periodicity condition of both the space and the internal spin variables, could be much more stringent than in the classical case thus amplifying an effect that appears otherwise to be quite a small one.

APPENDIX A: DYNAMICS OF THE SMALL OSCILLATING COMPONENTS

In the present appendix we consider a dipole moving in a magnetic field given by a static component with added small sinusoidal components. In particular, the present section

  axes are b x = b x , b y = b y , b z = -2 b z .

FIG. 1 :

 1 Fig.1

2 0 2 µb 2 mω 2 n

 2222 with b 0 standing for the order of magnitude of | b 1 | and | b 2 | of (A1). Now, for typical experimental parameters, b 0 = 3 ⊥ ≈ µb 2 mω 2 = 0.92 (here we have used b = 10.). Thus the amplitude of the corrections due to the second term in Eq. (A6) results of the order of 10 -8 . The first term of the same equation, is instead of the order of ω B

FIG. 2 :

 2 FIG.2:The shift obtained with the theoretical calculations (continuous line) and the zero order compensation shift (dotted line).

b.FIG. 3 :

 3 FIG.3: Comparison between the theoretical shift (continuous line) and the numerical one (crosses).

FIG. 4 :FIG. 5 :

 45 FIG. 4:The cosine of the angle between the magnetic field and the magnetic dipole at time t = 0 of the periodic orbit, as a function of the b parameter for both σ + (continuous line) and σ -(big dots) polarization of the bias field.

ACKNOWLEDGMENTS R. F. thanks B. Z., a right man, a good friend and a brilliant physicist. Bruno Zambon died on October 15 2006 in Pisa, before the publication of the present paper.

is aimed at derive the small deviations from the anti-alignment, due to the effect of the sinusoidal components like the ones of Eqs. (20). Let the total field acting on the magnetic dipole be

where n 1 and n 2 are small deviations from the anti-alignment. Here n 0 is the zero order solution which is anti-aligned with respect to the field B 0 . By inserting the above expressions in the torque equation for the dipole unit vector, and by neglecting the higher order terms coming from both the product between the small rotating components of the magnetic moment and the small rotating components of the magnetic field, one obtains the following equations

These equations can be solved for n 1 and n 2 , thus we obtain [START_REF]Here the condition n 1 • B 0 = 0 and n 2 • B 0 = 0, as it results from the Eqs[END_REF]

Here the subscript ⊥ means perpendicular to the static magnetic field. To the smallest order in the oscillating field amplitude, the constant term in the torque equations yields

This expression shows that the anti-alignment condition between n 0 and B 0 can no longer be satisfied. Thus we calculate the last two terms of the above equation by using the Eqs.

(A4), after a bit of algebra we obtain

This equation indicates how n 0 or B 0 must be modified in order to satisfy the torque equation.

We assume that the dipole and the magnetic field B 0 lies both in the x -y plane of the rotating frame and that the system moves in a quasi-periodic orbit. Let ϑ be the angle between the z axis and the field and ψ the angle formed by the dipole with the field, both measured with reference to the y axis rotations. It is not difficult to calculate the change in the vertical position of the orbit that occurs as a consequence of a change of the dipolefield angle ψ and of the bias pulsation frequency ω. Indeed from Eq. ( 11), noting that

which is the change that the pseudo-periodic orbit angle ϑ undergoes. Thus from Eq. ( 12)

we obtain that ∆z and ∆ω are correlated as follows

where the small term corresponding to the inhomogeneous field has been neglected with As a consequence of the instantaneous modification of the bias field rotation, the system performs oscillations around a position which is displaced by a quantity ∆z with respect to the initial position of the quasi-periodic orbit.


L Mazzola

S Maniscalco

J Piilo

K. -A Suominen

Interplay between entanglement and entropy in two-qubit systems

We study the exact entanglement and entropy dynamics of two qubits interacting with a common zero-temperature non-Markovian reservoir. It is a commonly held view that entanglement loss due to environmental decoherence is accompanied by loss of purity of the state of the system. We demonstrate that such intuitive picture does not always apply: the deterioration of entanglement and purity do not necessarily come together; i.e., revivals of entanglement can be accompanied by deterioration of purity. To complete our investigation on entanglement-mixedness interplay we consider the case of initial mixed states and study how the entanglement dynamics and its revivals are related to both the initial purity and the initial entanglement.

Introduction

Entanglement and purity are fundamental properties determining the usefulness of a quantum state in quantum information protocols. Entangled states play a central role in quantum key distribution, superdense coding, quantum teleportation and quantum error correction [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF]. Purity of quantum states crucially determines the efficiency of certain applications based on quantum theory, e.g., it characterizes the security of quantum key distribution protocols [START_REF] Usenko | [END_REF]. The determinant role of these quantities has motivated research in the field of quantum estimation theory with the aim of envisaging protocols for optimal estimation of purity and entanglement [3].

However, realistic quantum systems are never completely isolated from their surroundings. The inevitable interaction between a system and its environment leads to decoherence phenomena and loss of quantumness [START_REF] Breuer | The Theory of Open Quantum systems[END_REF]. So in realistic situation the state of a quantum system is mixed. That is why special status is accorded to the socalled maximally entangled mixed states (MEMS), i.e., states exhibiting the maximum amount of entanglement for a given degree of mixedness. The classification of bipartite MEMS and the study of their applications has received a lot of interest [START_REF] Wei | [END_REF], representing an important step in the comprehension of the trade-off between entanglement and mixedness. On the other hand, the dynamical relation between entanglement and mixedness for specific exemplary cases has been far-less investigated [6,7,8,9] especially in non-Markovian system [10].

It is a commonly held view that, whenever a small system interacts with the environment, the entanglement deteriorates with decreasing purity and the state goes toward a statistical mixture. However, this quite intuitive picture does not always apply. We demonstrate for two qubits interacting with a common zero temperature non-Markovian reservoir that the dynamical interplay between entanglement and mixedness depends dramatically on the initial state of the two qubits. The study of the dynamics of such a physical system, being exactly solvable, gives us the chance to investigate analytically the connection between entanglement and mixedness when sudden death/birth and revivals of entanglement occur.

The phenomenon of entanglement sudden death (ESD), consisting in the complete disappearance of bipartite entanglement in a finite time, is believed to set a limit on the life-time and usability of entanglement for practical purposes. Hence, lots of efforts have been done in order to understand the conditions under which ESD occurs [11,12,13,14,15,16,17,18,19]. In [19] we have studied the entanglement time evolution when two qubits are prepared in a Bell-like state with two excitations. Here, we extend our investigation to include the entanglement dynamics for initially mixed states. In particular we study a class of states having an "X"-structure density matrix, namely the extended Werner-like states (EWL). This class of states plays a crucial role in many applications of quantum information theory, such as teleportation [20] and quantum key distribution [21]. The extension of our analytical results to the case of initial mixed states also allows us to study how the entanglement dynamics and its revivals are related to the purity and the amount of entanglement of the initial state.

The paper is organized as follows. In Sec. II we review the exactly solvable model of two qubits interacting with a common Lorentzian structured reservoir. In Sec. III we focus on the entanglement-mixedness interplay comparing concurrence and von Neumann entropy dynamics for initial Bell-like states; we also prove that sudden death of entanglement can never occur if the qubits are initially in a mixed state having at most one excitation. In Sec. IV we present the entanglement time evolution of EWL states, and compare our results with those obtained in [22] for independent structured reservoirs. Finally, we summarize our results in Sec. V.

The model

In this section we describe the model used to study the dynamics of two two-level systems (qubits) interacting with a common zero-temperature bosonic reservoir. Our approach is non-Markovian and non-perturbative, i.e., it does not rely on either the Born or the Markov approximations [19].

The Hamiltonian of the system, in the rotating wave approximation, is given by

H = H 0 + H int , H 0 = ω 0 (σ A + σ A -+ σ B + σ B -) + k ω k a † k a k , (1)
H int = (σ A + + σ B + ) k g k a k + h.c., (2)
where σ A ± and σ B ± are the Pauli raising and lowering operators for qubit A and B respectively, ω 0 is the Bohr frequency of the two identical qubits, a k and a † k , ω k and g k are the annihilation and creation operators, the frequency and the coupling constant of the field mode k, respectively.

In Ref. [19] we observed that the dynamics of two qubits which are identical and equally coupled with the same reservoir can be conveniently described in the basis

{|00 , |+ , |-, |11 } with |+ = (|10 + |01 )/ √ 2 and |-= (|10 + |01 )/
√ 2 the superradiant and sub-radiant states, respectively. The advantage of this basis is that subradiant state turns out to be decoupled from the other states and the reservoir, as a consequence it does not decay. Therefore the two qubits can be effectively described by a three-state system {|00 , |+ , |11 } coupled to the reservoir in ladder configuration, and the decoupled sub-radiant state.

In order to solve the dynamics of the two qubits one needs to specify the properties of the environment. In Ref. [19] we assumed that the two qubits interact resonantly with a non-Markovian Lorentzian structured reservoir, such as the electromagnetic field inside a lossy cavity [START_REF] Haroche | Exploring the Quantum: Atoms, Cavities, and Photons[END_REF], having spectral distribution

J(ω) = Ω 2 2π Γ (ω -ω 0 ) 2 + (Γ/2) 2 , ( 3
)
where Γ is the width of the Lorentzian function and Ω the coupling strength. There we showed that the dynamics can be exactly solved by means of the pseudomode approach. The pseudomodes are auxiliary variables defined from the properties of the spectral distribution, allowing us to derive an exact Markovian master equation for the extended system comprised of the system of interest and the pseudomode degrees of freedom.

Only one pseudomode is associated to a Lorentzian spectral distribution and the master equation describing the extended system two-qubits-plus-pseudomode is

∂ ρ ∂t = -i[V, ρ] - Γ 2 (a † aρ + ρa † a -2aρa † ), (4)
with

V = √ 2Ω(a|+ 0| + a † |0 +| + a|2 +| + a † |+ 2|), ( 5
)
where ρ is the density matrix of the extended system in the interaction picture and V the effective coupling between qubits and pseudomode. The parameters Γ and Ω defining the spectral distribution in the pseudomode approach describe, respectively, the pseudomode decay rate and the coupling constant between the pseudomode and the ladder system. We solve the master equation in Eqs. (??) and ( 5) using the Laplace transform technique [19], (in Appendix we present the analytic solution for the qubits density matrix elements for a widely used class of initial states.)

In principle this approach can be suitably used to study the non-resonant case in which the spectral distribution of the environment is peaked at a different frequency than the transition frequency of the qubits. However, the price to pay is that the unitary part of the master equation, describing the coupling between the pseudomode and the qubits, becomes explicitly time-dependent. This of course greatly complicates the derivation of the analytic solution of the dynamics. The off-resonant case, indeed, does not present new qualitative features with respect to those we highlight in this paper, i.e. the relationship between loss of entanglement and purity.

In this work we are mainly interested in the time-evolution of entanglement and mixedness in the effects that the non-Markovian reservoir induces on the correlation between the two qubits. To quantify entanglement we use the Wootters concurrence [START_REF] Wootters | [END_REF], defined as

C(t) = max{0, √ λ 1 - √ λ 2 - √ λ 3 - √ λ 4 }, where {λ i } are the eigenvalues of the matrix R = ρ(σ A y ⊗ σ B y )ρ * (σ A y ⊗ σ B y
), with ρ * the complex conjugate of ρ and σ

A/B y are the Pauli matrices for atoms A and B. This quantity attains its maximum value of 1 for maximally entangled states and vanishes for separable states.

We focus now on the dynamics of initial "X" states. We use the method presented in Ref. [19] (we just briefly recalled) to calculate the time evolution of the elements of the pseudomode-plus-qubits density matrix. By tracing out the pseudomode degree of freedom we obtain the time evolution of the density matrix of the two qubits

ρ(t) =      a(t) 0 0 w(t) 0 b(t) z(t) 0 0 z * (t) c(t) 0 w * (t) 0 0 d(t)      , (6)
which is written in the basis {|00 , |10 , |01 , |11 }. Due to the structure of the differential equations derived from the master equation in ( 4) and (5), the "X" form of the state is preserved during the evolution. In the Appendix we present the analytical solution in the Laplace transform space for a particular type of "X" state.

For this class of states the concurrence assumes a simple analytic expression

C(t) = max{0, C 1 (t), C 2 (t)}, (7)
where

C 1 (t) = 2|w(t)| -2 b(t)c(t), (8)
C 2 (t) = 2|z(t)| -2 a(t)d(t). (9)
We notice that coherences give a positive contribution to C 1 (t) and C 2 (t) and so to concurrence, while the negative parts involve populations only.

In the next section we will also consider the evolution of the mixedness of the two qubit state, quantified by the von Neumann entropy

S(ρ) = -Tr{ρ(t) log N (ρ(t))}. ( 10
)
The von Neumann entropy is equal to zero for pure states, and attains its maximum value 1 (equal to log N N with N the dimension of the Hilbert space) for a maximally mixed state. Alternatively, to quantify the amount of mixedness in a quantum state one can use purity, defined us

P (ρ) = Tr{ρ 2 } (11)
The purity can range between unity, corresponding to a completely pure state, and 1/d (with d the dimension of the density matrix) for a maximally mixed state.

All the physical quantities we consider in this work (and the associated plots) have been calculated directly from the knowledge of the density matrix elements of Eq. ( 6) in the Laplace transform space. The initial conditions on the state are directly inserted in those analytic solutions, which are then anti-transformed and manipulated to calculate the quantity under investigation.

Entanglement-mixedness dynamical interplay

Here, we study and compare the dynamics of concurrence and von Neumann entropy of two qubits prepared in the Bell-like states

|Φ = α|10 + e iθ (1 -α 2 ) 1/2 |01 , (12)
and

|Ψ = α|00 + e iθ (1 -α 2 ) 1/2 |11 . ( 13
)
Our aim is to understand the interplay between these two different quantum properties, in particular when peculiar phenomena such as ESD or "entanglement sudden birth"(ESB) [18,19,25,26], and revivals of entanglement occur.

The entanglement dynamics of two qubits in a Lorentzian structured reservoir, prepared in a Bell-like state with two excitations as in (13), has been presented in [19]. The evolution of entanglement for a Bell-like state as in (12) has been studied also in [27]. For both the Bell-like states in ( 12) and ( 13) the entanglement dynamics is the result of two combined effects: the backaction of the non-Markovian reservoir and the reservoir-mediated interaction between the qubits. The memory effects due to the non-Markovianity of the reservoir causes oscillations in entanglement dynamics which in fact are not present in the Markovian case [15]. The sharing of the reservoir plays also a special role. Indeed, the common reservoir provides an effective coupling between the qubits, and consequently creates quantum correlations between them. As a result, qubits prepared in a factorized state can become entangled due to the interaction with the common reservoir, in contrast with the independent non-Markovian reservoirs case [16,17].

The results in [27] show that ESD does not occur for a Bell-like state with one excitation as in (12) for any value of α 2 , even when the dipolar interaction between the qubits is included [28]. Actually, a straightforward calculation shows that for every pure or mixed state of the qubits containing at most one excitation, ESD and ESB cannot take place. In fact, the density matrix describing a generic mixed state with maximum one excitation, written in the same basis of ( 6), has the form

ρ(t) =      a(t) j(t) k(t) 0 j * (t) b(t) z(t) 0 k * (t) z * (t) c(t) 0 0 0 0 0      . ( 14
)
The expression of the concurrence, for any value of the parameters, is

C(t) = max{0, 2|z(t)|}. ( 15
)
Here concurrence is directly given by the coherence between the |10 and |01 states. Since the coherence vanishes in asymptotic way, there cannot be ESD for any generic state with maximum one excitation. Analogously, entanglement can be smoothly generated but it cannot suddenly appear. This result does not depend on the degree of purity of the state. This is true as long as the form of the density matrix in ( 14) is maintained. On the other hand, if some population is transferred to the two excitations state then ESD can appear. This is the case of two qubits in a Bell state interacting with a non-RWA common reservoir [29].

Entanglement and von Neumann entropy time evolution

We will now try to understand how environmental decoherence affects the entanglement and the purity of the state. One could intuitively expect that the interaction with the environment not only deteriorates the entanglement content but also destroys purity.

Stated another way, one could think that the loss of entanglement is related to the transition from an entangled quantum superposition such as those given by equations (12) and ( 13), to the corresponding statistical mixture. In fact, as we will show in the following, this intuitive picture holds for two-photon Bell-like states and for those states containing the maximally excited state |11 . In such cases the concurrence and the von Neumann entropy dynamics are clearly connected. However, we have discovered that the dynamical relation between entanglement and entropy is far less trivial for Bell-like states with one excitation.

We mentioned in the previous section that the state of the two qubits is effectively equivalent to a four-state system, in which three states {|00 , |+ , |11 } interact in ladder configuration with the electromagnetic field, and the fourth state, the sub-radiant state |-, is completely decoupled from both the other states and the electromagnetic field.

We first look at the evolution of the von Neumann entropy of the qubit pair when the state is initially prepared in a Bell-like state of the form (12) and (13). We begin by considering the Ψ state dynamics as a function of the dimension-less time γ 0 t where the parameter γ 0 = 4Ω 2 /Γ is the Markovian decay rate of the atoms, i.e., the inverse of the atomic relaxation time in the Markovian limit. It is particularly interesting to see how the degree of mixedness of the state evolves when the qubits undergo ESD. Since ESD occurs when α 2 ≤ 1/4 [19], we choose α 2 = 1/20 and we set θ = 0 for the sake of convenience. In figure 1 we notice that ESD appears for high value of mixedness of the system, i.e., concurrence and von Neumann entropy oscillate out of phase. Stated another way, the revivals of entanglement appear roughly in correspondence of the minima of the von Neumann entropy, when the state becomes purer. Since the initial state does not contain any sub-radiant state component, the only states involved are |00 , |+ , |11 . For this particular initial conditions the super-radiant state population explains the link between entropy and concurrence. In fact, for the chosen initial state ρ ++ (t) = 2 b(t)c(t) and the time-dependent part of the concurrence is 2|w(t)|-ρ ++ (t). The von Neumann entropy has a more complicated analytical expression, however it follows closely the time evolution of ρ ++ (t). We can grasp the main features of the dynamics of mixedness by considering the analytical expression of the linear entropy,

S L (ρ) = 1 -Tr[ρ 2 ] = 2(ρ ++ (t) -ρ 2 ++ (t) + a(t)d(t) -|w(t)|).
When the initial state ρ(0) is such that 11|ρ|11 = 0 the flow of probability from the |11 component to the super-radiant component destroys entanglement, even if the super-radiant state is a maximally entangled state. Indeed, whenever the population ρ ++ (t) reaches its relative maxima, the state attains a maximum value of mixedness, the time-dependent part of the concurrence (2|w(t)| -ρ ++ (t)) becomes negative and entanglement disappears. On the other hand, whenever the population of the super-radiant state reaches a minimum, the population of the |11 excited state and the |00 ground state have their maxima (see inset in figure 1), and the system goes toward the initial Bell-like state |Ψ . As a consequence the system becomes purer and entanglement is partially recovered.

The same kind of dynamical relation between entanglement and entropy exists if the qubits are prepared in a state of the form

ρ(0) = (α 2 |0 A 0 A | + (1 -α 2 )|1 A 1 A |) ( 1 6
)

⊗ (α 2 |0 B 0 B | + (1 -α 2 )|1 B 1 B |). (17)
The dynamics of this state in a common non-Markovian reservoir is characterized by ESB and revivals of disentanglement, as we have shown in [19]. In figure 2 we see how these interesting features are related to the degree of purity of the state. As we have seen before, the positions of maxima and minima of the entropy and of the superradiant state population match. The sudden creation of entanglement happens roughly when the entropy hits its first minimum. Entanglement is again lost when the amount of mixedness increases, and it reappears again when the entropy reaches another minimum.

A completely different dynamics characterizes the system when the qubits are prepared in a Bell-like state as the one of (12). In fact in this case there is no ESD, and revivals of entanglement are not necessarily related to a reduction in mixedness. To start with consider as initial state the super-radiant Bell |+ state. As shown in figure 3 revivals of entanglement are always accompanied by increase of entropy, so the intuitive picture fails. Moreover the zeroes of entanglement and entropy coincide. For those times, in fact, the system goes into the ground state which is pure and factorized. When some population returns in the super-radiant state entanglement is recovered, and the state is again mixed. This examples certainly shows that loss of entanglement is not simply due to a transition from quantum superposition to statistical mixture.

A more complex situation appears for a generic one-photon Bell-like state. When α 2 = 1/2 and/or θ = 0 the initial state is a superposition of super-radiant and sub-radiant states. Although the sub-radiant component does not evolve in time, being decoupled from the super-radiant and ground states, its presence affects the entanglement and entropy dynamics. For the sake of convenience we choose α 2 = 1/5 and θ = 0. In figure 4 we see that entropy still follows closely the time evolution of the super-radiant state population, having its relative minima in the same positions of the zeroes of ρ ++ (t). On the contrary, entanglement has new relative maxima when the population of the super-radiant state is zero. This is due to the presence of the sub-radiant state. In fact, in this case both the super-radiant and the sub-radiant states contribute to the total entanglement. Thus there are two different sets of entanglement maxima, those associated with the maxima of the super-radiant state population, and those associated to the sub-radiant state. Entanglement is zero whenever the population of the super-radiant state ρ ++ (t), the population of the sub-radiant state ρ --(t), and the absolute value of the coherence between super-radiant and sub-radiant states ρ +-(t), are equal. This can be explained having in mind equations (22) in the Appendix, where b(t), c(t) and z(t) are written as a function of those quantities. Indeed when ρ ++ (t), ρ --(t) and |ρ +-(t)| have equal value, then z(t) is equal to zero, irrespective of the sign of ρ +-(t). Specifically, this happens in correspondence to the 1st, 4th and 5th zeroes of entanglement, where b(t) = 0 and c(t) = 2ρ --= k. In this case the state of the system becomes (1 -k)|00 00| + k|01 01|, which is clearly not entangled. Similarly, for the 2nd and 3rd zeroes, b(t) = k, c(t) = 0 and the non-entangled state is (1 -k)|00 00| + k|10 10|.

Dynamics in the concurrence-vs-purity phase space

To summarize and better understand the entanglement-mixedness dynamical interplay we look at our results from a different perspective. Figure 5 shows the time-evolution in the entanglement-purity phase diagram for different class of initial states of the two qubits. In particular, figure 5 a) shows the parametric plot of concurrence and purity for a Bell-like state with two excitation Ψ with α 2 = 1/20. In figure 5 b) the state of preparation of the system is a Bell-like state with one excitation Φ with α 2 = 1/5 and θ = 0. The two plots are qualitatively very different.

The dynamical feature we noticed from the analysis of figure 1 is even more apparent in figure 5 a): deterioration, and revivals as well, of entanglement and purity come together during the dynamics of a Ψ state. The initial moment of the time evolution corresponds to the upper corner of the C-P diagram, the system evolves toward smaller entanglement and purity, till sudden death of entanglement occurs; afterwards entanglement and purity revive accordingly. We plot the curve up to the time 30 γ 0 t, similar structures appear later on in the dynamics until the system reaches its pure and separable asymptotic state, the ground state (at the lower right corner of the CP diagram).

Analogously, in part b) of figure 5 the upper right corner of the C-P diagram corresponds to the starting time of the dynamics. Entanglement and purity decrease together until at (P = 0.5, C = 0.3) a change of behavior occurs: after that while concurrence is still decreasing purity starts to increase. Once more a sudden change in the dynamical interplay has place at the point when concurrence becomes zero and entanglement and purity behave again in the same way. As time passes the system starts travelling back on the same curve, therefore there is a period of time for which while purity is decreasing concurrence increases. This is a clear sign of non-Markovianity. In fact, as pointed out by Ziman and Buzek [7] the C-P line associated to a Markovian process must be non-increasing respect to concurrence, this prevents the presence of loops in the diagram even though allows the growth of purity when concurrence is decreasing. Here, revivals of concurrence can be accompanied by decrement of purity, and the system travels partially back and forth on the same curve till it reaches its asymptotic state which is a mixed non separable state, statistical mixture of the ground and sub-radiant state. To the best of our knowledge this type of dynamics was not yet been observed.

Werner states

In this section we study the two qubits entanglement dynamics of extended Wernerlike states in a common zero-temperature Lorentzian structured reservoir. Our goal here is to study how, starting from an initial state that is not perfectly pure, as in realistic experimental conditions, affects the entanglement dynamics, and in particular the occurrence of ESD and ESB phenomena.

Terminology and previous works

The standard two-qubit Werner state, introduced in 1989 by Werner [30], is defined as

ρ W = r|--| + 1 -r 4 I, (18)
where |is the singlet state. In [30] Werner demonstrated that while pure entangled states always violate the Bell inequality, mixed entangled states might not. The Werner state is the first entangled state to be proven to satisfy all Bell inequalities [30]. The generalized or Werner-like states are defined as

ρ W L = r|M M| + 1 -r 4 I ( 19
)
with |M one of the four maximally entangled Bell states. For a given r, Werner and Werner-like states exhibit the same entanglement. A further generalization is the extended Werner-like states (EWL), containing a non-maximally entangled state part, which are defined as

ρ Φ EW L = r|Φ Φ| + 1 -r 4 I, (20)
where |Φ = α|10 + e iθ (1 -α 2 ) 1/2 |01 , and

ρ Ψ EW L = r|Ψ Ψ| + 1 -r 4 I, (21)
with

|Ψ = α|00 + e iθ (1 -α 2 ) 1/2 |11 .
The dynamics of entanglement of two qubits prepared in Werner, Werner-like or extended Werner-like states has attracted a lot of attention. The appearance of ESD has been demonstrated for two qubits prepared in a Werner-like state and interacting with Markovian independent reservoirs [31,32] or with independent noisy channels [33]. Violation of Bell inequality has been examined in independent thermal reservoirs [34,35]. Entanglement has been studied in a Markovian common thermal reservoir [36], in the presence of collective dephasing [37] and in the maximal noise limit [38]. Finally, the dynamics of an EWL in two independent Lorentzian structured reservoirs has been investigated in [22].

Moreover, Werner and Werner-like states have been used so far in many applications in quantum information processes such as teleportation [20] and entanglement teleportation [39]. The experimental preparation and characterization of the Werner states have also been widely investigated. Werner states are prepared via spontaneous parametric down-conversion [40] or using a universal source of entanglement [41], and used in ancilla-assisted process tomography [42] and secure quantum key distribution [21].

Entanglement dynamics

We evaluate the dynamics of the entanglement of EWL states as a function of the initial amount of mixedness, controlled by the purity parameter r, and as a function of the initial degree of entanglement measured by α 2 .

The study of the evolution of this kind of states allows to infer how the degree of purity of the initial state influences the entanglement dynamics. We also compare our results with those obtained in [22] for two qubits in two independent Lorentzian structured reservoirs.

Figures 6, 7 and 8 show the entanglement dynamics as a function of the dimensionless quantity γ 0 t (with γ 0 = 4Ω 2 /Γ the Markovian decay rate of the atoms) and of the purity parameter r, for the two EWL states in (20) and (21). In figure 6 the qubits are initially in the state given by (20) with α 2 = 1/2 and θ = 0. When r = 1 the qubits are prepared in the super-radiant state, entanglement exhibits oscillations and ESD is never present. However, whenever a little amount of mixedness is added, periods of finite-time disentanglement immediately occur. This is due to the appearance of some population in the excited state |11 . As a consequence the oscillating part of the concurrence {2|z(t)| -2 a(t)d(t)} can become negative. Due to the non-Markovianity of the system and to the effective coupling provided by the common reservoir, ESD regions are followed by revivals, and for long times entanglement reaches a stationary value. When r < 1/3 the state is initially factorized, but as time passes, because of the reservoir-mediated interaction, entanglement between the qubits is suddenly created. The non-Markovianity of the reservoir enriches the dynamics causing revivals of disentanglement. Eventually the entanglement reaches a non-zero stationary value. The amount of entanglement that has been created depends on the population of the sub-radiant state, which carries the entanglement.

If the qubits are prepared in the state given by ( 20) with α 2 = 1/2 and θ = π, the dynamics is completely different compared to the previous θ = 0 case, as shown in figure 7. In fact when r = 1 the qubits are prepared in the sub-radiant state and concurrence does not evolve, being at any time equal to 1. When the mixed part is present little oscillations appear, but eventually the concurrence attains the stationary value (1 + 3r)/4. When r < 1/3 the initial state is factorized, entanglement is suddenly created, a revival of disentanglement is present, and again after some oscillations concurrence reaches its stationary value.

When α 2 = 1/2 and/or θ = 0, π the non-mixed part of the initial state is a superposition of super-radiant and sub-radiant states. Thus the entanglement dynamics depends on the weights of those two states. As a result there is a wide variety of entanglement dynamics in between the two asymptotic behaviors described above.

In figure 8 the qubits are prepared in the state given by ( 21) with α 2 = 1/2 and θ = 0. Note however that the results are independent of the choice of the relative phase θ. For r = 1 the initial state is the Bell state (|00 + |11 )/ √ 2. For this initial condition not only there is no ESD, but also concurrence vanishes only for infinite time. When an increasing amount of mixedness is present in the initial state, finitetime disentanglement appears. ESD is then followed by revivals and, as expected, a certain amount of entanglement is preserved. Analogously to the other EWL state, for r < 1/3 the initial state is factorized. For the reasons mentioned above, entanglement is suddenly created, momentarily deteriorated, and it finally goes to the stationary value r/4, coinciding with the population of the sub-radiant state. For this type of EWL state different choices of α lead to the same qualitative behavior of the entanglement dynamics.

The results we just described for two qubits in a common structured reservoir are quite different from those presented in [22] for two qubits in two independent Lorentzian structured reservoirs. First of all, if two uncoupled qubits interact with two independent reservoirs, entanglement cannot be created from a factorized state. Hence, the ESB region, characterizing the dynamics in a common reservoir, is not present in the results of [22], for any of the EWL states. Moreover, for qubits prepared in the initial state of ( 21), we notice that the reservoir-mediated interaction between the qubits keeps the value of concurrence higher compared to the two independent reservoirs case.

The crucial difference between the common and independent reservoirs cases is that for qubits in two independent reservoirs the decoherence-free sub-radiant state does not exist. When the qubits are in two independent reservoirs, the entanglement dynamics of the |+ and |states are the same. As a consequence the relative phase θ in ( 12) and ( 20) does not affect the results. This is definitely in contrast with the results we present in figures 6 and 7, showing two completely different asymptotic behavior when changing the relative phase θ of a π factor.

Conclusive Remarks

In this paper we have investigated the dynamical relation between entanglement and mixedness in a system of two qubits interacting with a common zero temperature non-Markovian reservoir. We have used the exactly solvable model presented in [19], in the Appendix we attach the analytical solution in the Laplace transform space in the case of extended-Werner-like states. All the plots we have presented in this work have been derived by manipulation of the exact analytic solution in the Laplace transform space.

Our results show that decoherence induced by the environment does not necessarily affect in the same way the entanglement content and the purity of the state. Nevertheless, for qubits prepared in two-photon Bell-like states or in states comprising the maximally excited state |11 , loss of entanglement is generally accompanied by increasing of mixedness. ESD seems to appear when the state becomes highly mixed, whereas revivals of entanglement are associated to minima of the entropy, where the state becomes purer. Moreover when starting from a factorized state, sudden birth of entanglement occurs for lower values of mixedness, while revivals of disentanglement are accompanied by peaks of the entropy.

The dynamical connection between entanglement and mixedness becomes much more complex for initial Bell-like state with one excitation. However, by looking at the dynamics of the population of the super-radiant and sub-radiant states, one sees two sets of maxima of entanglement can be identified. The super-radiant component mainly carries the entanglement when it is maximally populated; when its population is zero, the entanglement is associated to the sub-radiant component. When the populations of these states are equal, the entanglement vanishes.

The dynamical interplay between the entanglement content and the purity of the state of the two qubits can be more clearly seen in the concurrence-versus-purity phase diagram. There, it is apparent that when starting from a Bell-like state with two excitations entanglement and purity behave accordingly, decreasing or reviving at the same time. For an initial Bell-like state with one excitation their trade-off is far less trivial, the diagram shows loops (the system traveling back and forth over the same curve) and revivals of entanglement can take place in correspondence of decrement of purity. Those are further signs of the non-Markovianity of the system [7].

We have also considered the entanglement dynamics of a particular class of mixed states, the extended Werner-like states. We have demonstrated that the amount of purity of the initial state plays a key role in the entanglement dynamics, controlling the appearance of ESD and ESB phenomena. As expected, for stronger non-Markovian conditions, the dynamics exhibits stronger and longer lasting entanglement oscillations as well as an increasing number of dark periods and revivals as well. On the other hand, in the Markovian regime, no oscillations are present, however the basic features of the dynamics, and so the ESD and ESB regions, are still present.

Entanglement and purity are two different physical quantities characterizing the degree of non-classicality of a quantum state. Understanding their interplay is important not only from an applicative but also from a fundamental point of view. We believe that the present results contribute to the comprehension of their dynamical relations. d(t) = 1 -a(t) -b(t) -c(t) = ρff (t), w(t) = ρaf (t), ρ++ (t) = ρdd (t) + ρee (t).

The subscripts a, b and c refer to the states in which both the qubits are in the ground state, and the pseudomode has 0, 1 or 2 excitations respectively. In d and e the atomic system is in the super-radiant state, and the pseudomode has 0 or 1 excitations. f is the state having 0 excitations in the pseudomode and both the qubits in their excited states.

The analytic expressions are given in the Laplace transform space. In the following we provide the solutions for a EWL state of the form (20): ρaa (t) = L -1 -1 4s (-1 + r + 64Γ 2 Ω 4 (-1 + r)l(s) j(s)k(s)

- 8ΓΩ 2 (1 + r + 4α √ 1 -α 2 r cos θ) j(s) ) , (23)
ρbb (t) = L -1 2Ω 2 j(s) (1 + r - 8ΓΩ 2 (-1 + r)l(s) k(s) + 4α √ 1 -α 2 r cos θ) , (24) ρcc
(t) = L -1 - 48Ω 4 (-1 + r)(Γ + s) k(s) , (25) ρdd
ρ+-(t) = L -1 r(Γ + 2s)(-1 + 2α 2 + 2iα

√ 1 -α 2 sin θ)) 2(4Ω 2 + s(Γ + 2s) , (29)
ρ-+ (t) = ρ * +-(t), (30) ρ--(t) = r(1 -2α

√ 1 -α 2 cos θ) 2 + 1 -r 4 , (31)
ρaf (t) = 0. (

The solutions for both the EWL states are written as functions of k(s), j(s), l(s) defined as k(s) = (16ΓΩ 2 + 2Γ 2 s + 24Ω 2 s + 3Γs 2 + s 3 ) (3Γ 3 + 28ΓΩ 2 + 11Γ 2 s + 24Ω 2 s + 12Γs 2 + 4s 3 ), (33) j(s) = (Γ + 2s)(8Ω 2 + s(Γ + s)),

l(s) = (6Γ 3 + 31Γ 2 s + Γ(56Ω 2 + 45s 2 ) + 20(3Ω 2 s + s 3 )).

(

) 35
Here we present the solutions for a EWL state of the form (21): ρaa (t) = L -1 -1 4s (-1 + r -4α 2 r + 8ΓΩ 2 (-1 + r) j(s)

+ 64Γ 2 Ω 4 (-1 + (-3 + 4α 2 )r)l(s) j(s)k(s) ) , (36)
ρbb (t) = L -1 2Ω 2 j(s) (1 -r - 8ΓΩ 2 (-1 + (-3 + 4α 2 )r)l(s) k(s) ) , (37) ρcc
(t) = L -1 - 48Ω 4 (-1 + (-3 + 4α 2 )r)(Γ + s) k(s) , (38) ρdd

Figure 1 .

 1 Figure 1. (Color online) Dynamics in a common Lorentzian structured reservoir as a function of scaled time for two atoms prepared in the Bell-like state Ψ with α 2 = 1/20 and θ = 0. Solid red line is concurrence; dotted-dashed blue line is von Neumann entropy; dotted green line is the population of the super-radiant state ρ ++ (t). In the inset: dotted green line is the super-radiant state population; dotted-dashed black line is the ground state population; solid light blue line is the excited state |11 population.

Figure 2 .

 2 Figure 2. (Color online) Dynamics in a common Lorentzian structured reservoir as a function of scaled time for two atoms prepared in the mixed state (16) with α 2 = 0.75. Solid red line is concurrence; dotted-dashed blue line is von Neumann entropy; dotted green line is super-radiant state population.

Figure 3 .

 3 Figure 3. (Color online) Dynamics in a common Lorentzian structured reservoir as a function of scaled time for two atoms prepared in the Bell-like state Φ with α 2 = 1/2 and θ = 0. Solid red line is concurrence and dotted-dashed blue line is von Neumann entropy. The super-radiant state population overlaps with concurrence.

Figure 4 .

 4 Figure 4. (Color online) Dynamics in a common Lorentzian structured reservoir as a function of scaled time for two atoms prepared in the Bell-like state Φ with α 2 = 1/5 and θ = 0. Solid red line is concurrence; dotted-dashed blue line is von Neumann entropy; dotted green line is super-radiant state population. In the inset: dotted green line is super-radiant state population; dotted-dashed black line is ground state population; dashed pink line is sub-radiant state population; solid yellow line is absolute value of the coherence between super-radiant and sub-radiant states.

Figure 5 .

 5 Figure 5. (Color online) Parametric plot of the time evolution of concurrence and purity of two qubits when initially prepared in: a) a Bell-like state with two excitation Ψ with α 2 = 1/20; b) a Bell-like state with one excitation Φ with α 2 = 1/5 and θ = 0.

Figure 6 .

 6 Figure 6. (Color online) Concurrence as a function of scaled time and r for two atoms prepared in the Werner state (20) (Φ) with α 2 = 1/2 and θ = 0.

Figure 7 .

 7 Figure 7. (Color online) Concurrence as a function of scaled time and r for two atoms prepared in the Werner state (20) (Φ) with α 2 = 1/2 and θ = π.

Figure 8 .

 8 Figure 8. (Color online) Concurrence as a function of scaled time and r for two atoms prepared in the Werner state (21) (Ψ) with α 2 = 1/2 and θ = 0.

  2Ω 2 (-1 + r)(6Γ 3 + 8ΓΩ 2 + 13Γ 2 s + 12Ω 2 s + 9Γs 2 + 2s 3 ) + 31Γ 4 s + 4Γ 3 (38Ω 2 + 15s 2 ) +Γ 2 s(412Ω 2 + 55s 2 ) + 8Γ(40Ω 4 + 91 2 Ω 2 s 2 + 3s 4 ) +4(72Ω 4 s + 26Ω 2 s 3 + s 5 ))) ,

	(t) = L -1	1 4j(s)	(-	1 k(s)	(8ΓΩ 2 (-1 + r)(6Γ 5 + 31Γ 4 s + 20Γ 3 (-2Ω 2 + 3s 2 )
		+5Γ 2 s(-20Ω 2 + 11s 2 ) + 8Γ(56Ω 4 -	17 2	Ω 2 s 2 + 3s 4 )
		+4(120Ω k(s)	,	(27)
	ρff (t) = L -1 -	1 4k(s)	((-1 + r)(6Γ 5

4

s -

2Ω 2 s 3 + s 5 ))) + (8Ω 2 + (Γ + s)(Γ + 2s)) (1 + r + 4α √ 1 -α 2 r cos θ)) ,

(26)

ρee (t) = L -1 -

  (-1 + (-3 + 4α 2 )r)(6Γ 5 + 31Γ 4 s + 20Γ 3 (-2Ω 2 + 3s 2 ) +5Γ 2 s(-20Ω 2 + 11s 2 ) + 8Γ(56Ω 4 -17 2 Ω 2 s 2 + 3s 4 ) +4(120Ω 4 s -2Ω 2 s 3 + s 5 ))) , )r)(6Γ 3 + 8ΓΩ 2 + 13Γ 2 s + 12Ω 2 s + 9Γs 2 + 2s 3 ) ,(40) )r)(6Γ 5 + 31Γ 4 s + 4Γ 3 (38Ω 2 + 15s 2 ) + Γ 2 s + 12Ω 2 + 3Γs 2 + 2s 3 ,(42)

	(t) = L -1 -	1 4j(s)	(8Ω 2 + (Γ + s)(Γ + 2s))(-1 + r)
	+	1 k(s)	(8ΓΩ 2 (39)
	ρee (t) = L -1 -(-1 + (-3 + 4α 2 ρff (t) = L -1 -2Ω 2 k(s) 1 4k(s) ((-1 + (-3 + 4α 2 +Γ 2 s(412Ω 2 + 55s 2 ) + 8Γ(40Ω 4 +	91 2	Ω 2 s 2 + 3s 4 )
	+4(72Ω 4 s + 26Ω 2 s 3 + s 5 ))) ,	(41)
	√ 1 -α 2 ((8Ω 2 + (Γ + s)(Γ + 2s)) 4ΓΩ 2 ρaf (t) = L -1 αr

ρ--(t) = ρ+-(t) = ρ-+ (t) = 0. (

43)

Acknowledgments

We thank Barry Garraway for stimulating discussions. This work has been financially supported by M. Ehrnrooth Foundation, Väisälä Foundation, Turku University Foundation, Turun Collegium of Science and Medicine, and the Academy of Finland (projects 108699, 115682, 115982).

Appendix

Here we present the exact analytic solution for two qubits interacting with a Lorentzian structured reservoir, when the qubits are prepared in an EWL state as in (20) and (21). We provide the expressions for the density matrix element in (6) as a function of the solution ρij of the pseudomode master equation [43] in ( 5) and ( 6) of [19]