K Ziegler

Dynamical creation of entangled bosonic states in a double well

Keywords: numbers: 03.65.Aa, 03.65.Fd, 03.67.Bg

published or not. The documents may come

Dynamical creation of entangled bosonic states in a double well

I. INTRODUCTION

Recent experiments on ultracold gases in optical potentials [1][2][3] and experiments on photons in microwave cavities [4,5] have demonstrated that it is possible to prepare a Fock state as a pure state in a finite-dimensional system. After the preparation of the Fock state, the parameters of the system can suddenly be changed (performing a "quench") such that the Fock state is not an eigenstate of the new system Hamiltonian H. Then the evolution of the many-body state due to the evolution operator exp(-iHt) will lead to a random walk inside the available Hilbert space. The visited states include other Fock states as well as superpositions of Fock states. Typical questions in this context are: what is the probability for visiting different states and how is this affected by the interaction of the particles? A natural quantity for measuring this probability is the spectral density function of the Hamiltonian H with respect to the initial Fock state [6][7][8].

A classical candidate for modeling the evolution of a Fock state is the Hubbard model [9,10]. The corresponding discrete Hamiltonian describes the tunneling of a particle between neighboring potential wells and a local particle-particle interaction. The Hubbard model for bosons (Bose-Hubbard model) was realized as an atomic system in an optical lattice [11]. A possible realization of the Bose-Hubbard model by photons in coupled microwave cavities was proposed recently by Hartmann et al. [12]. An anharmonicity of the microwave cavities plays the role of the photon-photon interaction [8].

The simplest system for discussing the evolution of a Fock state within the Hubbard model is a double well, where particles can tunnel between the two wells. For N bosons the underlying Hilbert space is spanned by the (N + 1)-dimensional Fock base {|0, N , |1, N -1 , ..., |N, 0 }, where l bosons are in one well and N -l in the other well [13][14][15][16]. The initial state is prepared as a Fock state, where all the bosons are in one of the two wells (i.e. |0, N or |N, 0 ), while the tunneling between the wells is turned off. To start the evolution, a "quench" is provided by switching on the tunneling between the two wells. This is realized by a sudden reduction the potential barrier between the wells in an atomic system [1] or by connecting the two microwave cavities with an optical fiber [8,12,17]. A similar experiment was performed with two atomic clouds, subject to weak interaction and separated by an adjustable potential barrier [18,19].

On the theoretical side, mean-field descriptions of the Bose-Hubbard model, such as a Hartree approximation or the Gross-Pitaevskii equation, may work well for clouds with many bosons and weak boson-boson interaction [14]. However, they provide a rather poor approximation for the dynamics of small many-body systems (cf. Ref. [13]). This was also observed in a recent study by Streltsov et al. who compared the results of a simple Hartree (Gross-Pitaevskii) approximation with a sophisticated (multiconfigurational time-dependent) Hartree approximation for bosons [20]. The latter reveals that the bosonic clouds are related to superpositions of Fock states in the form of N00N states

|N 00N = 1 √ 2 |0, N + e iφN |N, 0 . (1)
In the following we will study the Hubbard dynamics of bosons in a double well in more detail. In particular, we are interested in the connection of spectral properties and the formation of N00N states, Confidential: not for distribution. Submitted to IOP Publishing for peer review 7 June 2011 based on a Fock state with all the particles in one well as the initial state. To avoid problems with uncontrolled approximations, we will rely on a full quantum calculation. An exact solution is available in a Fock-state base, as described previously in Refs. [6,8]. The paper is organized as follows: In Sect. II the model, based on the Bose-Hubbard Hamiltonian, is defined and in Sect. II A the dynamics of an isolated quantum system is explained. Then we discuss the dynamics of a noninteracting Bose gas in Sect. III and the dynamics of an interacting Bose gas in Sect. IV. The latter is divided into a study of a symmetric double well (Sect. IV A) and of a double well with mixed interaction, where one well has attractive the other repulsive interaction (Sect. IV B). Finally, we summarize the results of our calculation in Sect. V and discuss them in Sect. VI.

II. MODEL

The many-body Hamiltonian Ĥ of N bosons with mass m reads

Ĥ = N j=1 p 2 j 2m + V (r j ) + N j,k=1 U (r j , r k ) , (2)
where p j is the momentum of a boson, V (r j ) is the one-body potential of the double well and U r j , r k ) is the two-body interaction potential. For the latter we assume that it decays very quickly with the distance |r jr k | of the particles. This implies that particles located in different wells do not interact with each other. In a general situation, the two-body interaction can be attractive in one and repulsive in the other well. Then the many-body Hamiltonian is expressed in Fock-state representation as

• • • N -k, k|r 1 , r 2 , ..., r N r 1 , r 2 , ..., r N | Ĥ|r 1 , r 2 , ..., r N |r 1 , r 2 , ..., r N |N -k , k d 3 r 1 • • • d 3 r N = N -k, k|H|N -k , k . (3)
For the new Hamiltonian H, which acts in the Hilbert space spanned by the Fock base, we can use a generalized Bose-Hubbard Hamiltonian with different local interaction parameters U 1 , U 2 in the corresponding wells as a reasonable approximation

H = J(a † 1 a 2 + a † 2 a 1 ) + U 1 (a † 1 a 1 ) 2 + U 2 (a † 2 a 2 ) 2 , (4)
where a † j (a j ) are creation (annihilation) operators for bosons in the Fock states. H, which describes tunneling between the two wells and the local interaction inside the well with interaction strength U j , gives us a complete quantum description of the different Fock states and their superpositions. In particular, we can employ it to study the evolution of a Fock state to a N00N state of Eq. ( 1). This will be used subsequently for two special cases. The first case is a symmetric double well, where U 1 = U 2 . The second case has opposite signs of interactions U 1 = -U 2 which describes an attractive interaction in one well and a repulsive interaction in the other well.

A. Evolution of isolated systems

We consider a system which is isolated from the environment. Furthermore, we assume that the system lives in an (N +1)-dimensional Hilbert space. With the initial state |Ψ 0 we can get for the time evolution of the state

|Ψ t = e -iHt |Ψ 0 (5)
or the evolution of the return probability | Ψ 0 |Ψ t | 2 with the amplitude

Ψ 0 |Ψ t = Ψ 0 |e -iHt |Ψ 0 . (6)
In general, the amplitude Ψ 1 |Ψ t can be expressed via an integral transformation of the resolvent as

Ψ 1 |Ψ t = Ψ 1 |e -iHt |Ψ 0 = Γ Ψ 1 |(z -H) -1 |Ψ 0 e -izt dz , (7)
where the contour Γ encloses all the eigenvalues E j (j = 0, 1, ..., N) of H. With the corresponding eigenstates |E j the spectral representation of the resolvent is a rational function:

Ψ 1 |(z -H) -1 |Ψ 0 = N j=0 Ψ 1 |E j E j |Ψ 0 z -E j = P N (z) Q N +1 (z) , (8)
where P N (z), Q N +1 (z) are polynomials in z of order N , N +1, respectively, with the common denominator

Q N +1 (z) = N j=0 (z -E j ) .
These polynomials are readily evaluated by the recursive projection method (RPM) [6].

The expression in Eq. ( 8) for |Ψ 1 = |Ψ 0 can be interpreted as the bosonic spectral density ρ (E) with respect to the state |Ψ 0 :

ρ (E) = 1 π Im Ψ 0 |(E -i -H) -1 |Ψ 0 = π N j=0 | Ψ 0 |E j | 2 2 + (E -E j ) 2 . ( 9
)
The amplitude of the return probability then reads as the Fourier transform of the spectral density

Ψ 0 |Ψ t = lim →0 ρ (E)e -iEt dE . (10)
Analogously, the overlap Ψ 1 |Ψ t reads in terms of the resolvent

Ψ 1 |Ψ t = 1 π lim →0 Im Ψ 1 |(E -i -H) -1 |Ψ 0 e -iEt dE (11)
with

lim →0 Im Ψ 1 |(E -i -H) -1 |Ψ 0 = π j Ψ 1 |E j E j |Ψ 0 δ(E -E j ) , (12)
provided that the matrix elements are symmetric. The latter is the case for the Hubbard Hamiltonian.

The purpose of the subsequent calculation is to determine the evolution of the Fock state under the influence of the Bose-Hubbard Hamiltonian of Eq. ( 4). In general, this is expressed in the Fock base as

|Ψ t = N j=0 c j (t)|N -j, j (13)
with coefficients c j (t) = N -j, j|Ψ t . For the N00N state we only need to focus on the coefficients c 0 (t) and c N (t). Although we do not anticipate that the other coefficients vanish, for the existence of a N00N state it is sufficient to have simultaneously c 0 (t), c N (t) = 0. Moreover, the resulting entangled state may in general not have c 0 (t) = c N (t) but is still called a N00N state.

Comparing the result in Eq. ( 13) with the expressions in Eqs. ( 10), ( 11), (12), it turns out that the Fourier transform of c 0 (t) and c N (t) are just the imaginary parts of the matrix elements of the resolvent

c0 (E) = 1 π lim →0 Im N, 0|(E -i -H) -1 |N, 0 = j N, 0|E j E j |N, 0 δ(E -E j ) ( 1 4 )
and

cN (E) = 1 π lim →0 Im 0, N|(E -i -H) -1 |N, 0 = j 0, N|E j E j |N, 0 δ(E -E j ) . (15)
These two expressions will be called spectral coefficients, where c0 (E) measures the relative weight | N, 0|E j | 2 . Integration over the energy E gives 1 for this coefficient. The coefficient cN (E) measures the correlation between |N, 0 and |0, N due to the product 0, N|E j E j |N, 0 . The latter is real for a symmetric Hamiltonian. Integration over the energy E gives 0 for this coefficient.

III. DOUBLE WELL: NONINTERACTING BOSE GAS

The Bose-Hubbard Hamiltonian has two simple limits: The local limit J = 0 and the noninteracting limit U 1 = U 2 = 0. In the local limit for a symmetric double well with

U 1 = U 2 pairs Fock states |N -k, k , |k, N -k are doubly degenerate eigenstates with energy E k = U [(N -k) 2 + k 2 ].
A perturbation by a small tunneling term will break the degeneracy. This effect is stronger at lower energies because the parabolic spectrum is denser there. This agrees with a numerical study [13]. The fact that the states |N, 0 and |0, N are very close in energy may support the formation of a N00N state.

In the absence of particle-particle interaction the Bose-Hubbard Hamiltonian H t (i.e. the Hamiltonian in Eq. ( 4) with U 1 = U 2 = 0) describes only tunneling. A straightforward calculation shows that the eigenstate |N -k; k of H t with H t |N -k; k = J(N -2k)|N -k; k has an overlap with the Fock states |N, 0 and |0, N as

N, 0|N -k; k = 2 -N/2 N k , 0, N|N -k; k = (-1) k 2 -N/2 N k . ( 16
)
This implies that the spectral coefficients of Eqs. ( 14) and ( 15) have a binomial form

c0 (E) = 2 -N N k=0 N k δ(E + J(2k -N )) (17) cN (E) = 2 -N N k=0 N k (-1) k δ(E + J(2k -N )) . (18)
A Fourier transformation reveals a periodic behavior of the evolutionary coefficients as

c 0 (t) = N, 0|e -iHt |N, 0 = cos N (Jt), c N (t) = 0, N|e -iHt |N, 0 = (-i) N sin N (Jt) . (19)
Thus the evolution of the Fock state leads to a N00N state with a probability that decays exponentially with N . This is a consequence of the fact that for an increasing N the particles disappear in the (N + 1)dimensional Hilbert space because there is no constraint due to interaction.

IV. DOUBLE WELL: INTERACTING BOSE GAS

The double well with the two Fock states |N, 0 , |0, N as possible initial states can be treated within the RPM. This method is based on a systematic expansion of the resolvent Ψ 1 |(z -H) -1 |Ψ 0 , starting from the initial base {|N, 0 , |0, N }. The method can also be understood as a directed random walk in Hilbert space. This means that in comparison with the conventional random walk the directed random walk of the RPM visits a subspace H 2j only once and never returns to it. In terms of N bosons, distributed over the double well, the subspace H 2j is spanned by the base {|N -j, j , |j, N -j }. A step from H 2j to H 2j+2 is given by the Hamiltonian H in such a way that H 2j+2 is created by acting H on H 2j (cf. App. A). This step is provided by the tunneling of a single boson. Thus, the directed random walk follows a path with increasing numbers j. The directed random walk is the main advantage of the RPM which allows us to calculate the matrix elements

Ψ 0 |(z -H) -1 |Ψ 0 , Ψ 1 |(z -H) -1 |Ψ 0 of
the resolvent on a (N + 1)-dimensional Hilbert space exactly.

A. Symmetric double well

Now we choose U 1 = U 2 ≡ U for the Bose-Hubbard Hamiltonian. Assuming that N is even, all projected spaces H 2j are two-dimensional and spanned by {|N -j, j , |j, N -j } (j = 0, ..., N/2). This leads to a recurrence relation in the base of the two Fock states (|N, 0 , |0, N ) as initial states. The recurrence relation reads (App. A)

g k+1 = a k+1 b k+1 b k+1 a k+1 , g 0 = 1 z -UN 2 /2 1 0 0 1 (k = 0, 1, ..., N/2 -1) (20)
with coefficients

a k+1 = z -fk+1 -J 2 a k (N/2 + k + 1)(N/2 -k) z -fk+1 -J 2 a k (N/2 + k + 1)(N/2 -k) 2 -J 4 b 2 k (N/2 + k + 1) 2 (N/2 -k) 2 (21) b k+1 = J 2 b k (N/2 + k + 1)(N/2 -k) z -fk+1 -J 2 a k (N/2 + k + 1)(N/2 -k) 2 -J 4 b 2 k (N/2 + k + 1) 2 (N/2 -k) 2 (22)
and

fk+1 = U (N/2 + k + 1) 2 + U (N/2 -k -1) 2 .
The iteration terminates after N/2 steps with

g N/2 = a N/2 b N/2 b N/2 a N/2 , ( 23
)
where

a N/2 = N, 0|(z -H) -1 |N, 0 = 0, N|(z -H) -1 |0, N , (24)
and

b N/2 = 0, N|(z -H) -1 |N, 0 = N, 0|(z -H) -1 |0, N . (25)
There exists an invariance of the recurrence relation under the following simultaneous sign changes in Eqs. ( 21) and ( 22)

z → -z, U → -U, a j → -a j , b j → -b j . (26)
This implies that a change from a repulsive to an attractive Hubbard interaction results in a mirror image with respect to energy of the spectral coefficients

c0 (E, U ) = c0 (-E, -U ), cN (E, U ) = cN (-E, -U ) . (27)

B. Double well with mixed interaction

In the case U 1 = -U 2 ≡ U we have one more variable, namely a k , b k and c k with the following recurrence relations (App. A)

g k+1 = a k+1 c k+1 c k+1 b k+1 , g 0 = 1 z 1 0 0 1 (k = 0, 1, ..., N/2 -1) (28)
with matrix elements (n = N/2):

a k+1 = z + 4Un(k + 1) -J 2 (n + k + 1)(n -k)b k D k+1 (29) b k+1 = z -4Un(k + 1) -J 2 (n + k + 1)(n -k)a k D k+1 (30)
c k+1 = - J 2 (n + k + 1)(n -k)c k D k+1 (31)
and with

D k+1 = [z -4Un(k + 1) -J 2 a k (n + k + 1)(n -k)][z + 4Un(k + 1) -J 2 b k (n + k + 1)(n -k)] -J 4 c 2 k (n + k + 1) 2 (n -k) 2 .
The final result of the iteration is

g N/2 = a N/2 b N/2 b N/2 c N/2 , (32)
with

a N/2 = N, 0|(z -H) -1 |N, 0 , b N/2 = 0, N|(z -H) -1 |0, N , (33)
c N/2 = N, 0|(z -H) -1 |0, N = 0, N|(z -H) -1 |N, 0 . (34)
V. RESULTS

The iteration of Eqs. ( 21), ( 22) for a symmetric double well and the iteration of Eqs. ( 29)-(31) for a double well with mixed interaction gives us, according to Eqs. ( 24), ( 25) and ( 33), (34), the following four matrix elements of the resolvent

N, 0|(z -H) -1 |N, 0 , 0, N|(z -H) -1 |0, N , 0, N|(z -H) -1 |N, 0 = N, 0|(z -H) -1 |0, N .
These are rational functions of z, as shown in Eq. (7). For N bosons these are lengthy expressions with N + 1 poles. Therefore, it is convenient to present the results as plots with respect to energy. Examples of the spectral coefficients c0 (E) and cN (E) are shown for a symmetric double well with 100 bosons in Fig. 1 and with 20 bosons in Fig. 2, and for a double well with mixed interaction for 100 bosons in Fig. 3. A larger number of bosons shows a richer spectral structure. The diagonal coefficient c0 (E) in the case of 100 bosons is remarkably different from the off-diagonal coefficient cN (E) because the latter does not have spectral weight from eigenstates whose energy E j is larger than the energy of the initial Fock state Ē = UN 2 . The reason for this feature is the double degeneracy of the eigenvalues mentioned in Sect. III: The signs of the product 0, N|E j E j |N, 0 for adjacent eigenvalues are opposite to each other. Since the eigenvalues get closer pairwise as we increase their energy, the contribution of the two levels cancel each other for each pair inside the sum of Eq. ( 15). This interaction effect is also visible for 20 bosons (Fig. 2), although the cancellation is incomplete then due to a larger level distance. This can be considered as an effect of spectral fragmentation, where the spectrum has a nondegenerate low-energy and a degenerate high-energy part, caused by the competition of tunneling and interaction.

The contribution of the two Fock states |0, N , |N, 0 to the evolution in Eq. ( 13) is given by the coefficients c 0 (t) = N, 0|Ψ t , c N (t) = 0, N|Ψ t . In Fig. 4 For a double well with interaction strength ±U the spectrum is different because of the absence of double degeneracy of the eigenvalues (cf. Fig. 3). There are two "bands", one around E = UN 2 , the other around E = -UN 2 , where the widths of the bands is characterized by the tunneling rate J. Moreover, the off-diagonal part cN (E) appears closer to zero energy and its values are very small. This indicates that the off-diagonal part has overlaps with energy levels which are different from those of the diagonal part c0 (E) . For the evolution only the latter contribute substantially, preventing the system to create a N00N state.

VI. DISCUSSION AND CONCLUSIONS

In order to understand the evolution of an isolated many-body bosonic system, we start with noninteracting bosons (i.e. U 1 = U 2 = 0) of Sect. III. The spectral properties are characterized by (i) equidistant energy levels with distance J and (ii) a binomial weight distribution of the energy levels. The evolution of a Fock state is characterized by a periodic behavior with a single frequency ω = J/2π as a direct consequence of the equidistant energy levels. The amplitudes for visiting the initial Fock state |N, 0 or the complimentary Fock state |0, N vary with cos N (Jt) or (-i) N sin N (Jt), respectively. This implies for a large number N of bosons that (i) these states are visited only for a very short period of time and (ii) the two Fock states are visited at different times. Thus the formation of a N00N state is very unlikely for noninteracting bosons.

A simple qualitative picture for the general evolution of the Fock state is the random walk in Hilbert space. In case of noninteracting bosons the particles can walk independently of each other which enables them to explore the entire Hilbert space spanned by the Fock states without restriction. A simultaneous overlap of |Ψ t with both Fock states |N, 0 and |0, N is very unlikely then, as discussed above. Once we have turned on the boson-boson interaction the particles experience a mutual influence which restricts their individual random walks. This is related to the fact that the system stays much longer in the energetically (almost) degenerate Fock states |N, 0 and |0, N than it would for noninteracting bosons (cf. Fig. 5) and, what is even more important here, they can have a simultaneous overlap with both Fock states, such that they create a N00N state. In terms of the spectral properties the interaction modifies (i) the energy levels, which are not equally spaced, and (ii) the weight distribution of the levels, which are not binomial any longer (cf. Figs. 123). This, of course, affects also the evolution of the Fock state which is more complex now, since many different frequencies are involved. A particular feature is the spectral fragmentation (cf. Fig. 1), where only a part of the spectrum contributes to the off-diagonal coefficient cN (E). This is a kind of Hilbert-space localization, where transitions to the high-energy part of the Hilbert space are completely suppressed, similar to the self-trapping found in the Hartree approximation of the Bose-Hubbard model [13]. It should be noticed, however, that spectral fragmentation appears at a much weaker interaction strength than the self-trapping effect. For U ≈ J, which is the threshold for self-trapping [13], there is only one eigenvalue with significant weight | N, 0|E j | 2 [6]. Thus it is unlikely that the two effects are directly connected.

For the double well with mixed interaction the situation is different due to the existence of two "bands" and the absence of the double degeneracy. The main consequence is the absence of a support for the formation of N00N states because the off-diagonal coefficient c N (t) is strongly suppressed. From this observation we can conclude that the evolutionary entanglement is much more favorable in the symmetric double well. This is in qualitative agreement with the results of a multiconfigurational Hartree calculation for a one-dimensional Bose gas with Gaussian barrier in Ref. [20].

In conclusion, we have studied the evolution of a bosonic Fock state |N, 0 in a double well and found that a local particle-particle interaction supports the formation of a N00N state, provided that the interaction is not too strong. This is accompanied by a fragmentation of the spectrum. The latter is characterized by the fact that only eigenstates with energies less than the energy of the initial Fock state can be reached in the evolution. This interaction effect causes a Hilbert-space localization and prevents the evolution of the Fock state to disappear in the depth of the Hilbert space. This is the main reason for a favorable creation of a N00N state. The appearance of a N00N state is suppressed though for strong interaction because then the restriction of the Hilbert space is too severe and does not allow to reach the complementary Fock state |0, N . 14) and ( 15) for 100 bosons with U/J ≈ 0.023 and = 0.01. The energy of the initial Fock state is Ē = 1000. The spectral fragmentation appears around E ≈ 970, where the levels are nondegenerate at lower energies but almost degenerate for higher energies. This is a consequence of the competition between tunneling and interaction, which the latter wins at higher energies.

+ U 1 (N/2 -k) 2 + U 2 (N/2 + k) 2 |N/2 -k, N/2 + k N/2 -k, N/2 + k| .
The off-diagonal terms of the effective Hamiltonian in Eq. (A4) read

P 2(n-k) HP 2(n-k+1) = P 2N -2k HP N -2k+2 = -J N/2 + k N/2 -k + 1 (|N/2 + k, N/2 -k N/2 + k -1, N/2 -k + 1| + |N/2 -k, N/2 + k N/2 -k + 1, N/2 + k -1|) .
This leads for U 1 = U 2 to Eqs. (21), ( 22) and for U 1 = -U 2 to Eqs. (29), (30) and (31).

5 .√ 15 .

 515 the real parts of these coefficients are plotted for 100 bosons. Their evolution indicates a collapse and revival behavior. The latter is mainly due to a rapidly changing phase factor, since |c 0 (t)|, |c N (t)| has a much smoother behavior. This allows us to plot the evolution of the two-dimensional vector (|c 0 (t)|, |c N (t)|) to characterize the dynamics of the N00N state. Examples of the complex dynamical behavior are shown for 100 bosons and for 20 bosons in Fig. These results suggest a statistical description with a probability P (|c 0 (t)|, |c N (t)|) that measures how often certain values of |c 0 (t)|, |c N (t)| are visited during the evolution in a period of time. The result for 20 bosons is plotted in Fig. 6 for U = 0.1, J = It indicates that there is a strong correlation between the coefficients, where the most favored values are |c 0 (t)| ≈ |c N (t)| ≈ 0.35.

FIG. 1 :

 1 FIG.1: Spectral coefficients of Eqs. (14) and(15) for 100 bosons with U/J ≈ 0.023 and = 0.01. The energy of the initial Fock state is Ē = 1000. The spectral fragmentation appears around E ≈ 970, where the levels are nondegenerate at lower energies but almost degenerate for higher energies. This is a consequence of the competition between tunneling and interaction, which the latter wins at higher energies.

FIG. 2 :FIG. 3 :FIG. 4 :FIG. 5 :

 2345 FIG.2: Spectral coefficients for 20 bosons with U/J ≈ 0.1 and = 0.01. The energy of the initial Fock state is Ē = 40. The almost degenerate states appear above 42.

FIG. 6 :

 6 FIG. 6: DistributionP (|c0(t)|, |cN (t)|) of |c0(t)|, cN (t)|over a time period of 0.4h/U for 20 bosons with U/J ≈ 0.026. The axes are scaled by a factor 100 and the vertical axis is in arbitrary units. This plot indicates a strong correlation between the two spectral coefficients, supporting the formation of a N00N state.

Acknowledgments

I am grateful to A. Streltsov for discussing his work on the N00N state.

APPENDIX A: RECURSIVE PROJECTION METHOD

Given is a sequence of projectors P j (n ≥ j ≥ 0), defined by the recurrence relation

with initial conditions P -1 = 1, P 0 and by the Hamiltonian H through the properties

The projection of the resolvent (z -H) -1 defines

where (...)

-1

2(n-k) is the inverse on the P 2(n-k) -projected Hilbert space. Then g k satisfies the recurrence relation

Of interest is here only the case k = n, where we have from Eq. (A2)

For the specific case of the double well we choose n = N/2 and the projectors

With the Hubbard Hamiltonian of Eq. ( 4) the diagonal terms of the effective Hamiltonian in Eq. (A4) read


K Berrada

H Eleuch

Y Hassouni

Asymptotic Dynamics of Quantum Discord in Open Quantum Systems

Keywords: numbers: 03.67.-a, 03.67.Mn, 03.65.Ud

published or not. The documents may come

Asymptotic dynamics of quantum discord in open quantum systems

I. INTRODUCTION

In the field of quantum information it is important to distinguish between quantum and classical aspects of correlation in a composite quantum state. It is well known that many operations in various quantum information processing tasks depend largely on special kind of quantum correlation, that is, entanglement. Much work has been performed in order to subdivide quantum states into separable and entangled states [START_REF] Nielsen | Quantum Computation and Information[END_REF][START_REF] Alber | Quantum Information[END_REF][START_REF] Horodecki | [END_REF]. However, there are other nonclassical correlations apart from entanglement [4][5][6][7] that can be of great importance to this field and some of these have been verified experimentally [8,9]. These correlations are more general and more fundamental than entanglement. Therefore, it is important to study, characterize and quantify quantum and classical correlations. Several measures of these quantum correlations have been investigated in the literature [10][11][12], and among them the quantum discord [10] has recently received a great deal of attention [13][14][15][16][17][18][19]. Quantum discord has been proposed as the key resource present in certain quantum communication tasks and quantum computational models without containing much entanglement [7,20,21]. The quantum discord quantifies nonclassical correlations of more general and more fundamental type than entanglement because separable mixed states (without entanglement) can have nonzero quantum discord. This indicates that classical communication can give rise to quantum correlations due to the existence of nonorthogonal quantum states. For pure entangled states quantum discord coincides with the entanglement entropy. However, for two-qubit mixed states, the relation between the classical and quantum correlations is complicated and not yet clear.

For understanding and distinguishing between classical and quantum correlations, it is important to assume that the total correlation is a directly sum of the both correlations. For a bipartite system ρ ab , it is largely accepted that quantum information mutual measures its total correlation (classical and quantum correlation) defined as [22][23][24][25] I ρ ab := S (ρ a ) + S ρ b -S ρ ab ,

where ρ a and ρ b are the reduced density matrices of the bipartite state ρ ab and S (ρ) = -tr (ρ log 2 ρ) is the von Neumann entropy. However, the mutual quantum information may be written as a sum of classical correlation C ρ ab and quantum correlation, that is, quantum discord C ρ ab as

I ρ ab := C ρ ab + Q ρ ab . (2)
In order to quantify quantum discord contained in the bipartite state, Ollivier and Zurek [10] have proposed the use of von Neumann type measurements, which consists of a set of one-dimensional projectors that sum up to the identity. The projective measurements on subsystem remove all nonclassical correlations between the parts, i.e., after a measurement on particular subsystem, all quantum correlations are destroyed. If projector measurements {Λ j } are performed locally only on the subsystem b, then the conditional density operator ρ j with the measurement result j is

ρ ab j = 1 p j (I ⊗ Λ j ) ρ ab (I ⊗ Λ j ) ( 3 )
where p j = tr (I ⊗ Λ j ) ρ ab (I ⊗ Λ j ) is the probability for the measurement of the jth state and I is the identity operator of subsystem a. However, the density matrix Confidential: not for distribution. Submitted to IOP Publishing for peer review 30 May 2011 ρ ab j is considered as a conditional density operator and we can define the quantum conditional entropy with the respect of this measurement [25] S ρ ab |{Λ j } := j p j S (ρ j ) ,

and the related quantum mutual information of this measurement is defined as

I ρ ab |{Λ j } := S (ρ a ) -S ρ ab |{Λ j } . (5)
To derive classical correlation, since I ρ ab |{Λ j } depends on the projector operators {Λ j }, we take the maximum of I ρ ab |{Λ j } , taking into the all possible projectors. On other word, the classical part is defined as the maximal information about subsystem a that can be obtained by performing a measurement on the subsystem b. The classical correlation between the two subsystems a and b is defined as

C ρ ab := sup {Λj } I ρ ab |{Λ j } . ( 6
) = S (ρ a ) -min {Λj } S ρ ab |{Λ j } .
In this scenario, the quantum discord that measures the nonclassical correlations between two subsystem, can be introduced as a difference between the original quantum mutual information in Eq. ( 1) and the the classical correlation in Eq. ( 6) as

Q ρ ab := I ρ ab -C ρ ab , (7)
which is always non-negative quantity by expressing mutual information in terms of quantum relative entropy and invoking the monotonicity property of the latter [26,27]. The zero-discord states are relatively well study: Q ρ ab = 0 if and only if there exists a complete orthonormal basis {|k } for the subsystem a and some density operator ρ b for the subsystem b such that

ρ ab = k p k |k k| ⊗ ρ b .
Recently various methods to detect zero-discord [28][29][30] have been proposed for a given state as well as for an unknown state [31]. Moreover, it is found that vanishing quantum discord is related to the complete positivity of a map [32] and the local broadcasting of quantum correlations [15]. Besides the characterization and quantification of classical and quantum correlations, another important problem is the behavior of these correlations under the action of decoherence. Every natural object is in contact with its environment, so its dynamics is that of an open system. Thus, the interaction between a composite quantum system and its environment and understanding the dynamics of different kinds of correlations has attracted more interest. The problem of the formulation and characterization of the dynamics correlations of open systems in the quantum regime has a long and extensive history. The quantum entanglement dynamics in open quantum systems was broadly investigated in the literature. However few works dealt with the effect of the environment on quantum discord [33][34][35][36][37].

On the other hand, the definition of quantum discord in Eq. ( 7) exhibits an optimization problem. However, it is usually intractable to determine the quantum discord for generic cases. Even in the simplest bipartite quantum states, an analytical expression for general two-qubit state is still missing. The dynamics of the classical and quantum correlations in the presence of an environment has been studied for a certain class of high symmetrical states [33,34,38]. For this reason the relation between the dynamics of the different kinds of correlations is not known for more general quantum system.

More recently M. Ali et al. [39] have developed an analytic method to evaluate both classical correlations and quantum discord for the complete set of two-qubit X states depending on seven real-valued parameters in order to deepen the understanding the relation between different kinds of correlations. Using this result, we present a detailed analysis of the time evolution of classical and quantum correlations for a large class of two-qubit states in the case of a Lindblad-type master equation. This dynamics that exhibits a rich manifold of asymptotic states that may be more or less classical and quantum correlated with respect to the initial states they emerge from, including maximally entangled Bell states, separable, and nonseparable states. Analytical expressions are calculated and numerically represented, which exhibit the dependence of the classical correlation and quantum discord for asymptotic states. However, it is found that quantum discord can be larger than classical correlation and it can persist even in the asymptotic long-time regime. Furthermore, the comparison between the dynamics of quantum discord and quantum entanglement shown that quantum discord behave differently under the effect of environment. Our observations may have important implications in exploiting these correlations in quantum information processing and transmission.

The paper is organized as follows. In Sec. II, we study the open dynamics of classical and quantum correlations of two qubits in the present model. Analytic derivation of the expression of different kinds of correlations are presented. In Sec. III, we apply this result for a large family of two-qubit states, studying the relation between the dynamics of different kinds of correlations. Finally, we conclude this work in Sec. IV.

II. DYNAMICS OF THE CLASSICAL AND QUANTUM CORRELATION IN OPEN BIPARTITE SYSTEMS

In this section, we describe the asymptotic dynamics of two-qubit systems and we derive the classical correlation and quantum discord for them.

Let us consider a bipartite system composed of two qubits immersed in external environment via standard weak-coupling limit techniques [START_REF] Alicki | Quantum Dynamical Semigroups and Applications[END_REF]. The reduced irreversible dynamics is described by one-parameter semigroup of linear maps, which called quantum dynamical semigroups obtained from γ(t) = exp(tL)), and it is gov-erned by the master equation: ∂ t ρ(t) = L[ρ(t)] where the generator L takes care the effects of the environment through the elements of the Kossakowski matrix [START_REF] Alicki | Quantum Dynamical Semigroups and Applications[END_REF][START_REF] Breuer | The Theory of Open Quantum Systems[END_REF][START_REF] Benatti | [END_REF]. Formally, we have

∂ t ρ(t) = L[ρ(t)] = -i Ω 2 Σ 3 , ρ(t) + 3 i,j=1 A ij Σ i ρ(t)Σ j - 1 2 Σ j Σ i , ρ(t) (8)
where Ω define the system frequency, Σ i := σ i ⊗ I + I ⊗ σ i is the 2 × 2 identity matrix, σ i are the Pauli matrices and the matrix

A = [A ij ] =   1 iα 0 -iα 1 0 0 0 1   α ∈ R, α 2 ≤ 1, (9)
is positive semi-definite (elements of the Kossakowski matrix) [START_REF] Benatti | [END_REF]. This shows that the semigroup generated by the master equation consists of completely positive maps γ(t) for all t ≥ 0.

By means of the single qubit Pauli matrices σ

(1) i = σ i ⊗I and σ (2) i
= I⊗σ i , we can write the purely dissipative contribution to the generator

D[ρ(t)] = 3 i,j=1 A ij 2 a,b=1 σ (a) i ρ(t)σ (b) j - 1 2 {σ (b) j σ (a)
i , ρ(t)} .

(10) In this way there are six Kraus operators σ (a)

i , a = 1, 2; i = 1, 2, 3 and the 6 × 6 Kossakowski matrix reads

K = K (ab) ij = K (11) K (12) K (21) K (22) = A A A A . (11)
From the theory of open quantum systems [START_REF] Alicki | Quantum Dynamical Semigroups and Applications[END_REF][START_REF] Breuer | The Theory of Open Quantum Systems[END_REF]44], the coefficients

K (ab) ij
in the Kossakowski matrix relative to the i-th Pauli matrix of the a-th qubit, respectively the j-th Pauli matrix of the b-th qubit, a, b = 1, 2; i, j = 1, 2, 3, are determined by the Fourier transforms of the two-point time-correlation functions with respect to an environment equilibrium state ω, ω B

(a) i B (b) j (t) , of the environment operators B (a) i ap- pearing in the system-environment interaction H I = 3 i=1 σ (1) i ⊗ B (1) i + σ (2) i ⊗ B (2) i
. The symmetric form of (3) thus results when both qubits are linearly coupled to bath operators such that:

B (1) 1,2,3 = B (2) 1,2,3 = B 1,2,3 and ω (B 1,2 B 3 (t)) = 0.
In the following, we are interested to the time evolution of two-qubit system prepared initially in the form

ρ = a|1 1| + d|2 2| + b|3 3| + c|4 4| + e|1 2| + f |2 1|. ( 12
)
where a, b, c and d are real constants satisfying the positive conditions and unit trace, that is, a

+ b + c + d = 1,
|1 = |00 , |2 = |11 , |3 = |01 + |10 √ 2 , |4 = |01 -|10 √ 2 .
(13) The density matrix ( 12) is represented in standard basis as

ρ =     a 0 0 e 0 b+c 2 b-c 2 0 0 b-c 2 b+c 2 0 f 0 0 d     . (14)
Equation ( 14) describes a so-called quantum X state that is entangled if and only if either (a+b) 2

4

< e 2 or ad < (b-c) 2

4

. Both conditions cannot hold simultaneously. The density operator ( 14) includes a very large family of two-qubit states, including Bell states, Werner states, Horodecki states, etc.

The evolution of the density operator ρ(t) initially prepared in Eq. ( 12) can be obtained as

ρ(t) =a(t)|1 1| + d(t)|2 2| + b(t)|3 3| + c(t)|4 4| + e(t)|1 2| + f (t)|2 1|. ( 15
)
where

a(t) = (1 -α) 2 3 + α 2 R + 2(1 + α)a -(1 -α) 2 (b + d) 3 + α 2 E + (t) + 1 -α 2 (1 + α) 2 a -2 (1 -α) d + (1 + α) 2 b (1 + α) (3 + α 2 ) E -(t) d(t) = (1 + α) 2 3 + α 2 R - (1 + α) 2 a -2(1 + α)d + (1 + α) 2 b 3 + α 2 E + (t) -1 -α 2 2 (1 + α) a -(1 -α) 2 (b + d) (1 -α) (3 + α 2 ) E -(t) b(t) = (1 -α) 2 3 + α 2 R + 2(1 + α 2 )b -(1 -α 2 )(a + d) 3 + α 2 E + (t) + 1 -α 2 (1 + α) 3 a + (1 -α) 3 d -2 1 -α 2 b (1 -α 2 ) (3 + α 2 ) E -(t) c(t) = c e(t) = e exp(-12t) f (t) = f exp(-12t), (16)
with

R = a + b + d, E + (t) = e -8t cosh 4t 1 -α 2 , E -(t) = e -8t sinh 4t 1 -α 2 . ( 17
)
It is known that in certain specific situations, coupling an environment need not have destroy correlations between subsystems, but can create them, even create correlations that persist at arbitrary long times. This possibility depends on the form of the generator of the reduced dynamics and initial states. Furthermore, in Ref. [43] a necessary and sufficient condition has been obtained for environmental-induced entanglement in an initially separable states for two-qubit system. Starting from an initially separable states, the correlations generated at short times can persist asymptotically, and also starting from initially entangled states, its entanglement content asymptotically increase. This can happen in the present case and the correlation generation capability of the environment is due to the Kossakowski-generator of the reduced dynamics (8). Now, our primary aim is to evaluate the quantum discord dynamics Q (ρ(t)) defined by Eq. ( 7) corresponding to above-mentioned cases. This requires us to derive the quantum mutual information dynamics I (ρ(t)) defined by Eq. ( 1) and classical correlation C (ρ(t)) defined by Eq. ( 6).

Let us start by evaluating the total correlations, that is quantum mutual information I (ρ(t)), contained in the bipartite state (15). From Eq. ( 15), the eigenvalues of the density matrix ρ(t) are given by

λ 1 (t) = 1 2 (a(t) + d(t)) + (a(t) -d(t)) 2 + 4e(t)f (t) λ 2 (t) = 1 2 (a(t) + d(t)) -(a(t) -d(t)) 2 + 4e(t)f (t) λ 3 (t) = b(t) λ 4 (t) = c(t). (18)
The mutual quantum information is given by

I(ρ(t)) = S(ρ a (t)) + S(ρ b (t)) + 4 i=1 λ i (t) log 2 λ i (t), ( 19
)
where the quantities S(ρ a (t)) and S(ρ b (t)) are the marginal entropies of density operator ρ(t). They have the following form

S(ρ a (t)) = - 1 + a(t) -d(t) 2 log 2 (1 + a(t) -d(t)) - 1 + d(t) -a(t) 2 log 2 (1 + d(t) -a(t)) = S(ρ b (t)). ( 20
)
After evaluating the total correlations I (ρ(t)), we shall compute the evolution classical correlation C (ρ(t)) using von Neumann type measurements. It is known that any von Neumann measurement for subsystem b can be written as [25] Λ

j = V Π j V † j = 0, 1, (21)
where Π j = |j j| is the local measurement for qubit b along the computational base {|j } and V ∈ SU (2) is an unitary operator with unit determinant. After the measurement {Λ j }, the bipartite state ρ(t) will change to the ensemble {p j (t), ρ j (t)}, that is 2×2 density matrices, where

ρ j (t) := 1 p j (t) (I ⊗ Λ j ) ρ(t) (I ⊗ Λ j ) (22)
and p j (t) = tr ((I ⊗ Λ j ) ρ(t) (I ⊗ Λ j )) . The operator V can be written, up to a constant phase, as

V = tI + i - → y - → .σ. ( 23
)
With t, y 1 , y 2 , y 3 ∈ R and t 2 + y 2 1 + y 2 2 + y 2 3 = 1. This implies that these parameters, three among them independent, assuming their values in the interval

[-1, 1], that is, t, y i ∈ [-1, 1].
In order to obtain the evolution of classical correlation, We need to evaluate ρ j (t) and p j (t). In this way, Eq. (A1) can be rewritten as

p j (t)ρ j (t) = I ⊗ V Π j V † ρ(t) I ⊗ V Π j V † (24)
and the minimization over von Neumann measurements min {Λj } [S (ρ(t)|{Λ j })] is straight forwardly derived by the analytically method in Ref. [39] and is given as

min {Λj } [S (ρ(t)|{Λ j })] = min f 1 2 ,m (t) , f 1,m (t) . ( 25
)
This result, which we prove in Appendix A, helps to explicitly solve the minimization problem. Then, the evolution of the classical correlation can be expressed as

C(ρ(t)) = S(ρ a (t)) -min f 1 2 ,m (t) , f 1,m (t) . (26)
Therefore, the quantum discord dynamics can be written as

Q(ρ(t)) = I(ρ(t)) -S(ρ a (t)) + min f 1 2 ,m (t) , f 1,m (t) ( 27
)
where the amount of the total correlation I(ρ(t)) is given in Eq. (19).

The above expression exhibits the time evolution of the quantum discord for open bipartite quantum system characterized by a quantum dynamical one-parameter semigroup of completely positive linear maps.

In order to compare the quantum discord dynamics with the entanglement dynamics in the present model, we use the concurrence as a measure of entanglement. For two-qubit states, the concurrence dynamics can be given as [45] C (ρ(t)) = max{θ 1 (t) -θ 2 (t) -θ 3 (t) -θ 4 (t), 0}, (28) where θ i (i = 1, 2, 3, 4) are the square roots of the eigenvalues of the non-Hermitian matrix ρ(t)ρ(t) in descending order. Note that each θ i (t) is a positive real number and ρ(t) being the spin-flip operation of ρ(t) given by ρ(t) = (σ y ⊗ σ y )ρ * (t)(σ y ⊗ σ y ) where ρ * (t) denotes the complex conjugate of ρ(t). For a two-qubit with a structure defined as in Eq. ( 15) we have

ρ(t) =      a(t) 0 0 e(t) 0 (b(t)+c(t) 2) (b(t)-c(t) 2) 0 0 (b(t)-c(t) 2) (b(t)+c(t) 2) 0 f (t) 0 0 d(t)      , (29)
ρ(t)ρ(t) =     a(t)d(t) + e 2 (t) 0 0 a(t)f (t) + e(t)a(t) 0 b 2 (t)+c 2 (t) 2 b 2 (t)-c 2 (t) 2 0 0 b 2 (t)-c 2 (t) 2 b 2 (t)+c 2 (t) 2 0 f (t)d(t) + d(t)e(t) 0 0 a(t)d(t) + f 2 (t)     . (30)
and therefore the corresponding concurrence can be given as

C (ρ(t)) = max 0, Ξ 1 (ρ(t)) , Ξ 2 (ρ(t)) , (31)
with

Ξ 1 (ρ(t)) = 2e(t) -b(t) -c(t) and Ξ 2 (ρ(t)) = |b(t) -c(t)| -2 a(t)d(t).
In what follows one wants to inquire the relationship between the dynamics of the classical and quantum correlations. In particular, one wants to understand the changes in behaviour in the asymptotic dynamics of different kinds of correlations.

III. RELATION BETWEEN CLASSICAL CORRELATION, QUANTUM DISCORD AND ENTANGLEMENT IN OPEN QUANTUM SYSTEMS

In this section, we study in detail the relation between the asymptotic dynamics of the classical correlation, quantum discord, and entanglement for a large family of initial states.

To explore the influence of decoherence on the dynamical behavior of the classical and quantum correlations in the present model, we have plotted the time evolution of C(ρ(t)), Q(ρ(t)), and C(ρ(t)) as a function of time t for various values of parameter α in several classes of twoqubit states. The choice of the range of the values for the plot axes was made only for graphic reasons to make the plot clear. Changing the values of the parameter α does not alter the results of the evolution of the classical and quantum correlations.

Case 1. As a first example, we consider an initial pure separable state, that is, ρ = |1 1|. The state ρ is uncorrelated with zero classical and quantum correlations and it goes into a mixed separable state (32):

ρ(∞) = (1 -α) 2 3 + α 2 |1 1| + (1 + α) 2 3 + α 2 |2 2| + (1 -α 2 ) 3 + α 2 |3 3|. ( 32
) For this particular case, the different kinds of correlations are given as

C(ρ(∞)) = S(ρ a (∞)) -f 1 2 ,0 (∞) Q(ρ(∞)) = S(ρ a (∞)) + f 1 2 ,0 (∞) -S(ρ(∞)) C(ρ(∞)) = 0, (33)
where the entropies S(ρ a (∞)) and S(ρ(∞)) are given by Eqs. ( 19) and ( 20) respectively, and the minimal value of the function f reads

f 1 2 ,0 (∞) = - 1 -Γ(∞) 2 log 2 1 -Γ(∞) 2 - 1 + Γ(∞) 2 log 2 1 + Γ(∞) 2 (34)
where

Γ(∞) = (a(∞) -d(∞)) 2 + b 2 (∞) 1 2 . ( 35
)
Figure 1 displays the asymptotic dynamics of the classical correlation, quantum discord, and concurrence for the pure state ρ for α = 0.5. The solid line presents the dynamics of quantum discord, the dotted line is for the classical correlation. It can be seen that the dissipative time evolution is able to generate the quantum discord and classical correlation for the state ρ, contrariwise the quantum entanglement is zero at any time. Moreover, quantum discord is observed to be larger than classical correlation.

Case 2.

We consider an initial maximally entangled pure state, one of the four Bell states defined as

ρ + 1 = |φ + φ + | = 1/2 (|1 1| + |2 2| + |2 1| + |1 2|), ρ - 1 = |φ -φ -| = 1/2 (|1 1| + |2 2| -|2 1| -|1 2|), ρ + 2 = |ψ + ψ + | = |3 3|, and ρ - 2 |ψ -ψ -| = |4 4|, with |φ ± = 1/ √ 2 (|00 ± |11 ) and |ψ ± = 1/ √ 2 (|01 ± |10 ).
It is known that for any Bell states, we have

C(ρ(0)) = Q(ρ(0)) = C(ρ(0)) = 1. ( 36
)
For this example, classical correlation, quantum discord and any measure of entanglement coincide and are equal to the maximum value of the correlation between subsystems. Under the irreversible time evolution, the initial maximal pure states ρ ± 1 and ρ + 2 go into an asymptotic mixed separable state which is less classical and quantum correlated :

ρ(∞) = (1 -α) 2 3 + α 2 |1 1| + (1 + α) 2 3 + α 2 |2 2| + (1 -α 2 ) 3 + α 2 |3 3|. ( 37
)
In Fig. 2, we plot the dynamics of the correlations for the Bell states ρ ± 1 and ρ + 2 . From the figure, the dissipative time evolution shows a sudden death of quantum entanglement, that is, the concurrence vanishes at finite time. On the other hand, as we can see, contrary to what happens to entanglement, the quantum discord may remain constant for long times, that is the case when quantum coherence is alive. This implies that quantum discord is more resistant to noise than entanglement. Moreover, we see that the discord is always greater than classical correlation.

In contrast to previous maximal pure states, the classical and quantum correlations of the state ρ - 2 are unaffected under decoherence and we have

ρ - 2 (∞) = |4 4|, ( 38
)
with the corresponding correlations

C(ρ - 2 (t)) = Q(ρ - 2 (t)) = C(ρ - 2 (t)) = 1, ∀t ≥ 0. ( 39
)
This is due to rotational invariance of |4 . Case 3. We take an initial mixed separable state, that is, ρ = 1/2 (|3 3| + |4 4|), with vanishing values of quantum correlations (Q(ρ(0)) = C(ρ(0) = 0) and maximal classical correlation (C(ρ(0)) = 1, all correlations are classical). During the time evolution, the state ρ goes into an asymptotic state with more quantum correlations and less classical correlation :

ρ(∞) = 1 2 (1 -α) 2 (3 + α 2 ) |1 1| + 1 2 (1 + α) 2 (3 + α 2 ) |2 2| + 1 2 (1 -α 2 ) (3 + α 2 ) |3 3| + 1 2 |4 4|, ( 40
)
with the corresponding correlations

C(ρ(∞)) = S(ρ a (∞)) -f 1 2 ,0 (∞) Q(ρ(∞)) = S(ρ a (∞)) + f 1 2 ,0 (∞) -S(ρ(∞)) C(ρ(∞)) = 2α 2 3 + α 2 (41)
where the minimal value of the function f is obtained as

f 1 2 ,0 (∞) = - 1 -γ(∞) 2 log 2 1 -γ(∞) 2 - 1 + γ(∞) 2 log 2 1 + γ(∞) 2 (42)
with

γ(∞) = (a(∞) -d(∞)) 2 + (b(∞) -c) 2 1 2 . ( 43
)
The dynamics of the classical correlation, quantum discord, and concurrence have been plotted in Fig. 3 for ρ against t. The solid line is for quantum discord, the dashed line is for concurrence, and dotted line is for classical correlation. Interestingly, in this particular case, we note that the dissipative time evolution is able to generate the quantum correlations. However, we see that the correlations have a different order as functions of t, with classical correlation initially larger than discord and concurrence, but after some time, discord becomes larger than concurrence and classical correlation. Moreover, the all correlations remain constant as time become large.

Case 4.

Let us consider an initial mixed entangled state, that is, ρ = 1/10 (3|2 2| + |3 3| + 6|4 4|) with concurrence grater than quantum discord and classical correlation (C(ρ(0)) = 0.5, Q(ρ(0)) = 0.285829, and C(ρ(0)) = 0.286845). During the process of decoherence, the state ρ goes into an asymptotic mixed entangled state: with the corresponding classical and quantum correlations

ρ(∞) = 2 5 (1 -α) 2 (3 + α 2 ) |1 1| + 2 5 (1 + α) 2 (3 + α 2 ) |2 2| + 2 5 (1 -α 2 ) (3 + α 2 ) |3 3| + 3 5 |4 4|, (44)
C(ρ(∞)) = S(ρ a (∞)) -f 1 2 ,0 (∞) Q(ρ(∞)) = S(ρ a (∞)) + f 1 2 ,0 (∞) -S(ρ(∞)) C(ρ(∞)) = 3(1 + 3α 2 ) 5(3 + α 2 ) (45)
where ( 46)

f 1 2 ,0 (∞) = - 1 -γ(∞) 2 log 2 1 -γ(∞) 2 - 1 + γ(∞) 2 log 2 1 + γ(∞) 2 (47)
with

γ(∞) = (a(∞) -d(∞)) 2 + (b(∞) -c) 2 1 2 . (48)
In Fig. 4 we depict the graphs of the dynamics of the classical and quantum correlations for the state ρ. In contrast to previous examples, the dynamic behavior of the different kinds of correlations depends on the parameter α. As we can see in Fig. 4, the correlations of the asymptotic state can be larger or smaller than initial one depending on the parameter α. If we take for example α = 0.5, then C(ρ(∞)) < 0.286845, Q(ρ(∞)) < 0.285829, and C(ρ(∞)) < 0.5, i.e. the asymptotic state has less correlation than the initial state (Fig. 4(a)). Whereas for α = 0.9, the initial entanglement and classical correlation first diminish and then increase again, leading to an asymptotic state with more correlations than initial one. Moreover, it can be seen that for this particular initial state, quantum discord is always less than concurrence but always greater than the classical correlation (Fig. 4(b)). Case 5. We take an initial two-qubit mixed entangled state, that is, ρ = 1/10 (5|2 2| + |3 3| + 4|4 4|), with entanglement grater than quantum discord and classical correlation (C(ρ(0)) = 0.3, Q(ρ(0)) = 0.139036, and C(ρ(0)) = 0.122556). During the time evolution, the state ρ goes into an asymptotic mixed state with less or more quantum correlations depending on the parameter α :

ρ(∞) = 3 5 (1 -α) 2 (3 + α 2 ) |1 1| + 3 5 (1 + α) 2 (3 + α 2 ) |2 2| + 3 5 (1 -α 2 ) (3 + α 2 ) |3 3| + 2 5 |4 4|. ( 49
)
with the corresponding correlations

C(ρ(∞)) = S(ρ a (∞)) -f 1 2 ,0 (∞) Q(ρ(∞)) = S(ρ a (∞)) + f 1 2 ,0 (∞) -S(ρ(∞))    C(ρ(∞)) = 0 for α 2 ≤ 3 11 C(ρ(∞)) = (11α 2 -3) 5(3+α 2 ) for 3 11 < α 2 ≤ 1, (50)
where

f 1 2 ,0 (∞) = - 1 -γ(∞) 2 log 2 1 -γ(∞) 2 - 1 + γ(∞) 2 log 2 1 + γ(∞) 2 (51) with γ(∞) = (a(∞) -d(∞)) 2 + (b(∞) -c) 2 1 2 . (52)
In Fig. 5, we display the dynamics of the classical correlation, quantum discord, and concurrence versus t for different values of α. Depending on the choice of α, we can see that the correlations go into a state with less or more quantum correlations than initial one. On other hand, the dissipative time evolution shows a sudden death of the quantum entanglement, that is, the concurrence vanishes at finite time, and also it exhibits a sudden birth (Fig. 5(b)). Contrary to what happens to entanglement, the discord remains constant for a long time. More interestingly, the quantum discord dominates the classical correlation in this limit.

Case 6. we consider an initial state, that is, Werner state [46]

ρ = x|Ψ -Ψ -| + (1 -x) 4 I (53)
where

|Ψ -= 1/ √ 2 (|01 -|10
) is a maximally entangled state and 0 ≤ x ≤ 1.

The parameters defined in the density matrix ( 14) are given as

a = b = d = 1 -x 4 c = 1 + 3x 4 e = f = 0. ( 54
)
The initial Werner state, has a particular property is that the minimization of the function f does not depend on the parameters k and m, and we have

C(ρ(0)) = 1 -x 2 log 2 (1 -x) + 1 -x 2 log 2 (1 -x) Q(ρ(0)) = 1 4 [(1 -x) log 2 (1 -x) + (1 + 3x) log 2 (1 + 3x) -2(1 + x) log 2 (1 + x)] C(ρ(0)) = max{0, 3x -1 2 }. (55)
We display these correlations versus x for this state in Fig. 6, we can see that the correlations have a different order depending the parameter x with quantum discord initially larger than concurrence and the classical correlation when x 0.5234, but for x > 0.5234 the concurrence becomes larger than quantum discord and the classical correlation. During the process of decoherence, the initial state goes into an asymptotic mixed state with more correlations depending on the parameters x and α:

ρ(∞) = 3 4 (1 + x)(1 -α) 2 (3 + α 2 ) |1 1| + 3 4 (1 + x)(1 + α) 2 (3 + α 2 ) |2 2| + 3 4 (1 + x)(1 -α 2 ) (3 + α 2 ) |3 3| + (1 + 3x) 4 |4 4|, (56)
with the corresponding classical and quantum correlations

C(ρ(∞)) = S(ρ a (∞)) -f 1 2 ,0 (∞) Q(ρ(∞)) = S(ρ a (∞)) + f 1 2 ,0 (∞) -S(ρ(∞))
     C(ρ(∞)) = 0 for α 2 ≤ 3-9x 5-3x C(ρ(∞)) = (5-3x)α 2 +9x-3 2(3+α 2 ) for 3-9x 5-3x < α 2 ≤ 1, (57)
where the minimal value of the function f is given by

f 1 2 ,0 (∞) = - 1 -γ(∞) 2 log 2 1 -γ(∞) 2 - 1 + γ(∞) 2 log 2 1 + γ(∞) 2 (58)
with

γ(∞) = (a(∞) -d(∞)) 2 + (b(∞) -c) 2 1 2 . ( 59
)
In Fig. 7, we display the variation of different types of correlations against x for different values of t. From the figure it can be seen that for this particular initial state, the behavior of the correlations is the same as the initial ones. The correlations have different order depending on the values of the parameter α and x.

In Fig. 8, we show the time evolution of classical and quantum correlations for the initial Werner states with different values of the parameters α and x. We can see that the quantum discord can be smaller or larger than quantum entanglement in the asymptotic limit. Moreover, the quantum discord is always greater than the classical correlation.

Case 7. Finally we consider an initial Horodecki state [START_REF] Horodecki | Quantum Information: An Introduction to Basic Theoretical Concepts and Experiments[END_REF] defined as a mixture of a maximally entangled state, say the Bell state |Ψ + = 1/ √ 2 (|00 + |11 ), and a separable state orthogonal to it, see |00 : the matrix ( 14) are given as

ρ = x|Ψ + Ψ + | + (1 -x)|00 00|, ( 60
a = 1 -x b = x c = d = e = f = 0. (61)
In this case, the minimum of the function f is attained

for k = 1 2 , that is, f 1 2 ,0 (0) = f 1 2 , 1 4 (0)
. Therefore, The different correlations are given as

C(ρ(0)) = S(ρ a (0)) -f 1 2 ,0 (0) Q(ρ(0)) = S(ρ a (0)) + f 1 2 ,0 (0) -S(ρ(0)) C(ρ(0)) = x (62)
where

S(ρ a (0)) = - x 2 log 2 x 2 - 2 -x 2 log 2 2 -x 2 f 1 2 ,0 (0) = - 1 -x 2 + (1 -x) 2 2 log 2 1 -x 2 + (1 -x) 2 2 - 1 + x 2 + (1 -x) 2 2 log 2 1 + x 2 + (1 -x) 2 2 S(ρ(0)) = -x log 2 x -(1 -x) log(1 -x). (63)
In Fig. 9, we plot these correlations for various values of x. We can see that for the initial Horodecki state, the quantum discord is always greater than classical correlation but always small than the concurrence. During the time evolution, the initial Horodecki state goes into an asymptotic mixed state with less correlations than initial one: We display the classical correlation, quantum discord, and concurrence versus t for this state in Fig. 10. Similarly to the case 5, we can see that under decoherence the concurrence exhibits sudden death, whereas the quantum discord remains constant for asymptotic time. Moreover, the discord is always larger than the classical correlation.

ρ(∞) = (1 -α) 2 3 + α 2 |1 1| + (1 + α) 2 3 + α 2 |2 2| + (1 -α 2 ) 3 + α 2 |3 3|. ( 64

IV. CONCLUDING REMARKS

In conclusion, we have studied the classical and quantum correlation rates in open quantum systems. We have focused on the dynamics of classical correlation, quantum discord, and entanglement of two-qubit system under dynamical semigroup generators of Kossakowskitype, which exhibits a rich manifold of asymptotic states that may be more or less correlated. Despite decoherence, the presence of an environment need not have only destructive effects in relation to classical and quantum correlations; the correlations can even be asymptotically increased with respect to the initial amount. This can be happen in the present case and the correlations generation capability of the environment is due to the non-Hamiltonian contribution (10) to the generator (8). Indeed, the two-qubit Hamiltonian does not contain coupling terms and can not be a source of correlations. this can be true for (10) because the off-diagonal contributions in the Kossakowski matrix (11) couple the two qubits. This is only necessary, but not sufficient to ensure correlations generation and their asymptotic persistence. They indeed depend on a trade-off between the off-diagonal couplings and the purely decohering diagonal terms in (11).

We have evaluated the different kinds of correlations and obtained analytical expressions. The results are used to illustrate the relation between the dynamics of correlations in two-qubit sates including maximally states or partially entangled states and mixed states. We have shown that the quantum discord is more resistant to the action of the environment than the quantum entanglement and it can persist even in the asymptotic longtime regime. On the other hand, the discord is always larger than the classical correlation in this limit. Hence, the results provide the description of the behavior of the classical and quantum correlations and their robustness under the decoherence process. In fact, the quantum discord behavior between the two-qubit depends on the minimization of the classical correlations during the time evolution. This means that for initial conditions and non-Hamiltonian contribution, we expect that von Neumann measurements quantifying the classical correlations are subjected to an initial change before persisting unchanged later.

We have also noticed that, during the dissipative time evolution even without quantum entanglement, the correlations introduced by the environment are transferred to the two-qubit, producing a finite quantum discord. In fact, the quantum states with nonzero-discord are much more common than entangled states and therefore they will be easier to generate that entanglement. Furthermore, we have observed that, the classical and quantum correlations can remain unaffected for certain classes of states.

  and the set {|1 , |2 , |3 , |4 } form an orthonormal basis related to two-qubit standard basis {|00 , |01 , |10 , |11 } by

FIG. 1 .

 1 FIG. 1. Asymptotic classical and quantum correlations. Quantum discord (solid line) , Classical correlation (dotted line), and C = 0 with α = 0.5.

FIG. 2 .

 2 FIG. 2. Asymptotic classical and quantum correlations for Bell states. Quantum discord (solid line), classical correlation (dotted line), and concurrence (dashed line) with α = 0.5. (a) for ρ ± 1 and (b) for ρ + 2 .

FIG. 3 .

 3 FIG.3. Asymptotic classical and quantum correlations. Quantum discord (solid line), classical correlation (dotted line), and concurrence (dashed line) with α = 0.5.

FIG. 4 .

 4 FIG. 4. Asymptotic classical and quantum correlations. Quantum discord (solid line), classical correlation (dotted line), and concurrence (dashed line). (a) for α = 0.5 and (b) for α = 0.9.

FIG. 5 .

 5 FIG. 5. Asymptotic classical and quantum correlations. Quantum discord (solid line), classical correlation (dotted line), and concurrence (dashed line). (a) for α = 0.5 and (b) for α = -0.9.

FIG. 6 .

 6 FIG.6. Classical and quantum correlations for the initial Werner state. Quantum discord (solid line), classical correlation (dotted line), and concurrence (dashed line).

FIG. 7 .

 7 FIG. 7. Asymptotic classical and quantum correlations for the initialWerner state with (a) α = 0.5; t = 0.2, (b) α = 0.5; t = 1, (c) α = 0.9; t = 0.2, and (d) α = 0.9; t = 1. Quantum discord (solid line), classical correlation (dotted line), and concurrence (dashed line).

FIG. 8 .

 8 FIG. 8. Asymptotic classical and quantum correlations for the initialWerner state with (a) α = 0.5; x = 0.3, (b) α = 0.5; x = 0.5, (c) α = 0.9; x = 0.3, and (d) α = 0.9; x = 0.5. Quantum discord (solid line), classical correlation (dotted line), and concurrence (dashed line).

FIG. 9 .

 9 FIG. 9. Classical and quantum correlations for the initial Horodecki state. Quantum discord (solid line) , classical correlation (dotted line), and concurrence (dashed line).

FIG. 10 .

 10 FIG. 10. Asymptotic classical and quantum correlations for the initialHorodecki state. Quantum discord (solid line) , classical correlation (dotted line), and concurrence (dashed line) with α = 0.5 and x = 0.5.

Acknowledgments

The authors would like to thank Prof. Dr. Gernot Alber, Prof. Fabio Benatti and Dr. Joseph Renes for helpful discussions. K.B. pleased to express his sincere gratitude for the hospitality at the Institut für Angewandte Physik at TU Darmstadt. The authors wish to thank the referees for their valuable comments that resulted in improvements of the paper in many aspects.

APPENDIX A: CALCULATIONAL DETAILS OF THE CLASSICAL CORRELATION

After the measurement {Λ j }, the bipartite state ρ(t) will change to the ensemble {p j (t), ρ j (t)}, that is 2 × 2 density matrices, where ρ j (t) := 1 p j (t) (I ⊗ Λ j ) ρ(t) (I ⊗ Λ j ) (A1) with p j (t) = tr ((I ⊗ Λ j ) ρ(t) (I ⊗ Λ j )) .

(A2)

In order to obtain the evolution of classical correlation, We need to evaluate ρ j (t) and p j (t). Using the von Neumann measurement Λ j = V Π j V † , we write

Through some calculations, we can obtain these density matrices

where

According to Eqs. (A4) and (A5), it is straightforward to obtain the eigenvalues characterizing the ensembles {p j (t), ρ j (t)} as follows

where we have defined γ(t) and δ(t) as

and we have introduced the coefficients k, l and m as

such that

The corresponding probabilities are given as

The entropies of the ensembles {p k (t), ρ k (t)} are written as

(A16) Thus substituting Eqs. (A15) and (A16) into the quantum conditional entropy expression in Eq. ( 4), we obtain S (ρ(t)|{Λ j }) = p 0 (t)S(ρ 0 (t)) + p 1 (t)S(ρ 1 (t)), (A17) which leads to the classical correlation

(A18) The minimal value of the f (k, m) = S (ρ(t)|{Λ j }) can be easily done analytically by appropriate choice of the parameter k, l, and m:

then the evolution of the classical correlation can be expressed as

APPENDIX B: CLASSICAL AND QUANTUM CORRELATIONS FOR ASYMPTOTIC STATES

In this Appendix, we give some details of the calculations of classical and quantum correlations for asymptotic states.

When the time tends to infinity, the coefficients of the asymptotic states resulting from the inial density matrix (12) are given by lim t→∞

which leads to the following asymptotic states

There is a one-parameter set {ρ(∞)} 0≤c≤1 of asymptotic states, where the initial states of the form (12) with the same c go into the same ρ(∞). The corresponding concurrence is given by

(B6) Evaluation of the classical correlation given by Eq. (A18) can obtained by taking the minimal value of the function f , that is, min

, and S (ρ 0 (∞)) = S (ρ 1 (∞)) and the minimization of f is equal to the minimization of either S (ρ 0 (∞)) or S (ρ 1 (∞)) and is given by

where

For k = 1, we have

where

and

(B14) Based on this result, we are able to evaluate the classical correlation and quantum discord. We note that f 1 2 ,0 (∞) ≤ f 1,0 (∞), that is, min f 1 2 ,0 (∞) , f 1,0 (∞) = f 1 2 ,0 (∞). The classical correlation is given by

and the quantum discord is Q(ρ(∞)) = S(ρ a (∞)) + f 1 2 ,0 (∞) -S(ρ(∞)). (B16)


F Buda

P A Glaude

F Battin-Leclerc
email: frederique.battin-leclerc@ensic.inpl-nancy.fr

R Porter

K J Hughes

J F Griffiths

USE OF DETAILED KINETIC MECHANISMS FOR THE PREDICTION OF AUTOIGNITIONS

Keywords: Alkanes, Gas-phase explosion, Autoignition, Cool flames, Detailed kinetic mechanisms

This paper describes how automatically generated detailed kinetic mechanisms are obtained for the oxidation of alkanes and how these models could lead to a better understanding of autoignition and cool flame risks at elevated conditions. Examples of prediction of the occurrence of different autoignition phenomena, such as cool flames or two-stage ignitions are presented depending on the condition of pressure, temperature and mixture composition. Three compounds are treated, a light alkane, propane, and two heavier ones, n-heptane and n-decane.

Introduction

Large-scale hydrocarbon (partial) oxidation processes form the basis of much of the (petro-) chemical process industry and in particular of the production of feedstock for polymer production (such as ethylene and propylene oxide, maleic and terephthalic acids). To prevent explosions, the European industry is in need of basic safety information, safety tools and predictive safety software, for instance, for design of equipment under all anticipated operating conditions. It is common knowledge that gas and vapour explosions are specific risks of the chemical industry. Much time and money are therefore spent in preventing such explosions, for example by eliminating all potential ignition sources. Commonly one thinks of the usual ignition sources such as flames, hot surfaces and sparks, but forgets that sufficiently high temperatures can lead to autoignition of the fuel-air-mixture.

An important feature of the oxidation of hydrocarbon is the possible occurrence of cool flame phenomenon at temperatures several hundred degrees below the minimum autoignition temperature. During a cool flame, or multiple cool flames, the temperature and the pressure increase strongly over a limited temperature range (typically up to 200°C), but the reaction stops before combustion is complete. Cool flames in themselves do not present a significant hazard, but they change the composition of the initial mixture, making it more reactive, with the accumulation of peroxide intermediates, which play an important part in causing normal ignition, as a two-stage process [START_REF] Pekalski | The relation of cool flames and auto-ignition phenomena to process safety at elevated pressure and temperature[END_REF].

A better understanding of autoignition and cool flame risks at elevated conditions and of the fundamentals of hydrocarbon oxidation leading to new and better tools for process engineering is therefore a very important challenge. The chemistry of hydrocarbon oxidation is very complex, as a mixture of hydrocarbons and air can react either in a slow combustion process or through an uncontrolled exponential increase in rate which leads to ignition, and can only be well described by using detailed kinetic mechanisms.

Detailed kinetic mechanisms are based on elementary steps, the rate constants of which depend only of temperature and pressure and can, therefore, be used in a predictive way.

Detailed reaction mechanisms of the oxidation of hydrocarbons in the gas phase have been developed world-wide over several decades. Most studies have been confined to oxidation chemistry at temperatures above about 1000 K with more limited attention of detailed kinetic models to the important lower temperature regime. A major problem in constructing a detailed chemical kinetic model, especially at lower temperatures, is the very large number of possible reactions, products, and reaction intermediates involved.

Because manual assembly of a comprehensive kinetic model is extremely difficult and prone to error, the only practical way to construct and use large models lies in the use of formal computer-based methods. The work presented here is based on the use of EXGAS, a fully automatic software able to generate detailed reaction mechanisms for gas-phase oxidation reactions [START_REF] Warth | Computed aided design of gas-phase oxidation mechanisms: Application to the modelling of normal-butane oxidation[END_REF][START_REF] Warth | Computer based generation of reactions mechanisms for gas-phase oxidation[END_REF].

The purpose of this paper is to show how automatically generated detailed kinetic mechanisms can model the oxidation of alkanes and how they can predict the occurence of different autoignition phenomena, such as cool flames or two-stage ignitions, depending on the condition of pressure, temperature and mixture composition. We will first describe the case of a light alkane, propane, and then, that of two heavier compounds, n-heptane and n-decane.

Automatically generated detailed mechanism of the oxidation of alkanes

Although the global reaction for the complete combustion of an alkane and its exothermicity, e.g. propane,

C 3 H 8 + 5 O 2 = 3 CO 2 + 4 H 2 O , ∆ r H(1000K) = -2045.63 kJ/mole
can be used to estimate the maximum explosion pressure, it is of no help to determine the temperature and the pressure where autoignition can take place, to explain the happening of cool flames (see further in the text) or to assess the by-products, such as carbon monoxide (CO) or aldehydes, which are also formed. For that is necessary to investigate in more details the reaction mechanism by considering the all possibilities of breaking and creation of chemical bonds, that is to say, by writing all the possible elementary reactions.

For the oxidation of alkanes containing at least three atoms of carbon, the possible elementary reactions of the reactant molecules and the derived radicals can be grouped according to the 7 following main types (taking as examples reactions related to the oxidation of propane) :

1) Unimolecular and bimolecular initiation steps, e.g.

C 3 H 8 ↔ •CH 3 + •C 2 H 5 ; C 3 H 8 + O 2 ↔ •C 3 H 7 + •OOH,
2) Reactions leading to alkenes from alkyl and hydroperoxyalkyl radicals :

•

Decompositions by beta-scission, e.g.

OO• + •OOH ↔ C 3 H 7 OOH + O 2 ,
As the rate constants of these elementary steps depend only on temperature and pressure, a model built in that way can be used to predict the autoignition conditions or the formation of by products. But to do so, this model has to be comprehensive, so that no important channel has been omitted. For that reason, a software EXGAS, has been developed in order to write systematically all the possible elementary reaction, by using an algorithm to ensure the comprehensiveness of the obtained mechanism [START_REF] Warth | Computer based generation of reactions mechanisms for gas-phase oxidation[END_REF]. Details about the way the reactions of compounds containing less than two atoms of carbon and the decomposition of stable products obtained by the above-mentioned elementary steps are also taken into account can be found in previous papers [START_REF] Warth | Computed aided design of gas-phase oxidation mechanisms: Application to the modelling of normal-butane oxidation[END_REF].

For a given alkane, or mixture of alkanes, EXGAS produces a file compatible with the CHEMKIN II softwares library [START_REF] Kee | Chemkin II. A fortran chemical kinetics package for the analysis of gas-phase chemical kinetics[END_REF]) and containing all the relevant elementary reactions, as well as all the thermodynamic and kinetic data related [START_REF] Warth | Computed aided design of gas-phase oxidation mechanisms: Application to the modelling of normal-butane oxidation[END_REF]. While the number of reactions involved in the mechanism can rapidly be important, the number of parameters needed to estimate the rate constants remains much more limited, as they are estimated using correlations mainly based on the type of elementary step and on the type of bond which is broken. For instance, the rate constant of the abstraction of an H-atom from the reactant (6) depends only on the type of H-atom which is abstracted, as shown in table 1. The same correlation is used for each linear alkane.

FIGURE 1

The two next parts show how models generated by EXGAS for the oxidation of propane, n-heptane and n-decane between 600 and 1500 K can reproduce observed experimental results.

Prediction of the phenomena observed during the oxidation of propane in a closed vessel

The processes described above may be exemplified from the modelling of ignition in small vessels at sub-atmospheric pressures. A definitive pressure-temperature (p-T a )

ignition diagram for the combustion phenomena in propane was published by [START_REF] Newitt | The oxidation of propane[END_REF] and there have been many subsequent studies [START_REF] Griffiths | Experimental and numerical studies of oxidation chemistry and spontaneous ignition[END_REF] but it is only very recently that the full richness of the behaviour has been modelled successfully [START_REF] Griffiths | The role and rate of hydrogen peroxide decomposition during hydrocarbon two-stage autoignition[END_REF].

The simulation of multiple cool flames, and ignitions in a propane + oxygen mixture at a total pressure of 300 Torr (0.4 bar) is shown in Figure 2. The transition from cool flames through a three-stage ignition to two-stage ignition is brought about by a reduction in temperature from 610 K to 570 K and, as such, is a dramatic example of the consequence of the negative temperature dependence of reaction rate on temperature, as discussed below. The complete ignition diagram obtained by simulation is shown in Figure 3.

FIGURES 2 AND 3

Modelling of the autoignition delay times of n-heptane and n-decane

The autoignition of n-heptane have been recently investigated in two complementary apparatuses : [START_REF] Minetti | A rapid compression machine investigation of oxidation and autoignition of n-heptane : measurements and modeling[END_REF] has measured auto-ignition delays in a rapid compression for mixtures n-heptane/oxygen/argon/nitrogen for temperatures after compression from 600 to 900 K and pressures from 3.0 to 4.6 bar. [START_REF] Ciezki | Shock-tube investigation of self-ignition of n-heptane-air mixtures under engine relevant conditions[END_REF] have studied the autoignition behind reflected shock wave of mixtures n-heptane/oxygen/argon for temperatures from 660 and to 1350 K and pressures from 3 to 42 bar.

Figure 4 shows the simulated temporal profiles of pressure and concentration of hydroxyl radicals vs. which can be obtained after the compression. This simulated profile of pressure exhibits the same behaviour as that of the experiments. After an induction period around 8 ms, a small rise of pressure is observed and is followed 25 ms later by a sharp increase corresponding to autoignition. That is characteristic of the two-stage autoignition, which starts first by a cool flame putted in evidence by the small pressure rise around 8 ms and is followed later by the actual autoignition. The time elapsed between the end of the compression and the sharp rise of pressure is the ignition delay time, which is a good characteristic of the reactivity of a compound.

As hydroxyl radicals are the main chain carriers, their concentration is a good indicator of the global rate of the reaction. Figure 4 presents the simulated profile of the mole fraction of these radicals and shows how a first acceleration of the reaction around 8 ms is followed by a slow down, while the second acceleration leads to autoignition.

Figure 4 shows also how during the cool flame an important concentration of hydrogen peroxide, which is a very explosive species used in rocket engines, is built on and is consumed only during the ignition itself. Figures 5 and6 shows that our model reproduce well the fact that the higher the pressure, the lower the ignition delay time. Modelling shows that an increase of the pressure from 50 to 100 bar induce only a 30% decrease of the ignition delay time.

FIGURES 5 AND 6

These figures illustrate also clearly the existence of a zone of negative temperature coefficient regime between 700 and 900 K, which is characteristic of the oxidation of alkanes and is displaced towards the higher temperatures when the pressure increases, and show that our detailed kinetic models are able to correctly capture these phenomenena.

Figure 7 displays a sensitivity analysis performed for n-heptane between 650 and 700 K, i.e. in the lowest part of the studied temperatures range. This figure shows that the primary reactions with the highest promoting effect are additions to oxygen molecule, which is the most sensitive reaction, metatheses and isomerizations. A secondary reaction, the decomposition of hydroperoxide species has also an important promoting effect. The primary reactions with the highest inhibiting effect are those directly competing with the additions to oxygen, i.e. the oxidations (to give alkenes and HO 2 • radicals), the beta-scissions of alkyl radicals and the formations of cyclic ethers from hydroperoxyalkyl radicals.

FIGURE 7

The formation of hydroperoxides is extremely important, because they include a O-OH bond, which can easily be broken and lead to the formation of two radicals, which can in their turn react with alkane molecules to give alkyl radicals. This multiplication of the numbers of radicals in a chain reaction induces an exponential acceleration of reaction rates, leading in some conditions to spontaneous auto-ignition.

The reversibility of the addition of alkyl radicals to oxygen molecules when the temperature increases to the benefit of the formation of alkenes leads to an overall reduction of the reaction rate and induces the negative temperature dependent regime.

When the temperature increase further, as a result of relatively weak exothermic propagation, other reactions (such as

H 2 O 2 → 2OH• and H• + O 2 → •OH + •O•) ensure
the multiplication of the number of radicals and the reactivity accelerates again.

The decrease of reactivity with temperature explains the limited increase of temperature observed during a cool flame, since even a small elevation of temperature (a few ten of degrees) can lead to an important reduction of the reactivity, which, in some conditions almost stops the course of the reaction. In some cases, several accelerations/slow downs of the reaction can be observed, which are then observed as multiple cool flames.

Conclusion

Figures captions:

Figure 1 : Rules governing the generation of the primary mechanism for the oxidation of alkanes. Numbers refer to the types of elementary steps described in the text.

Peroxide species are written in italics. .

Figure 1

 1 Figure1presents the global structure of a mechanism generated by EXGAS for an

FIGURE 4 Figure 5

 45 FIGURE 4

  Detailed kinetic mechanisms, automatically generated based on generic types of elementary reactions, have been successfully used to model autoignition and cool flame characteristics for propane, n-heptane and n-decane. The experimental and modelled pressure-temperature ignition diagrams mapping the transition between the different observed oxidation phenomena (slow reaction, cool flames, multiple stages ignition) are qualitatively in good agreement for propane. The variations of autoignition delay times with temperature are well reproduced for n-heptane and n-decane, including the existence of a "negative temperature coefficient" zone, where the global reactivity of the mixture decreases when temperature increases. A sensitivity analysis has allowed the reactions of influence to explain these behaviours to be determined.

Figure 2 :

 2 Figure 2: Simulated temperature -time profiles for multiple cool flames and two-stage

Figure 3 :

 3 Figure 3: Simulated pressure -temperature ignition diagram for propane combustion

Figure 4 :

 4 Figure 4: Simulated profiles of pressure and concentration of hydroxyl radicals after

Figure 5 :

 5 Figure 5: Modelling of ignition delay times obtained in a rapid compression machine

Figure 6 :

 6 Figure 6: Modelling of ignition delay times versus temperature for n-decane/air

Figure 7 :

 7 Figure 7: Sensitivity analysis computed in the case of the autoignition of n-heptane at a

  Figure 7

  3 H 7 + O 2 ↔ C 3 H 6 + •OOH, 3) Additions of alkyl and hydroperoxyalkyl radicals on a molecule of oxygen, e.g. •C 3 H 7 + O 2 ↔ C 3 H 7 OO•, •C 3 H 6 OOH + O 2 ↔ •OOC 3 H 6 OOH, 4) Isomerizations of alkyl and peroxy radicals involving cyclic transition state, e.g. CH 3 CH 2 CH 2 OO• ↔ •CH 2 CH 2 CH 2 OOH, 5) Decompositions of hydroperoxyalkyl and di-hydroperoxyal alkylradicals to form cyclic ethers, alkenes, aldehydes or ketones oxohydroperoxyalkanes, e.g. •CH 2 CH 2 CH 2 OOH ↔ cyclic-C 3 H 6 O + •OH,

•C 3 H 7 ↔ •CH 3 + C 2 H 4 , • Oxidations to form the conjugated alkene and •OOH, e.g. •C e.g. CH 2 OOHC•OOHCH 3 ↔ CH 2 OOHC=OCH 3 + •OH, 6) Metathesis reactions to abstract an H atom from the initial reactant, e.g. •OH + C 3 H 8 ↔ H 2 O + •C 3 H 7 , 7) Termination steps : • Combinations of two free radicals, e.g. CH 3 • + CH 3 • ↔ C 2 H 6 , • Disproportionations of peroxy radicals with •OOH, e.g. C 3 H 7

Table 1 :

 1 Examples of quantitative structure-reactivity relationships used for the estimation of rate constants in the case of the abstraction (reaction 6) of an H-atom from linear alkane (per H-atom abstractable). Rate constants are assumed to follow the modified Arrhenius law: k = A T b exp(-E/RT), with the units cm 3 , mol, s, kJ.

		Primary H-atom	Secondary H-atom
		example for propane :	example for propane :
						H H
		H CH 2 -CH 2 -CH 3 H CH 2 -CH 2 -CH 3	CH 3 -CH-CH 3 CH 3 -CH-CH 3
	Radicals	A	n	Ea	A	n	Ea
	•O•	13.23	0	32.8 13.11	0	21.7
	•H	6.98	2	32.2	6.65	2	20.9
	•OH	5.95	2	1.88	6.11	2	-3.20
	•OOH	11.30	0	71.1 11.30	0	64.8

Acknowledgement

This work has been supported by the European Commission within the "Safekinex" Project EVG1-CT-2002-00072.


T J Li

L Bellon

Dissipation of micro-cantilevers as a function of air pressure and metallic coating

Keywords: 46.35.+z -Viscoelasticity, plasticity, viscoplasticity 05.40.Jc -Brownian motion 07.10.Cm -Micromechanical devices and systems

In this letter, we characterize the internal dissipation of coated micro-cantilevers through their mechanical thermal noise. Using a home-made interferometric setup, we achieve a resolution down to 10 -14 m/ √ Hz in the measurement of their deflection. With the use of the fluctuation dissipation theorem and of the Kramers-Kronig relations, we rebuilt the full mechanical response function from the measured noise spectrum, and investigate frequency dependent dissipation as a function of the air pressure and of the nature of the metallic coatings. Using different thicknesses of gold coatings, we discuss the source of the internal viscoelastic damping.

Microsized cantilevers are present in many applications, ranging from chemical and biological sensors [1] to atomic force microscopy (AFM) [START_REF] Meyer | Scanning probe microscopy: The lab on a tip[END_REF]. They are also used as a basic brick of microelectromechanical systems (MEMS). As mass detectors for example, they can be functionalized for the adsorption of specific chemical compounds by proper coating of their surface: the mass increment is detected as a resonance frequency shift. Gold coating is widely used due to the large selection of materials that can be adsorbed via thiol chemistry [1,[START_REF] Ziegler | [END_REF]. Besides, from a gold layer, electrically conducting layers can be patterned and integrated into the cantilever for local heating [4], magnetomotive actuation [START_REF] Lange | Magnetic actuation and mos-transistor sensing for cmos-integrated resonators in 15th[END_REF] or piezoresistive readout of the deflection [START_REF] Rasmussen | [END_REF].

The sensitivity of the resonant cantilever when used as a mass sensor depends on the spectral resolution, thus of its quality factor Q defined as the ratio of stored vibrational energy over energy lost per cycle of vibration [7]. The greatest sensitivity will thus be reached with the smallest dissipation. The functionality of MEMS, AFM probe or mass sensors is based on the deformation of the cantilever. Thermally induced mechanical fluctuations determine the ultimate deflection sensitivity of these sensors and represent one of the most important noise sources. They are linked to the damping of the system, as shown by the fluctuation dissipation theorem (FDT) [8]. Smaller damping (a) Email: Ludovic.Bellon@ens-lyon.fr will thus lead to reduced thermal noise, thus increasing the sensitivity and usability of the probes. It is therefore of prime importance to understand and characterize the dissipation sources in these systems.

Great effort has been put into describing the effect of ambient pressure and cantilever coating on their resonant behavior, both experimentally [9][10][11][12][13][14] and theoretically [15][16][17]. For example, Sandberg [13] investigated the effect of gold coating on the quality factor of a resonant cantilever and showed that in vacuum Q is severely reduced by the deposition of even a thin gold film (100nm), especially for higher order modes. Considering only structural damping, Saulson [15] proposed a viscoelastic model, in which the power spectrum density (PSD) of thermal induced deflection presents a characteristic 1/f like trend. In a previous work [14], we introduced a simple power law to describe the frequency dependence of this viscoelasticity on a gold coated cantilever, and a model that includes Sader's approach to describe the coupling with the surrounding atmosphere [17,18]. We showed that the damping is only due to the coating when viscous dissipation vanishes in vacuum. Understanding the source of this coating induced viscoelasticity will undoubtedly help in designing more sensitive and accurate cantilever based sensors and microprobes. However, measuring the thermal noise or small damping over a wide range of frequency is a great challenge, and very few experiments [14,[18][19][20][21] have succeeded so far in directly measuring fluctuations out of res-onance, notably at low frequency.

In this work, we measure the thermal noise of AFM cantilevers, with a highly sensitive interferometric technique [14,18], allowing access to the whole spectrum from 1 Hz to 20 kHz with a background noise as low as 10 -28 m 2 /Hz (see fig. 1). With the use of the FDT and of the Kramers-Kronig relations, we rebuild the full mechanical response function from the measured PSD [14]. We can therefore characterize dissipation of micro-cantilevers as a function of the air pressure and of the metallic coatings. Let us first recall briefly the theoretical background for our analysis of thermal noise, before describing the experiments and concluding with a discussion on the results.

Since our cantilevers are in equilibrium at temperature T , the thermal fluctuations of their deflection d are described by the FDT, relating the PSD to the mechanical response function G of the system:

S d (f ) = - 4k B T ω Im 1 G(ω) (1)
where k B is the Boltzmann constant, ω = 2πf the pulsation corresponding to frequency f , and Im stands for the imaginary part of its argument. G(ω) is defined in the Fourier space as G(ω) = F (ω)/d(ω), where F is the force coupled to d in the Hamiltonian of the system. Considering both viscoelastic and viscous damping [14,15], we have:

G(ω) = k 1 - ω 2 ω 2 0 + i ω ω 0 Q a + φ (2)
with k * = k(1 + iφ) the complex spring constant, ω 0 the resonant pulsation and Q a the quality factor linked to viscous damping in air. At resonance, we observe an effective quality factor Q eff defined as:

1 Q eff = 1 Q a + φ (3)
In the experiments, it is hard to characterize dissipation outside resonance: the imaginary part of G much smaller than its real part, and the forcing method has to be perfectly controlled to perform such measurement. We therefore use the FDT (eq. 1) to infer Im[1/G(ω)] from the measurement of the thermal noise S d , then Kramers-Kronig integral relations to rebuilt the full mechanical response function G from the knowledge of this imaginary part. Using this method, we have shown that both Q a and φ are in fact weakly frequency dependent for the model to match the experimental observations [14]. Q a (ω) is described by Sader's approach for viscous damping [17], and the viscoelasticity is following a simple power law with a small exponent: φ(ω) ∝ ω α , with α ≈ -0.11 for a commercial gold coated cantilever [14]. In this letter, we demonstrate how viscous damping vanishes when the pressure decreases, and characterize the behavior of φ for different coatings.

We use silicon AFM cantilevers (Budget Sensors CONTx) with the following geometry: length 450±10 µm, width The inset in (a) illustrates our high precision interferometric method: interferences between beams reflecting on the base and free extremity of the cantilever directly measure its deflection d [14,22]. Inset in (b): effective quality factor Q eff at resonance as a function of pressure.

50 ± 5 µm, thickness 2 ± 1 µm. Their resonant frequency is around 13 kHz. We present here data corresponding to three samples : the coating of cantilever cAu is 5 nm of titanium plus 70 nm of gold (both sides), that of cantilever cAl is 70 nm of aluminum, and that of cantilever cPt is 70 nm of PtIr5. In addition, we performed successive gold layers deposition on three other probes: cantilever dAu is initially equivalent to cAu, cantilever dPt is initially equivalent to cPt, and cantilever dSi is initially a raw silicon cantilever, prepared with a 2 nm Cr adhesion layer. Layer thickness increments were around 10 nm for the first 3 evaporations, then around 20 nm.

We first study the influence of ambient pressure P on cantilever cAu. In fig. 1, we plot the PSD of the deflection for different pressures from ambient to 10 -3 mbar, and the corresponding reconstructed response function (imaginary part only, normalized by spring constant k). Dissipation is clearly the sum of two contribution: pressure independent viscoelasticity φ(ω), always dominant at low frequency, and viscous damping term ω/ω 0 Q a (ω), dominant at high frequency and high pressure, but vanishing in vacuum. At resonance, the effective quality factor thus increases when the pressure drops, but saturate to a finite value due to the internal dissipation. For this cantilever and the frequency range probed here, viscous damping becomes negligible when P < 10 -2 mbar. In vacuum, the dissipative part of the response function directly leads to the viscoelaticity: φ(ω) = Im[G(ω)]/k when Q a → ∞. A simple power law matches well the observed frequency dependence [14], with an exponent α Au ≈ -0.15.

To check the generic validity of this behavior, we measure in vacuum the viscoelastic damping of different metallic coatings, and report the results for cantilevers cAu, cPt and cAl on fig. 2. The magnitude of this internal dissipation is smaller for the 2 other coatings, but the weak frequency dependence is a common characteristic. At high frequency (above 1 kHz), the viscous damping is still observable for cantilever cPt at the lowest pressure achievable in our system (10 -3 mbar). If we limit the frequency range to 1 kHz for this sample, we observe a similar power law dependance with an exponent α P t ≈ -0.12. The data of the last cantilever, cAl, is a a little bit different with a very flat measurement, leading to α Al ≈ 0. The power law dependence thus works generically, though the exponent changes slightly for each cantilever.

To investigate the origin of this internal damping, we finally measured the viscoelasticity of cantilevers dSi, dPt and dAu as a function of the thickness of an added gold coating : the same 3 cantilevers where characterized in vacuum between successive layers deposition, and we plot in fig. 3 resonance (value of 1/Q eff , deduced from a lorenzian fit of the PSD). Viscoelasticity remains weakly frequency dependent in these observations, and is roughly proportional to the thickness. This behavior suggests that the main contribution to the internal damping of the cantilever originates in the bulk of the coating, rather than from a surface or interface effect. Indeed, if such phenomenon was dominant, it would be present even for the lowest thickness, whereas we observe a vanishing dissipation when the coating layer gets thinner. Moreover, the value of the dissipation is independent of the adhesion layer for the gold coating : we get the same results for 2 nm of Cr (cantilever dSi) and for 70 nm of PtIr5 (cantilever dPt).

As a conclusion, let us summarize the main points illustrated in this letter. We measure the thermal noise of AFM cantilevers with a high resolution interferometer, and reconstruct from this data their dissipation. The damping of coated cantilevers is shown to arise from two distinct sources: viscous damping due to the external atmosphere, dominant at high frequency and ambient pressure, and internal viscoelastic damping, dominant at low frequency or in vacuum. This viscoelasticity is linked to the presence of a coating on the cantilever, and can be modeled by a complex spring constant. This dissipation is accurately fit by a simple power law ω α on up to 4 decades in frequency, with a small exponent α depending on coating material and procedure: for the cantilevers probed here, α Au ≈ -0.15, α P t ≈ -0.12, and α Al ≈ 0. Damping increases by one order of magnitude between PtIr5 and Au coating, Al falling in between those 2 materials. Finally, by changing the coating thickness, we show that dissipation is due to the bulk of the coating rather than to some interface friction between layers.

These results demonstrate that the choice of the coating material is critical with respect to internal dissipation in micro-cantilevers, and that gold is actually a bad choice with respect to this criterium. If this material is however required for its chemical properties, using the lowest functional thickness is advisable to minimize damping. Our characterization procedure features an excellent resolution with measurement of mechanical loss tangents down to 10 -4 . This way, the viscoelasticity due to the coating can be accurately quantified and our measurements should be useful in the perspective of testing models of internal friction, eventually leading to improved coating procedures and better performance of cantilever based sensors. Our method would also be suited to study other type of coatings, such as those implied in chemical or biological sensors, alone or linked to the target molecules. * * *

Fig. 1 :

 1 Fig. 1: (a) Thermal noise spectrum and (b) reconstructed dissipative part of the response function G of cantilever cAu at different pressures.Background noise due to the electronics has been subtracted from the measured PSD. Viscous damping, dominant at high frequency and high pressure over viscoelastic dissipation, vanishes in vacuum, leading to a sharper resonance. The inset in (a) illustrates our high precision interferometric method: interferences between beams reflecting on the base and free extremity of the cantilever directly measure its deflection d[14,22]. Inset in (b): effective quality factor Q eff at resonance as a function of pressure.

Fig. 2 :

 2 Fig. 2: Viscoelasticity φ(ω) of cantilevers cAu, cAl and cPt reconstructed from the noise spectrums measured in vacuum.

Fig. 3 :

 3 Fig. 3: Viscoelasticity at low and high frequency (average of φ(ω) around 150 Hz and 1/Q eff at resonance) as a function of the total gold layer thickness for cantilevers dSi, dAu and dPt. Typical error bars (dispersion over a few acquisitions) are shown for two cantilevers.

We thank F. Vittoz and F. Ropars for technical support, and S. Ciliberto, A. Petrosyan and C. Fretigny for stimulating discussions.