Open quantum systems in Heisenberg picture

TL;DR

This paper derives exact formulas for the reduced density matrix, characteristic function, and Wigner distribution of a driven dissipative quantum harmonic oscillator using the Heisenberg picture, with potential generalizations.

Contribution

It provides a novel derivation of the reduced density matrix and related functions for open quantum systems in the Heisenberg picture, including generalizations via Magnus expansion.

Findings

Exact expressions for the reduced density matrix are obtained.

Formulas for the characteristic function and Wigner distribution are derived.

A possible generalization using Magnus expansion is proposed.

Abstract

In the framework of the Heisenberg picture, an alternative derivation of the reduced density matrix of a driven dissipative quantum harmonic oscillator as the prototype of an open quantum system is investigated. The reduced density matrix for different initial states of the combined system is obtained from a general formula, and different limiting cases are studied. Exact expressions for the corresponding characteristic function in quantum thermodynamics and Wigner quasi distribution function are found. A possible generalization based on the Magnus expansion of the evolution operator is presented.

Equations

$$\hat{H}$$

$$\hat{a}(t)$$

$${\hat{a}}^{†}(t)$$

$$\displaystyle G(t)=\mathcal{L}^{-1}\bigg{[}\frac{1}{s+\tilde{\chi}(s)+i\omega_{0}}\bigg{]},$$

$$\tilde{\chi }(s)=\sum _j\frac{|f_j{|}^2}{s+i{\omega }_j},$$

$$M_j(t)=f_j{\int }_0^{t}d{t}^{\prime }{e}^{-i{\omega }_j(t-t^{\prime })}G(t^{\prime }),$$

$$\zeta (t)={\int }_0^{t}d{t}^{\prime }G(t-t^{\prime })K(t^{\prime }).$$

$$|G(t){|}^2+\sum _j|M_j(t){|}^2=1.$$

$${\hat{b}}_j(t)$$

$${\hat{b}}_j^{†}(t)$$

$${\Lambda }_{jk}(t)$$

$${\Gamma }_j(t)$$

$${\Omega }_j(t)$$

$$(\Lambda (t){\Lambda }^{†}(t){)}_{jj}+|{\Gamma }_j(t){|}^2=1.$$

$$\hat{\rho }(t)=\hat{U}(t)\hat{\rho }(0){\hat{U}}^{†}(t).$$

$$⟨n|{\hat{\rho }}_S(t)|m⟩$$

$$\displaystyle\hat{Q}_{mn}(t)={\hat{U}}^{\dagger}(t)\big{(}|m\rangle\langle n|\otimes I_{R}\big{)}\,\hat{U}(t),$$

$$\sum _{n=0}^{\infty }{\hat{Q}}_{nn}(t)=I,$$

$$\sum _{n=0}^{\infty }{n}^s{\hat{Q}}_{nn}(t)=[{\hat{a}}^{†}(t)\hat{a}(t){]}^s.$$

$${\hat{Q}}_{mn}(t)=\frac{1}{\sqrt{m!n!}}\sum _{s=0}^{\infty }\frac{(-1{)}^s}{s!}[{\hat{a}}^{†}(t){]}^{s+m}[\hat{a}(t){]}^{s+n},$$

$$⟨n|{\hat{\rho }}_S(t)|m⟩$$

$$[\hat{a}(t){]}^{s+n}$$

$$[{\hat{a}}^{†}(t){]}^{s+m}$$

$$\hat{B}$$

$${\hat{B}}^{†}$$

$$⟨n|{\hat{\rho }}_S(t)|m⟩$$

$${\hat{\rho }}_R(0)=\frac{1}{Z_R}\prod _j\otimes e^{-\beta \hslash {\omega }_j{\hat{b}}_j^{†}{\hat{b}}_j},$$

$$Z_R=\prod _j{z}_j,$$

$$z_j=\text{mbox}{Tr}_j[e^{-\beta \hslash {\omega }_j{\hat{b}}_j^{†}{\hat{b}}_j}].$$

$$\mbox{Tr}_{R}\Big{(}[{\hat{B}}^{\dagger}]^{u}[\hat{B}]^{v}\,\hat{\rho}_{R}(0)\Big{)}=\delta_{uv}\,u!\,[\eta(t)]^{u},$$

$$\eta (t)=\sum _j\frac{|M_j(t){|}^2}{e^{\beta \hslash {\omega }_j}-1}.$$

$$⟨n|{\hat{\rho }}_S(t)|m⟩$$

$$⟨n|{\hat{\rho }}_S(t)|n⟩$$

$$⟨n⟩=|\alpha G(t)-i\zeta (t){|}^2.$$

Full text

††journal: European Physical Journal Plus

1 Introduction

Experimental techniques in design and manufacturing nanoscale devices to be applied in quantum technology have considerably improved and reached a high level of accuracy in recent years. These devices due to working in the quantum domain are so sensitive to external sources and noises that a precise theoretical understanding of their function is vital to controlling and correcting unwanted behaviors. These devices or systems belong to a wider class of quantum systems interacting with their environment known as open quantum systems Breuer . The subject of open quantum systems covers a wide range of applications in quantum physics due to the fact that no quantum system can be completely isolated from its environment. An important paradigm of an open quantum system is the quantum Brownian motion A1 ; A2 which has been investigated extensively by different approaches B1 ; B2 ; B3 ; B4 ; B5 and appears in miscellaneous problems in physics and chemistry D1 ; D2 ; D3 ; D4 ; D5 ; D6 ; D7 ; D8 ; D9 ; D10 ; D11 ; D12 . Generally, the Hamiltonian of an open quantum system consists of three parts, Hamiltonian of the main system $H_{S}$, Hamiltonian of the environment $H_{R}$ and the interaction Hamiltonian $H_{SR}$. The main ingredient in the context of open quantum system theory is the reduced density operator $\hat{\rho}_{S}(t)$ describing the main system under consideration. This operator is obtained by tracing out the environmental degrees of freedom of total density matrix $\hat{\rho}(t)$ describing the combined system. Knowing the explicit form of the reduced density matrix, a complete description of the quantum dynamics of the main subsystem is achievable.

Another issue that should be taken into account while investigating the thermodynamical properties of nanoscale devices working at the quantum regime is the significance of quantum fluctuations that may lead to a reconsideration of thermodynamical laws. Therefore, a precise investigation of the relation between thermodynamics and quantum mechanics is unavoidable and has opened a new issue nowadays referred to as quantum thermodynamics C1 ; C2 ; C3 ; C4 ; C5 ; C6 ; C7 ; C8 ; C9 ; C10 .

Our aim in the present letter is to introduce an alternative derivation of the reduced density matrix of a driven dissipative quantum harmonic oscillator. Although, the method applied in the present work is based on the existence of exact expressions for the time-evolution of dynamical observables, as a generalization, a perturbative method based on the Magnus expansion (Sec. (11)) can also be developed leading to approximate expressions for the time-evolution of the main dynamical variables. The main result of this investigation is the general formula Eq. (17) and its particular form Eq. (21) from which exact results can be extracted for different initial states of the combined system or for different limiting cases. Due to the important role played by the characteristic function in quantum thermodynamics, an exact expression for this function is given in Eq. (49) and its limiting cases are also considered. In the following, the exact Wigner quasi distribution function corresponding to the reduced density matrix is found.

2 Basics

Here the complex conjugation of an arbitrary c-valued quantity $z$ is denoted by $\bar{z}$, its complex norm by $|z|$ and Laplace transform of an arbitrary function $\varphi(t)$ is denoted by $\mathcal{L}[\varphi(t)]=\tilde{\varphi}(s)$, $(\varphi(t)=\mathcal{L}^{-1}[\tilde{\varphi}(s)]$). The Hamiltonian that we have considered here is the Hamiltonian of a dissipative quantum harmonic oscillator under the influence of a classical external source given by

[TABLE]

where $f_{j}$’s are coupling constants coupling the oscillator to its environment and $K(t))$ is an arbitrary time-dependent classical external source. Our aim is to find the exact reduced density matrix of the oscillator as the main subsystem in an alternative simple and efficient way. For this purpose, we first find the time-evolution of the oscillator ladder operators in the Heisenberg picture as (App. A)

[TABLE]

where we have defined

[TABLE]

From Eqs. (2) and $[\hat{a}(t),\hat{a}^{\dagger}(t)]=1$ we deduce

[TABLE]

From Heisenberg equations for reservoir operators we find (App. A)

[TABLE]

where

[TABLE]

From $[\hat{b}_{j}(t),\hat{b}^{\dagger}_{j}(t)]=1$ one easily finds

[TABLE]

3 Matrix elements of the reduced density matrix

Let us assume that the initial density matrix ($\hat{\rho}(t=0)$), of the oscillator and its environment is a separable state $\hat{\rho}(0)=\hat{\rho}_{S}(0)\otimes\hat{\rho}_{R}(0)$. Then, the total density matrix at an arbitrary time $t$ is

[TABLE]

where the unitary operator $\hat{U}(t)$ is the evolution operator of the combined system. The reduced density matrix of the oscillator can be obtained by tracing out the environment degrees of freedom of $\hat{\rho}(t)$ as $\hat{\rho}_{s}(t)=\mbox{Tr}_{R}\{\hat{\rho}(t)\}$. For the matrix elements of the reduced density matrix in the basis of number states we have

[TABLE]

where $\mbox{Tr}_{S}$ denotes trace over oscillator degrees of freedom, Tr denotes the total trace and $I_{R}$ is identity operator over the environment Hilbert space. In Eq. (9) we have defined

[TABLE]

with the following properties (App. B)

[TABLE]

Eq. (10) can be rewritten in terms of the Heisenberg operators as (App. C)

[TABLE]

and plays a fundamental role in what follows. Now we can rewrite Eq. (9) as

[TABLE]

on the other hand, from Eq. (2) we have

[TABLE]

and similarly

[TABLE]

where for convenience we defined

[TABLE]

By inserting Eqs. (14) and (15) into Eq. (13), we find

[TABLE]

Eq. (17) is the most general formula representing the reduced density matrix in the basis of number states. In order to proceed, we assume that the initial state of the environment is a thermal Maxwell-Boltzmann state defined by

[TABLE]

where $\beta=1/k_{B}T$, $k_{B}$ is the Boltzmann constant, $T$ is the temperature of the environment and $\mbox{Tr}_{j}$ is the trace operator in the Hilbert space of the $j$th oscillator in the environment. By making use of the generating function method it can be proved that (App. D)

[TABLE]

where we have defined

[TABLE]

Therefore,

[TABLE]

Eq. (21) is the main result of this section. In the last line of this equation, the trace operator acts on subsystem operators at initial time $t=0$ that can be achieved straightforwardly. In the rest of this letter, we will extract different physical results from Eq. (21) by considering different initial states or conditions on the total Hamiltonian given in Eq. (1).

4 The oscillator is initially in a coherent state

In this case we have $\hat{\rho}_{s}(0)=|\alpha\rangle\langle\alpha|$, therefore,

[TABLE]

which is an exact expression suitable for both analytical and numerical calculations. Let us consider some limiting cases. In low temperature limit ($\beta\rightarrow\infty$), using Eq. (20) we have $\eta(t)\rightarrow 0$, therefore, by setting $\eta(t)=0$ in Eq.(22) we obtain

[TABLE]

By setting $n=m$ in Eq. (23) we find the probability of finding the system in number state $|n\rangle$

[TABLE]

which is a Poisson distribution with mean number parameter

[TABLE]

5 The oscillator is initially in a number state

In this case we have $\hat{\rho}_{s}(0)=|k\rangle\langle k|$, since $\hat{\rho}_{s}(t)$ is a hermitian operator we have $\langle n|\hat{\rho}_{S}(t)|m\rangle=\overline{\langle m|\hat{\rho}_{S}(t)|n\rangle}$, so with no lose of generality we can assume $m\geq n$ and find (App. E)

[TABLE]

where $r_{min}=\max(s+n-k,0)$. By setting $n=m$, we can find the probability of finding the oscillator in state $|n\rangle$ at time $t$ knowing that it was initially prepared in the state $|k\rangle$

[TABLE]

In zero temperature limit, Eq. (27) reduces to

[TABLE]

For $k=0$, the initial state is the vacuum state $|0\rangle$ that is a coherent state, in this case $r_{min}=s+n$ and Eq. (28) reduces to

[TABLE]

which is a Poisson distribution that could also be found from Eq.(24) by setting $\alpha=0$. Note that in large-time limit $|G(t)|\rightarrow 0$ and the non-vanishing terms in Eq. (28) are obtained by setting $r=s+n$ leading to the same equation Eq. (29) due to the fact that in large-time limit oscillator will decay to its vacuum state.

6 Expectation values

Let $F(\hat{a}^{\dagger}(0),\hat{a}(0))$ be an arbitrary function in terms of the oscillator ladder operators at the initial time $t=0$. The expectation value $\langle F\rangle_{t}$ at arbitrary time $t$ is (App. F)

[TABLE]

As an example, let $F=\hbar\omega_{0}(\hat{a}^{\dagger}\hat{a}+\hat{a}\hat{a}^{\dagger})/2$ then from Eqs. (2) and the initial state

[TABLE]

we find the energy of the oscillator at time $t$ as

[TABLE]

leading to a probabilistic interpretation of Eq. (4).

7 Special limiting cases

7.1 The dissipation can be ignored and coupling to the external source is strong

In this case the coupling constants $f_{j}$ are zero and from Eqs. (2) we have

[TABLE]

To proceed let the initial state of the oscillator be a coherent state $\hat{\rho}_{S}(0)=|\alpha\rangle\langle\alpha|$, then from Eq. (21) we have (App. G)

[TABLE]

and the diagonal elements are given by

[TABLE]

which is a Poisson distribution function with mean value

[TABLE]

7.2 The external source is switched off

Now let the initial state be the number state $\hat{\rho}_{s}(0)=|k\rangle\langle k|$ with zero coherency. In the absence of a driving force ($\zeta(t)=0$), we set $u=r,\,\,n=m$ in Eq. (21) and find

[TABLE]

Therefore, the evolved reduced density matrix has remained diagonal with zero coherency and diagonal elements are given by

[TABLE]

In low temperature limit, we set $r=0$ in Eq. (38) and find

[TABLE]

which is a binomial distribution with parameter $|G(t)|^{2}$. The Heaviside step function in Eq. (39) is defined by

[TABLE]

and shows that only transitions to lower energy levels ($k\geq n$) are possible.

Now let the initial state be a coherent state $\hat{\rho}_{s}(0)=|\alpha\rangle\langle\alpha|$, and the preferred basis be the number states, then in the absence of external source ($\zeta(t)=0$), we set $p=r=u$ in Eq. (22) and with no loss of generality we can assume $m\geq n$, therefore,

[TABLE]

From Eq. (43) it is seen that the coherency of the density matrix vanishes at large-time limit since the non-diagonal elements tend to zero in this limit (App. H).

In low-temperature limit Eq. (43) reduces to

[TABLE]

and the probability of finding the oscillator in the number state $|n\rangle$ is a Poisson distribution

[TABLE]

8 Quantum thermodynamics. The Characteristic function

Let $P(Q,t)$ be the probability distribution for the heat amount $Q$ to be transferred to the environment between times $t=0$ and $t$. Then Talkner ; Esposito

[TABLE]

where $P[e_{1}]$ is the probability of obtaining $e_{1}$ when measuring environment energy $\hat{H}_{R}$ at $t=0$ and $P[e_{1}\rightarrow e_{2};t]$ is the conditional probability that a measurement of $\hat{H}_{R}$ gives $e_{2}$ at time $t$ when it gave $e_{1}$ at $t=0$. The characteristic function is defined by the Fourier transform

[TABLE]

Let the initial density matrix of the total system be a factorized state $\rho(0)=\rho_{S}(0)\otimes\rho_{R}(0)$, where $\rho_{S}(0)$ is the initial density matrix of the oscillator and $\rho_{R}(0)$ is the initial density matrix of the reservoir which we assume to be a thermalized state

[TABLE]

Now we have

[TABLE]

therefore,

[TABLE]

In Eq. (50), $\mbox{TR}_{S(R)}$ means taking trace over system(reservoir) degrees of freedom and $\hat{H}_{R}(t)$ is the reservoir Hamiltonian in Heisenberg picture defined by

[TABLE]

where $\hat{b}_{j}(t)$ and $\hat{b}^{\dagger}_{j}(t)$ are given by Eq. (2). From the characteristic function $G(\nu,t)$, the moments of $Q$ can be found as

[TABLE]

As an example, let the oscillator be initially prepared in the number state $\rho_{S}(0)=|n\rangle\langle n|$, using Eq. (50), the average heat transferred to the reservoir is

[TABLE]

By inserting Eqs. (2) into Eq. (53) we find

[TABLE]

In zero temperature ($\beta\rightarrow\infty$) and in the absence of external source ($\Omega_{j}(t)=0$), we have

[TABLE]

9 Thermal equilibrium

In the absence of an external source ($\zeta=0$) the oscillator tends to an equilibrium state in large-time limit for arbitrary initial state. This can be easily deduced by setting $r=p=u,\,\,(m\geq n)$ in Eq. (21), we have

[TABLE]

in large-time limit we have $|G(t)|\rightarrow 0$, leading to (App. I)

[TABLE]

which is a thermal state with mean number $\langle n\rangle=\eta(t)$.

10 Wigner function

An important quasi distribution function on phase space is the Wigner function. In this section we find an expression for the Wigner function corresponding to a driven dissipative harmonic oscillator. The components of the reduced density matrix in the continuous position basis $\{|x\rangle\}$ are

[TABLE]

where $\psi_{n}(x)$ is the $n$th eigenvector of the free oscillator Hamiltonian. The Wigner quasi distribution function is defined by

[TABLE]

As a special case, let the initial state of the oscillator be the number state $\hat{\rho}_{S}(0)=|k\rangle\langle k|$, and the external source be switched off, then

[TABLE]

where

[TABLE]

is the Wigner function corresponding to the pure state $|n\rangle\langle n|$ and $P^{\zeta=o,T\neq 0}_{k\rightarrow n}(t)$ is given by Eq. (38). From Eq. (60) it is seen that the quasi distribution function $W$ is the average of quasi distributions $W_{n}$ with respect to the probability distribution $P^{\zeta=o,T\neq 0}_{k\rightarrow n}(t)$. In low temperature limit, using Eq. (39), we find

[TABLE]

As another case, let the initial state of the oscillator be a coherent state $\hat{\rho}_{S}(0)=|\alpha\rangle\langle\alpha|$, then in the zero temperature ($\beta\rightarrow\infty$) we have from Eq. (23)

[TABLE]

which is the Wigner function corresponding to the following evolved coherent state

[TABLE]

11 Generalisation

In the present work we have concentrated on the Hamiltonian Eq. (1) which is the Hamiltonian of a driven dissipative harmonic oscillator. However, there are few systems that their Heisenberg equations are integrable. So a perturbative approach is unavoidable. Among the perturbative methods, the Magnus expansion Magnus has some preferences. The main preference is the preservation of the unitarity of the evolution operator at any order of approximation. Having an approximate evolution operator $\hat{U}_{appr.}(t)$ we find an approximation expression for an arbitrary dynamical variable $\hat{A}$ in Heisenberg picture as

[TABLE]

Let the total Hamiltonian be given by

[TABLE]

where $[\hat{H}_{0},\hat{H}_{1}]\neq 0$. The exact evolution operator $\hat{U}(t)$ satisfies the Schrödinger equation

[TABLE]

To find an approximate unitary solution, following Magnus we assume

[TABLE]

then

[TABLE]

therefore,

[TABLE]

and the time-evolution of an arbitrary dynamical variable $\hat{A}$ in Heisenberg picture is approximately given by

[TABLE]

and similar steps can be followed to find an approximate reduced density matrix in a preferred basis.

12 Conclusions

An alternative derivation of the reduced density matrix of a driven dissipative quantum harmonic oscillator as the prototype of an open quantum system was introduced in the Heisenberg picture. The reduced density matrix for different initial states of the combined system was obtained from a general formula, and different limiting cases were studied. Exact expressions for the corresponding characteristic function in quantum thermodynamics and Wigner quasi distribution function were found. Though the method introduced here was based on the existence of exact expressions for the time-evolution of dynamical observables, as a generalization, a perturbative method based on the Magnus expansion could be developed leading to approximate expressions for the time-evolution of the main dynamical variables in Heisenberg picture.

Appendix A Derivation of Eqs. (2) and (2)

From Hamiltonian Eq. (1) and Heisenberg equations of motion for $\hat{a}$ and $\hat{b}_{j}$ we have

[TABLE]

[TABLE]

The solution of Eq. (72) is

[TABLE]

by inserting this solution into Eq. (71), we find

[TABLE]

where the response function of the medium is defined by

[TABLE]

By taking the Laplace transform of both sides of Eq. (74) we obtain

[TABLE]

and by taking inverse Laplace transform we finally find

[TABLE]

with the hermitian conjugation

[TABLE]

where

[TABLE]

By inserting Eq. (77) into Eq. (73) we find

[TABLE]

Therefore,

[TABLE]

Appendix B Derivation of Eqs. (3)

[TABLE]

[TABLE]

Appendix C Derivation of Eq. (12)

[TABLE]

On the other hand Louisell

[TABLE]

by inserting Eq. (85) into Eq. (84) we deduce

[TABLE]

Appendix D Derivation of Eq. (19)

We have

[TABLE]

By inserting the definitions

[TABLE]

into the generating function defined by Eq. (87) we find

[TABLE]

where $\mbox{Tr}_{k}$ means taking trace over the Hilbert space of the kth oscillator in the environment and $z_{k}$ is the corresponding partition function

[TABLE]

Now we have

[TABLE]

Therefore, the generating function is

[TABLE]

By making use of Eq. (87) we finally find

[TABLE]

Appendix E Derivation of Eq. (26)

We have

[TABLE]

where $r_{min}=\max(s+n-k,0)$.

Appendix F Derivation of Eq. (30)

[TABLE]

Appendix G Derivation of Eq. (34)

In the absence of dissipation ($f_{k}=0,\,\,\,k=0,1,2,\cdots$), we have from Eq. (21)

[TABLE]

Therefore,

[TABLE]

Appendix H Vanishing of coherency in Eq. (43)

In large-time limit, we have $G(t)\rightarrow 0$ and for the non diagonal matrix elements of the reduced density matrix ($n\neq m$) we easily see that the maximum degree of $G(t)$ in denominator is smaller than its degree in numerator, since

[TABLE]

therefore, the non diagonal elements tend to zero and accordingly there is no coherency in the preferred number basis $|n\rangle$.

Appendix I Derivation of Eq. (57)

In large-time limit $G(t)\rightarrow 0$, so by setting $s+n-r=0$ and $s+m-r=0$ we find $n=m$, leading to

[TABLE]