Bath dynamics in an exactly solvable qubit model with initial qubit-environment correlations

We study the bath dynamics in the dephasing model of a two-state quantum system (qubit) coupled to an environment of harmonic oscillators. This model was shown [Morozov et al., Phys. Rev. A, 2012, 85, 022101] to admit the analytic solution for the qubit and environment dynamics. Using this solution, we derive the exact expression for the bath reduced density matrix in the presence of initial qubit-environment correlations. We obtain the non-equilibrium phonon distribution function and discuss in detail the time behavior of the bath energy. It is shown that only the inclusion of dynamic correlations between the qubit and the bath ensures the proper time behavior of the quantity which may be interpreted as the"environment energy".


Introduction
The dynamics of open quantum systems has attracted a great deal of interest over the last few decades. Its importance is dictated by the prospects of applications in quantum optics, quantum computing, quantum measurements and control [1,2], as well as by the necessity of a deeper understanding of the theory itself [3][4][5][6][7][8]. The dynamics of open quantum systems was studied in several aspects: (i) the effect of initial correlations between an open system and its environment has been investigated in [9,10]; (ii) a new viewpoint concerning the nature and the measure of non-Markovianity has been presented in [11,12]; (iii) the effect of non-equilibrium environment on quantum coherence and the level populations has been considered in [13,14].
The latter problem is of a particular interest. Usually, when constructing a master equation for the reduced density matrix of an open system, one considers the bath to be at thermal equilibrium [1], even though there is a build-up of dynamical correlations [4,9] caused by the entanglement of quantum states. Moreover, even if the effect of non-equilibrium environment on the system behavior is not neglected [13], the intrinsic bath dynamics is beyond consideration, since the corresponding bath variables are always integrated out from the equations of motion. Although this approach seems to be quite natural as long as one studies solely the open system dynamics, the investigation of the bath evolution itself can undoubtedly be an interesting problem, yielding some useful hints about how to deal with more realistic systems, especially with those which possess slow relaxation to equilibrium and do not admit exact solutions.
In this brief report, we show the relevance of proper (or intrinsic) bath dynamics. The paper is structured as follows. In section 2 we derive an exact expression for the bath reduced density matrix in the so-called dephasing model, describing a two-state system (qubit) coupled to a bosonic bath [5][6][7]. In section 3 the non-equilibrium distribution function for the bath modes (phonons) is calculated. Special at-tention is paid to the description of time evolution of the phonon energy and the correlation energy in the "qubit-bath" system. In the last section we discuss the results and draw final conclusions.

Bath density matrix in the dephasing model
We consider a simple version of the spin-boson model describing a two-state system (qubit) (S) coupled to the bath (B) of harmonic oscillators [5][6][7][8]. In the "spin" representation for the qubit, the total Hamiltonian of the model is written as follows (in our units ħ = 1) where ω 0 is the energy difference between the excited |1〉 and the ground |0〉 states of the qubit, and σ 3 is one of the Pauli matrices σ 1 , σ 2 , σ 3 . Note that the operator σ 3 can also be written in the basis of the ground |0〉 and the excited |1〉 states of the two-level system as σ 3 = |1〉〈1| − |0〉〈0|. Bosonic creation and annihilation operators b † k and b k correspond to the kth bath mode with frequency ω k , and g k are the coupling constants. The distinctive feature of the dephasing model (2.1) is that the average populations of the qubit states do not depend on time, and hence there is no relaxation to the complete equilibrium between the qubit and the environment. In other words, the model is non-ergodic 1 . Note, however, that the Heisenberg picture operators σ ± (t ) = [σ 1 (t )±iσ 2 (t )]/2 evolve in time, leading to a non-trivial decay of the coherences 〈σ ± (t )〉. Thus, we have a unique situation, where the system relaxation may be interpreted physically as "pure" decoherence and the entropy exchange [8] rather than the energy dissipation.
It was shown in [7] that equations of motion for the Heisenberg picture operators σ ± (t ), b † k (t ), and b k (t ) can be solved exactly with the results We shall use the exact expressions (2.2)-(2.3) to evaluate the reduced density matrix of the bath. We start with the obvious relation for the non-equilibrium density matrix of the composite system (S +B) and assume that the initial density matrix has the form [7] where the projector P ψ = |ψ〉〈ψ| can be expressed in terms of the Bloch vector v = 〈 σ〉 [7,8]: where 〈σ 3 〉 = 〈ψ|σ 3 |ψ〉, 〈σ ± 〉 = 〈ψ|σ ± |ψ〉, and I denotes the 2 × 2 identity matrix. In equation (2.6) |ψ〉 = a 0 |0〉 + a 1 |1〉 with |a 0 | 2 + |a 1 | 2 = 1 is the state vector of the qubit. The constrained initial density matrix of the bath, ρ B (ψ), is given by 2 Physically, the density matrix (2.6) corresponds to a situation where at times t < 0 the open system S is in thermal equilibrium with its environment B, and at time zero one makes a perfect (selective) measurement on the system S only. As a result [1], the system S is prepared in some pure state |ψ〉. Making use of equations (2.6)-(2.7) and noting that σ 3 is an integral of motion, we can recast equation (2.5) into the form where we have introduced the notation (2.11) The bath density matrix is obtained from (2.10) by taking the trace over the qubit states, Since H (±) B (−t ) does not contain the spin operators σ ± , the terms with σ ± (−t ) do not contribute to (2.12). Expressing the probabilities |a i | 2 in terms of 〈σ 3 〉 in a usual way, it is straightforward to manipulate the time-dependent bath density matrix (2.12) to denotes a phonon part of H int in the interaction picture (with the replacement t → −t , see equation (2.5)), whereas the quantity is the non-equilibrium correction to the phonons energy (see also the next section for discussion). The expression for the bath density matrix becomes much simpler if one neglects the initial correlations in the system by taking a direct product instead of (2.6). Proceeding in a similar way, after some algebra one obtains where ∆H I (t ) = H I (t ) − H I (0) denotes the non-equilibrium contribution to the correlation energy in the interaction picture. It is seen from equation (2.17) that even in this case there is a dynamical build-up of correlations in the system. We touch upon this point in the next section, when analyzing the relevance of non-equilibrium correlations.

Phonon non-equilibrium distribution function and energy
To gain some insight into the time behavior of the bath modes (phonons), let us first calculate the phonon distribution function n k (t ). This can be done in two equivalent ways: either using the exact expression (2.13) for the bath density matrix, or averaging the bosonic Heisenberg picture operators (2.3) over the initial state of the composite system. Here, we shall follow the latter procedure which is simpler.
Taking the initial density matrix in the form (2.6) and then applying the unitary transformation technique [7] (or the method of the displaced harmonic oscillator modes [15]), after some straightforward algebra one obtains is the initial phonon distribution function, and the function represents the contribution of initial correlations. With equations (3.1) and (3.2), it is easy to calculate the time evolution of the non-equilibrium phonon energy ε(t ) = k ω k n k (t ). The final result is conveniently written in terms of the bath spectral density J (ω) which is introduced by the well-known rule [1,2,7] After simple manipulations we arrive at (3.5) where the initial phonon energy is given by Here the last term occurs due to initial correlations in the system. Usually [1,2,7], the spectral density function is chosen in the form where s > 0 and λ s is a dimensionless coupling constant. This formula ensures both a proper lowfrequency behavior of J (ω) and a cut-off at high frequencies (ω ≫ Ω). The case s = 1 is usually called the "Ohmic" case, the case s > 1 "super-Ohmic", and the case 0 < s < 1 "sub-Ohmic". Using expression (3.7), it is possible to analyze the time behavior of ε(t ) for different s, but in this brief report we would like to discuss only one physically interesting point related to the result (3.5) for the phonon energy. At first glance, the fact that the phonon energy (3.5) depends on time may appear as an apparent paradox. Indeed, on the one hand, one may conclude that there is an energy exchange between the qubit and the bath. On the other hand, the qubit Hamiltonian H S commutes with the total Hamiltonian (2.1), and hence the qubit energy 〈H S 〉 does not depend on time. To explain this paradox, let us calculate the nonequilibrium correlation energy ε cor (t ) = 〈H int (t )〉, where H int (t ) is the interaction term in equation (2.1) (in the Heisenberg picture) and the average is taken over the initial state (2.6). After some algebra, which we omit, we obtain ε cor (t ) ≡ Tr S,B ρ(t = 0)H int (t ) = ε 0 + 2A(ψ)〈σ 3 〉 + 1 ∆ε − ε(t ). (3.8) Combining this expression with equation (3.5), it is easy to check that the sum ε(t ) + ε cor (t ) is a timeindependent quantity (see also a footnote on page 2). We see that the dynamics of the correlation energy ε cor (t ) exactly compensates the time dependence of the non-equilibrium phonon energy ε(t ), ensuring the energy conservation law. Physically, the sum ε(t ) + ε cor (t ) is precisely the quantity which should be interpreted as the environment energy.
One more remark is to the point. It can be seen from equations (3.3), (3.5), and (3.8) that both the non-equilibrium phonon energy ε(t ) and the correlation energy ε cor (t ) do not depend on time under conditions 〈σ 3 〉 = ±1. Note in this connection that the correlational contribution γ cor (t ) to the generalized decoherence function [7] vanishes for the same values of the mean inversion population of the levels, manifesting a close relationship between the essentially non-equilibrium behavior of the correlation energy and the onset of the additional channel of decoherence in the system.

Conclusions
Here, we present a summary of the results and discuss their relation to some problems in the dynamics of open quantum systems.
We have derived exact expressions (2.13) and (2.17) for the bath density matrix in the model (2.1) which describes the dephasing mechanism of decoherence in a qubit interacting with a bosonic environment. To the best of our knowledge, a derivation of a bath density matrix has never been performed in the theory of open quantum systems. The explicit form of ρ B (t ) could be essential, for instance, when constructing master equations (especially non-Markovian) and taking into account the intrinsic dynamics of the environment along with the equation of motion for ρ S (t ). Such an approach would modify the well-known Zwanzig-Nakajima projection technique [1,3,9] where the bath degrees of freedom are "eliminated". We believe that this modification is quite natural in the case of the finite size of the bath, when all the environmental modes are involved in the composite system dynamics, and a back-flow of energy (information) from the bath to the open system is essential. Thus, the exact solutions (2.13) and (2.17) can give a valuable insight into general properties of the dynamics of decoherence and can serve as a step toward consistent derivation of master equations ensuring regular behavior of composite systems on all timescales and for strong coupling regimes.
Our analysis of the phonon energy in section 3 illustrates the special role of dynamic correlations between an open system and its environment. We have seen that a "naive" picture with the "energy exchange between the qubit and the bath" is inadequate (even in the case of weak coupling), and only the proper inclusion of non-equilibrium correlations ensures the conservation of the total energy.