No thermalization without correlations

TL;DR

This paper proves that quantum systems require persistent correlations with their environment to thermalize, challenging previous assumptions and proposing a new method to model baths that produce thermal states.

Contribution

It demonstrates that correlations are essential for thermalization in quantum systems and introduces a systematic way to construct translation invariant bath models.

Findings

Correlations are necessary for quantum thermalization.

A systematic procedure for bath model construction is proposed.

A quantum optical scheme for experimental validation is outlined.

Abstract

The proof of the long-standing conjecture is presented that Markovian quantum master equations are at odds with quantum thermodynamics under conventional assumptions of fluctuation-dissipation theorems (implying a translation invariant dissipation). Specifically, except for identified systems, persistent system-bath correlations of at least one kind, spatial or temporal, are obligatory for thermalization. A systematic procedure is proposed to construct translation invariant bath models producing steady states that well-approximate thermal states. A quantum optical scheme for the laboratory assessment of the developed procedure is outlined.

Equations

$$\frac{d}{dt}⟨{\hat{p}}_n⟩=-⟨\frac{\partial }{\partial {\hat{x}}_n}U(\hat{\mathbf{x}})⟩+⟨{\hat{F}}_n^{\text{mbox}fr}⟩,$$

$$\frac{d}{dt}⟨{\hat{x}}_n⟩=\frac{1}{m_n}⟨{\hat{p}}_n⟩,$$

$$\frac{\partial }{\partial t}\hat{\rho }=\mathcal{L}[\hat{\rho }],\text{leavevmode}\mathcal{L}={\mathcal{L}}_0+{\mathcal{L}}_{\text{mbox}rel},$$

$${\mathcal{L}}_0[\odot ]=\frac{i}{\hslash }[\odot ,\hat{H}],\text{leavevmode}\hat{H}=H(\hat{\mathbf{p}},\hat{\mathbf{x}})={\sum }_{n=1}^N\frac{{\hat{p}}_n^2}{2m_n}+U(\hat{\mathbf{x}}),$$

$${\mathcal{L}}_{\text{mbox}rel}=\sum _{k=1}^K{\mathcal{L}}_{{\hat{L}}_k}^{\text{mbox}lbd},\text{leavevmode}\text{leavevmode}{\mathcal{L}}_{\hat{L}}^{\text{mbox}lbd}[\hat{\rho }]\stackrel{\text{mbox}def}{=}\hat{L}\hat{\rho }{\hat{L}}^{†}-\frac{1}{2}({\hat{L}}^{†}\hat{L}\hat{\rho }+\hat{\rho }{L}^{†}\hat{L}),$$

$${\mathcal{L}}_{\text{mbox}rel}={\sum }_k{\mathcal{L}}_{{\hat{A}}_k}^{\text{mbox}lbd}+{\mathcal{L}}_{\text{mbox}aux}\text{mbox}with$$

$${\hat{A}}_k\stackrel{\text{mbox}def}{=}{e}^{-i{\mathbf{\kappa }}_k\hat{\mathbf{x}}}{\tilde{f}}_k(\hat{\mathbf{p}}),\text{leavevmode}\text{leavevmode}{\mathcal{L}}_{\text{mbox}aux}=-i[{\mathbf{\mu }}_{\text{mbox}aux}\hat{\mathbf{x}}+f_{\text{mbox}aux}(\hat{\mathbf{p}}),\odot ].$$

$${\mathcal{L}}_{\text{mbox}rel}[{\Psi }_0(\hat{\mathbf{p}})g(\hat{\mathbf{x}}){\Psi }_0(\hat{\mathbf{p}}{)}^{†}]=0$$

$$B(\mathbf{p},\mathbf{\lambda })={\int }_{-\infty }^{\infty }\dots {\int }_{-\infty }^{\infty }{e}^{i\mathbf{\lambda }\mathbf{x}}W(\mathbf{p},\mathbf{x}){\mathrm{d}}^N\mathbf{x}.$$

$$\forall \mathbf{p},\mathbf{\lambda }:B_{\theta }(\mathbf{p},\mathbf{\lambda })>0,\text{leavevmode}\text{leavevmode}{B}_{\theta }(\mathbf{p},-\mathbf{\lambda })=B_{\theta }(\mathbf{p},\mathbf{\lambda }),$$

$$\forall \mathbf{p}\ne 0,\mathbf{\lambda }\ne 0:B_{\theta }(\mathbf{p},\mathbf{\lambda })<B_{\theta }(0,0).$$

$$\frac{\partial }{\partial t}\varpi (\mathbf{p})\simeq \mathrm{Tr}[\delta (\mathbf{p}-\hat{\mathbf{p}}){\mathcal{L}}_0[\hat{\rho }]]+$$

$$\sum _{n,l}\frac{{\partial }^2{D}_{n,l}(\mathbf{p})\varpi (\mathbf{p})}{\partial p_n\partial p_l}-\sum _n\frac{\partial F_n^{\text{mbox}fr}(\mathbf{p})\varpi (\mathbf{p})}{\partial p_n}$$

$${\mathbf{F}}^{\text{mbox}fr}(\hat{\mathbf{p}})=-{\sum }_k\hslash {\mathbf{\kappa }}_k|{\tilde{f}}_k(\hat{\mathbf{p}}){|}^2,$$

$$D_{n,l}(\hat{\mathbf{p}})=\frac{{\hslash }^2}{2}{\sum }_k|{\tilde{f}}_k(\hat{\mathbf{p}}){|}^2{\kappa }_{k,n}{\kappa }_{k,l}.$$

$$\displaystyle\mathopen{}\mathclose{{}\left.\mathop{\tfrac{d}{dt}}\mathop{\langle e^{-\alpha\hat{H}}\rangle}\nolimits_{\theta}}\right|_{t{=}0}{=}0;\leavevmode\nobreak\ \mathopen{}\mathclose{{}\left|\mathop{\tfrac{d^{2}}{dt^{2}}}\mathop{\langle e^{-\alpha\hat{H}}\rangle}\nolimits_{\theta}}\right|_{t{=}0}{\to}\min\mbox{ for all }\alpha,$$

$${\tilde{f}}_k(p)=c_k{e}^{p\beta \hslash {\lambda }_k},\text{leavevmode}\text{leavevmode}{\lambda }_k={\kappa }_k\tanh (\frac{\hslash \omega }{4\theta }),$$

$$\displaystyle\mathop{\tfrac{d}{dt}}{\mathop{\langle\hat{H}\rangle}\nolimits_{\theta^{\prime}}}\big{|}_{t{=}0}{=}\tfrac{c_{k}^{2}}{\omega}\tilde{\gamma}^{\mbox{\scriptsize en}}_{k}(\theta^{\prime},\theta)(\mathop{\langle\hat{H}\rangle}\nolimits_{\theta}{-}\mathop{\langle\hat{H}\rangle}\nolimits_{\theta^{\prime}})\big{|}_{t{=}0},$$

$${\mathcal{B}}_{\kappa ,{\tilde{f}}^{\text{mbox}iso}}\stackrel{\text{mbox}def}{=}{\mathcal{L}}_{{\hat{A}}^{+}}^{\text{mbox}lbd}+{\mathcal{L}}_{{\hat{A}}^{-}}^{\text{mbox}lbd},\text{leavevmode}\text{leavevmode}{\hat{A}}^{\pm }=e^{\mp i\kappa \hat{x}}{\tilde{f}}^{\text{mbox}iso}(\pm \hat{p}).$$

$${\tilde{f}}_k(\hat{\mathbf{p}})=0\text{mbox}when\mathbf{p}{\mathbf{\kappa }}_k<0\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{mbox}(classicalmechanics).$$

$$\displaystyle\mathop{\tfrac{d}{dt}}{\mathop{\langle\hat{x}_{n}^{2}\rangle}\nolimits_{\theta}}\bigg{|}_{t{=}0}{=}\hbar^{2}\sum_{k}\mathop{\langle\bigl{|}\mathop{\tfrac{\partial}{\partial\hat{p}_{n}}}{\tilde{f}_{k}(\hat{\bm{p}})}\bigr{|}^{2}\rangle}\nolimits_{\theta}\bigg{|}_{t{=}0}{>}0$$

$$\kappa =\frac{{\omega }_{\text{mbox}l}}{c},\text{leavevmode}\text{leavevmode}{\tilde{g}}^{\text{mbox}iso}(p)={\tilde{c}}_1({\tilde{c}}_2^2+(p-{\tilde{c}}_3{)}^2{)}^{-\frac{1}{2}},\text{leavevmode}\text{leavevmode}{\tilde{c}}_k\in \mathbb{R}$$

$$\forall \mathbf{\delta x}:{\mathcal{R}}_{\mathbf{\delta x}}{\mathcal{L}}_{\text{mbox}rel}{\mathcal{R}}_{-\mathbf{\delta x}}={\mathcal{L}}_{\text{mbox}rel},$$

$${\mathcal{R}}_{\mathbf{\delta x}}=e^{-\frac{i}{\hslash }\mathbf{\delta x}\hat{\mathbf{p}}}\odot e^{\frac{i}{\hslash }\mathbf{\delta x}\hat{\mathbf{p}}}$$

$${\mathcal{L}}_{\text{mbox}rel}=\sum _{k,l_1,m_1,l_2,m_2}{c}_{k,l_1,m_1}{c}_{k,l_2,m_2}^{\ast }{\tilde{\mathcal{L}}}_{{\hat{B}}_{l_1,m_1},{\hat{B}}_{l_2,m_2}}^{\text{mbox}lbd},$$

$${\tilde{\mathcal{L}}}_{{\hat{A}}_1,{\hat{A}}_2}^{\text{mbox}lbd}\stackrel{\text{mbox}def}{=}{\hat{A}}_1\odot {\hat{A}}_2^{†}-\frac{1}{2}({\hat{A}}_1{\hat{A}}_2^{†}\odot +\odot {\hat{A}}_1{\hat{A}}_2^{†}).$$

$$\displaystyle{\cal L}_{\mbox{\scriptsize rel}}{=}\mathopen{}\mathclose{{}\left.\frac{1}{(2L)^{N}}\int_{-L}^{L}...\int_{-L}^{L}\,{\cal R}_{\bm{\delta x}}{\cal L}_{\mbox{\scriptsize rel}}{\cal R}_{-\bm{\delta x}}\mathrm{d}^{N}\mathop{\mathrm{\mathstrut d}}\nolimits\!\bm{\delta x}}\right|_{L{\to}\infty}{=}$$

$$\sum _l\sum _{m_1,m_2}{\tilde{c}}_{m_1,m_2}^{(l)}{\tilde{\mathcal{L}}}_{{\hat{B}}_{l,m_1},{\hat{B}}_{l,m_2}}^{\text{mbox}lbd},$$

$${\tilde{c}}_{m_1,m_2}^{(l)}=\sum _k{c}_{k,l,m_1}{c}_{k,l,m_2}^{\ast }.$$

$${\mathcal{L}}_{\text{mbox}rel}=\sum _{l,m}{\mathcal{L}}_{{\hat{A}}_{l,m}}^{\text{mbox}lbd},$$

$${\mathcal{L}}^{\text{mbox}G}=\sum _k{\mathcal{L}}_{{\hat{A}}_k^{\text{mbox}G}}^{\text{mbox}lbd},\text{leavevmode}\text{leavevmode}{\hat{A}}_k^{\text{mbox}G}={\mathbf{\mu }}_k\hat{\mathbf{x}}+{\tilde{f}}_k^G(\hat{\mathbf{p}})\text{leavevmode}\text{leavevmode}({\mathbf{\mu }}_k\in {\mathbb{R}}^N)$$

$$\displaystyle{\cal L}^{\mbox{\scriptsize lbd}}_{\hat{A}^{\mbox{\scriptsize G}}_{k}}{=}\mathopen{}\mathclose{{}\left.{\cal L}^{\mbox{\scriptsize lbd}}_{\hat{A}^{\mbox{\scriptsize G}}_{k,+}}{+}{\cal L}^{\mbox{\scriptsize lbd}}_{\hat{A}^{\mbox{\scriptsize G}}_{k,-}}}\right|_{\epsilon{\to}0}{-}i\frac{\hbar}{2}\bm{\mu}_{k}\mathopen{}\mathclose{{}\left[\mathop{\frac{\partial\tilde{{f}}^{G}_{k}(\hat{\bm{p}})}{\partial\hat{\bm{p}}}},\odot}\right],$$

$$\displaystyle\hat{A}^{\mbox{\scriptsize G}}_{k,\pm}{=}\tfrac{1}{\sqrt{2}}\mathopen{}\mathclose{{}\left(\tfrac{i}{\epsilon}{\pm}\tilde{{f}}^{G}_{k}(\hat{\bm{p}})}\right)e^{{\mp}i\epsilon\bm{\mu}_{k}\hat{\bm{x}}}.$$

$$\forall g({\mathbf{x}}^{\prime }):{\mathcal{L}}_{\text{mbox}rel}[\int g({\mathbf{x}}^{\prime }){\mathcal{R}}_{\mathbf{\delta x}}[{\hat{\rho }}_0]d^N{\mathbf{x}}^{\prime }]=0.$$

Full text

I Introduction

A stochastic interaction of a quantum system with a bath brings up the term $\hat{\bm{F}}^{\mbox{\scriptsize fr}}$ in the relations for time-dependent expectation values of system momenta $\hat{\bm{p}}{=}\{\hat{p}_{1},\ldots,\hat{p}_{N}\}$ and positions $\hat{\bm{x}}{=}\{\hat{x}_{1},\ldots,\hat{x}_{N}\}$:

[TABLE]

[TABLE]

where $U(\hat{\bm{x}})$ is a potential energy operator and $m_{k}$ are effective masses. In this Letter, we study the case where $\hat{\bm{F}}^{\mbox{\scriptsize fr}}{=}\hat{\bm{F}}^{\mbox{\scriptsize fr}}(\hat{\bm{p}})$ is position-independent. In this form, Eqs. (1) apply to many quantum phenomena including the translational motion of an excited atom in vacuum Sonnleitner et al. (2017), Brownian motion in a dilute background gas Accardi et al. (2003), light-driven processes in semiconductor, nanoplasmonic and optomechanical systems Auffèves et al. (2013); Milburn and Walls (2008); Barchielli and Vacchini (2015), superconducting currents O’Connell (1996), quantum ratchets Reimann (2002), energy transport in low-dimensional systems Dhar (2008), dynamics of chemical reactions Kuhn et al. , two-dimensional vibrational spectroscopy and NMR signals Steffen and Tanimura (2000); Tanimura (2006) as well as more exotic entirely quantum dissipative effects Rezek (2010); Kardar and Golestanian .

The term $\hat{\bm{F}}^{\mbox{\scriptsize fr}}(\hat{\bm{p}})$ in Eqs. (1) admits a simple classical interpretation as friction acting on effective particles moving in a potential $U(\bm{x})$. Such classical dynamics are described by the familiar Langevin, Drude and Fokker-Plank models when the system-bath interactions are treated as

(i) memoryless (Markovian) and

(ii) translation invariant (position-independent).

However, we will show that these two assumptions are at odds with quantum thermodynamics. Specifically, we will prove a long-standing no-go conjecture that completely positive111Positivity of quantum evolution guarantees satisfaction of the Heisenberg uncertainty principle at all times. It was shown that the requirements for positivity and complete positivity coincide for some quantum systems including a harmonic oscillator Kohen et al. (1997). Markovian translation-invariant quantum dynamics obeying Eqs. (1) cannot thermalize.

The no-go conjecture was demonstrated by Lindblad as early as in 1976 Lindblad (1976a) for a quantum harmonic oscillator with a Gaussian damping222The Gaussian damping corresponds to ${\cal L}_{\mbox{\scriptsize rel}}{=}{\cal L}^{\mbox{\scriptsize lbd}}_{\bm{\mu}\hat{\bm{x}}{+}\bm{\eta}\hat{\bm{p}}}$ ($\bm{\mu},\bm{\eta}\in\mathbb{C}^{N}$) in Eq. (2a) and can be cast to form (3), as shown in Appendix A). The original paper Lindblad (1976a) deals with one-dimensional case. The multidimensional extension can be found e.g. in Dodonov and Manko (1985).. Subsequently his particular result was extended to a general quantum system under the weight of mounting numerical evidence, however without proof. The no-go conjecture is de-facto incorporated in all popular models such as the Redfield theory Redfield (1957), the Gaussian phase space ansatz of Yan and Mukamel Yi J. Y. and S. Mukamel (1988), the master equations of Agarwal Agarwal (1971), Caldeira-Leggett Caldeira and Leggett (1983), Hu-Paz-Zhang Ford and O’Connell (2005), and Louisell/Lax Louisell and Walker (1965), and the semigroup theory of Lindblad Lindblad (1976b) along with specialized extensions in different areas of physics and chemistry. These models break either one of assumptions (i) and (ii) or the complete positivity of quantum evolution (see Kohen et al. (1997); Yan et al. (2000); Bondar et al. (2016) for detailed reviews, note errata Bondar et al. ). This circumstance is a persistent source of controversies (see e.g. the discussions Wiseman and Munro (1998); O’Connell (2001); Vacchini (2001) of original works Gao (1997); Vacchini (2000)). The matters were further complicated by the discovery that the free Brownian motion $U(\hat{\bm{x}}){=}0$ circumvents the conjecture Vacchini (2002) (we will identify the full scope of possible exceptions below).

The no-go result challenges studies of the long-time dynamics of open systems. On the one hand, model’s thermodynamic consistency is undermined by assumptions (i) and (ii). On other hand, the same assumptions open opportunities to simulate large systems that are otherwise beyond the reach. Specifically, the abandonment of Markovianity entails a substantial overhead to store and process the evolution history. The value of assumption (ii) can be clarified by the following example. Consider the re-thermalization of a harmonic oscillator coupled to a bath (represented by a collection of harmonic oscillators) after displacement from equilibrium by, e.g., an added external field, a varied system-bath coupling, or interactions between parts of a compound system. To account for such a displacement without assumption (ii), one needs to self-consistently identify the equilibrium position for each bath oscillator, re-thermalize the bath and modify the system-bath couplings accordingly. In practice, this procedure is intractable without gross approximations that lead to either numerical instabilities or physical inaccuracies. Choosing among a polaron-transformation-based method, Redfield, and Förster (hopping) models of quantum transfer epitomizes this dilemma Jang et al. (2013).

Remarkably, assumption (ii) enables to model the displaced state equilibrium by simply adjusting the potential energy $\hat{U}$. Fig. 1a shows that without this assumption the potential adjustment yields steady state $\hat{\rho}_{\mbox{\scriptsize st}}$ significantly different from the canonical equilibrium $\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}{\propto}e^{{-}\frac{\hat{H}}{\theta}}$, where $\theta{=}k_{\mbox{\scriptsize B}}T$ and $\hat{H}$ is system Hamiltonian.

Motivated by these arguments, we propose in this Letter a general recipe to construct approximately thermalizable bath models under assumptions (i) and (ii). Fig. 1 illustrates this recipe in application to the above example. The resulting mismatch between $\hat{\rho}_{\mbox{\scriptsize st}}$ and $\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}$ is small, especially at high temperatures and in the weak system-bath coupling limit. (The calculations details will be explained below.)

It will be shown elsewhere that the proposed recipe is capable of accurately accounting for electronic and spin degrees of freedom. We found it helpful in reservoir engineering and optimal control problems. Moreover, the resulting bath models are realizable in the laboratory and can be used for coupling atoms and molecules nonreciprocally Zhdanov et al. . However, the scope of our recipe is limited by the applicability of assumptions (i) and (ii) and, therefore, cannot encompass strongly correlated systems (as in the case of Anderson localization Nandkishore and Huse (2015)).

II The key results

Starting by formalizing the problem, we write the general master equation that accounts for memoryless system-bath interactions and ensures positivity of the system density matrix $\hat{\rho}$ at all times Lindblad (1976b):

[TABLE]

[TABLE]

where $\odot$ is the substitution symbol defined, e.g., in Grishanin and Zadkov (2000). The superoperator $\cal{L}_{\mbox{\scriptsize rel}}$ accounts for system-bath couplings responsible for the friction term $\hat{\bm{F}}^{\mbox{\scriptsize fr}}$ in Eq. (1a) and depends on a set of generally non-Hermitian operators $\hat{L}_{k}$. Based on theorems by A. Holevo Holevo (1995, 1996), B. Vacchini Petruccione and Vacchini (2005); Vacchini (2009); Vacchini and Hornberger (2009) has identified the following criterion of translational invariance for the $\cal{L}_{\mbox{\scriptsize rel}}$:

Lemma 1** **(The justification is in

Appendix A).

Any translationally invariant superoperator ${\cal L}_{\mbox{\scriptsize rel}}$ of the Lindblad form (2c) can be represented as

[TABLE]

[TABLE]

*where $\bm{\kappa}_{k}$ and $\bm{\mu}_{\mbox{\scriptsize aux}}$ are $N$-dimensional real vectors, $\tilde{f}_{k}$ are complex-valued functions and $f_{\mbox{\scriptsize aux}}$ is real-valued 333 The Gaussian dissipators ${\cal L}^{\mbox{\scriptsize lbd}}_{\bm{\mu}_{k}\hat{\bm{x}}{+}\tilde{{f}}^{G}_{k}(\hat{\bm{p}})}$ $(\bm{\mu}_{k}{\in}\mathbb{R}^{N})$ can be reduced to the form Eq. (3) as a limiting case $\bm{\kappa}_{k}{\to}0$, as shown in Appendix A. The generalized unitary drift term ${\cal L}_{\mbox{\scriptsize aux}}$ accounts for ambiguity of the separation of the quantum Liouvillian ${\cal L}$ in Eq. (2a) into Hamiltonian and relaxation parts. * .The converse holds as well.

The primary findings of this work are summarized in the following two no-go theorems.

No-go theorem 1**.**

Let $\mathop{|{\Psi_{0}}\rangle}\nolimits$ be the ground state (or any other eigenstate of $\hat{H}$), such that ${\mathop{\langle{\Psi_{0}}|\hat{\bm{p}}\,|{\Psi_{0}}\rangle}\nolimits}{=}0$, and which momentum-space wavefunction $\Psi_{0}(\bm{p}){=}{\mathop{\langle{\bm{p}}|{\Psi_{0}}\rangle}\nolimits}$ is nonzero almost everywhere, except for some isolated points. Then, no translationally invariant Markovian process of form (2) and (3) can steer the system to $\mathop{|{\Psi_{0}}\rangle}\nolimits$.

The idea of the proof, whose details are given in Appendix B, is to show that the state $\hat{\rho}_{0}{=}\mathop{|{\Psi_{0}}\rangle\langle{\Psi_{0}}|}\nolimits$ can be the fixed point of superoperator $e^{t\cal L}$ only if ${\cal L}_{\mbox{\scriptsize rel}}{\equiv}0$. First, note that the linearity and translation invariance of the dissipator (3) imply that ${\cal L}_{\mbox{\scriptsize rel}}[\int g(\bm{x}^{\prime})e^{{-}\frac{i}{\hbar}\bm{x}^{\prime}\hat{\bm{p}}}\hat{\rho}_{0}e^{\frac{i}{\hbar}\bm{x}^{\prime}\hat{\bm{p}}}d^{N}\bm{x}^{\prime}]{=}0$ for any function $g(\bm{x}^{\prime})$. This equation can be equivalently rewritten as

[TABLE]

using the identities $e^{{-}\frac{i}{\hbar}\bm{x}^{\prime}\hat{\bm{p}}}\mathop{|{\Psi_{0}}\rangle}\nolimits{=}\sqrt{2\pi\hbar}\Psi_{0}(\hat{\bm{p}})\mathop{|{\bm{x}^{\prime}}\rangle}\nolimits$ and $\int g(\bm{x}^{\prime})\mathop{|{\bm{x}^{\prime}}\rangle\langle{\bm{x}^{\prime}}|}\nolimits d^{N}\bm{x}^{\prime}{=}g(\hat{\bm{x}})$, where $\mathop{|{\bm{x}^{\prime}}\rangle}\nolimits$ is the eigenstate of position operator: $\hat{x}_{k}\mathop{|{\bm{x}^{\prime}}\rangle}\nolimits{=}x^{\prime}_{k}\mathop{|{\bm{x}^{\prime}}\rangle}\nolimits$. Let us choose $g(\bm{x}){=}e^{-i\bm{\lambda}\bm{x}}$, where $\bm{\lambda}$ is an arbitrary real vector, and move to the right the $\hat{\bm{x}}$-dependent terms in the lhs of Eq. (4) using the commutation relations $e^{-i\tilde{\bm{\lambda}}\hat{\bm{x}}}\hat{\bm{p}}=(\hat{\bm{p}}{+}\hbar\tilde{\lambda})e^{-i\tilde{\bm{\lambda}}\hat{\bm{x}}}$ with $\tilde{\bm{\lambda}}{=}\bm{\lambda},\pm\bm{\kappa}_{k}$. This rearrangement brings Eq. (4) to the form $\tilde{G}_{\bm{\lambda}}(\hat{\bm{p}})e^{-i\bm{\lambda}\hat{\bm{x}}}{=}0$ (note that all the operators of form $e^{\pm i\tilde{\bm{\kappa}}_{k}\hat{\bm{x}}}$ expectedly cancel out owing to translation invariance of ${\cal L}_{\mbox{\scriptsize rel}}$). The last equality can be satisfied only if the function $\tilde{G}_{\bm{\lambda}}(\bm{p})$ vanishes identically for all $\bm{p}$ and $\bm{\lambda}$. However, careful inspection of Appendix B shows that the latter happens only if ${\cal L}_{\mbox{\scriptsize rel}}{=}0$.

The statement of the 1-st no-go theorem can be strengthened for a special class of quantum systems. Let $B(\bm{p},\bm{\lambda})$ be the Blokhintsev function Blokhintzev (1940), which is related to Wigner quasiprobability distribution $W(\bm{p},\bm{x})$ as

[TABLE]

No-go theorem 2**.**

Suppose that the Blokhintsev function $B_{\theta}(\bm{p},\bm{\lambda})$ of the thermal state $\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}{\propto}e^{-\frac{\hat{H}}{\theta}}$ characterized by temperature $k_{\mbox{\scriptsize B}}T{=}\theta$ is such that

[TABLE]

[TABLE]

Then, no translationally invariant Markovian process (2) and (3) can asymptotically steer the system to $\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}$.

The proof of this theorem is given in Appendix C and generally follows the same logic as the outlined proof of the 1-st no-go theorem. Using Eq. (5) and the familiar formula for the thermal state Wigner function Bartlett and Moyal (1949), it is easy to check that the criteria (6) are satisfied for any $\theta$ in the case of a quadratic potential $U$. This means that the Lindblad’s original conclusion on inability to thermalize the damped harmonic oscillator using the Gaussian friction term ${\cal L}_{\mbox{\scriptsize rel}}{=}{\cal L}^{\mbox{\scriptsize lbd}}_{\bm{\mu}\hat{\bm{x}}{+}\bm{\eta}\hat{\bm{p}}}$ is equally valid for all Markovian translationally invariant dissipators.

Corollary 2.1**.**

No translationally invariant Markovian process of form (2) and (3) can steer the quantum harmonic oscillator into a thermal state of form $\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}{\propto}e^{{-}\frac{\hat{H}}{\theta}}$.

III Practical implications of the no-go theorems

In classical thermodynamics, the bath is understood as a constant-temperature heat tank “unaware” of a system of interest. However, the no-go theorems indicate that system-bath correlations of at least one kind – spatial or temporal – become obligatory for thermalization once quantum mechanical effects are taken into account. These correlations break the bath translation invariance or Markovianity assumptions, respectively.

Nevertheless, in the view of computational advantages outlined above, it is desirable to incorporate these assumptions into the master equations (2) and (3). Now we are going to introduce the recipe to construct such models with a minimal error in the thermal state. In order to proceed, note that in the limit $(\hbar\bm{\kappa}_{k})^{2}{\ll}\mathop{\langle\hat{\bm{p}}^{2}\rangle}\nolimits$ Eqs. (2) and (3) reduce to the familiar Fokker-Planck equation

[TABLE]

for the momentum probability distribution $\varpi(\bm{p}){=}\mathop{\rm{Tr}}[\delta(\bm{p}{-}\hat{\bm{p}})\hat{\rho}]$. The friction forces ${\bm{F}}^{\mbox{\scriptsize fr}}$ in Eq. (7) as well as Eq. (1a) have the form

[TABLE]

whereas the momentum-dependent diffusion operator is

[TABLE]

Equations (8) and (9) can be satisfied by different sets of $\bm{\kappa}_{k}$ and $\tilde{f}_{k}(p)$. We will exploit this non-uniqueness to reduce the system-bath correlation errors. Our strategy is reminiscent to the familiar way of making density functional calculations practical via error cancellation in approximated exchange-correlation functionals. We shall demonstrate the generic procedure by considering a one-dimensional oscillator with the Hamiltonian $\hat{H}{=}\tfrac{m}{2}\hat{p}^{2}{+}\frac{m\omega^{2}}{2}\hat{x}^{2}$ (here the dimension subscript $n$ is omitted for brevity). Corollary 2.1 implies that ${\cal L}_{\mbox{\scriptsize rel}}[\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}]{\neq}0$ and $\hat{\rho}_{\mbox{\scriptsize st}}{\neq}\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}$ for any $\theta$, where \hat{\rho}_{\mbox{\scriptsize st}}{=}\mathopen{}\mathclose{{}\left.\hat{\rho}}\right|_{t\to\infty} is the actual fixed point of the evolution operator $e^{t{\cal L}}$. However, the net discrepancies can be reduced by imposing the following thermal population conserving constraint:

[TABLE]

where $\mathop{\langle\odot\rangle}\nolimits_{\theta}(t){=}\mathop{\rm{Tr}}[\odot e^{t{\cal L}}[\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}]]$. This constraint can be intuitively justified when the characteristic decay rates are much smaller than the typical transition frequencies, such that the dissipation can be treated perturbatively. Since the term ${\cal L}_{\mbox{\scriptsize rel}}[\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}]$ generates only rapidly oscillating off-diagonal elements in the basis of $\hat{H}$, Eq. (10) ensures that the first-order perturbation vanishes on average for the exact thermal state: $\lim_{t\to{\infty}}\tfrac{1}{t}\int_{0}^{t}e^{\tau{\cal L}_{0}}{\cal L}_{\mbox{\scriptsize rel}}e^{(t-\tau){\cal L}_{0}}[\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}]\mathop{\mathrm{\mathstrut d}}\nolimits\!\tau{=}0$.

In the case of the driftless dissipation ${\cal L}_{\mbox{\scriptsize aux}}{=}0$, Eq. (10) is satisfied by the following functions $\tilde{f}_{k}(p)$ in Eq. (3):

[TABLE]

where $\beta{=}(m\hbar\omega)^{-1}$ and the constants $c_{k}$ should be chosen to satisfy Eq. (8). The corresponding dissipator (3) reproduces the familiar microphysical model of quantum Brownian motion (see e.g. Eq. (16) in Ref. Vacchini (2002)) in the limit $\bm{\kappa}{\to}0$, $\omega{\to}0$. Furthermore, the resulting dynamics tends to decrease (increase) the average system energy $\mathop{\langle\hat{H}\rangle}\nolimits_{\theta}$ if its initial temperature $\theta^{\prime}$ is higher (lower) than $\theta$:

[TABLE]

where \tilde{\gamma}^{\mbox{\scriptsize en}}_{k}(\theta^{\prime},\theta){=}2\omega{\beta\hbar^{2}\kappa_{k}\lambda_{k}\exp\mathopen{}\mathclose{{}\left({{\beta\hbar^{2}\lambda_{k}^{2}}{\coth(\frac{\hbar\omega}{2\theta^{\prime}})}}}\right)}{>}0.

Equation (12) suggests that $\hat{\rho}_{\mbox{\scriptsize st}}$ is close to $\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}$. This conclusion is supported by the simulations presented in Fig. 2a for the isotropic dissipator ${\cal L}_{\mbox{\scriptsize rel}}{=}{\cal B}_{\kappa,\tilde{f}^{\mbox{\scriptsize iso}}}$,

[TABLE]

One can see that the high-quality thermalization is readily achieved by tuning the free parameters $c_{k}$ and $\kappa_{k}$ even in the strong dissipation regime.

To understand the result (11), note that the terms ${\cal L}^{\mbox{\scriptsize lbd}}_{\hat{A}_{k}}$ in Eq. (3) represent independent statistical forces $\mathop{\langle{-}\hbar\bm{\kappa}_{k}|\tilde{f}_{k}(\hat{\bm{p}})|^{2}\rangle}\nolimits$ contributing to the net friction $\mathop{\langle\hat{\bm{F}}^{\mbox{\scriptsize fr}}\rangle}\nolimits$. In classical mechanics, such forces at $\theta{=}0$ steer the system to the state of rest by acting against the particles’ momenta, hence

[TABLE]

However, clipping the functions (11) according to Eq. (14) introduces significant errors, as displayed by dotted curves in Fig. 2a. Thus, the “endothermic” tails of $\tilde{f}_{k}(\hat{\bm{p}})$ at $\bm{p}\bm{\kappa}_{k}{>}0$ break the thermalization in the classical case, but reduce errors in the quantum mechanical treatment. To clarify this counterintuitive observation, note that the physical requirement $\mathop{\frac{d}{dt}}\mathop{\langle\hat{O}\rangle}\nolimits_{\theta}{=}0$ for any observable $\hat{O}$ in the thermodynamic equilibrium $\hat{\rho}_{\mbox{\scriptsize st}}{=}\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}$ is violated by the master equations (2) and (3) due to the no-go theorems, i.e.,

[TABLE]

in the driftless case ${\cal L}_{\mbox{\scriptsize aux}}{=}0$. The inequality (15) provides further evidence for the no-go theorems and is the hallmark of the “position diffusion”, a known artifact in the quantum theory of Brownian motion Vacchini and Hornberger (2009).

According to Eq. (15), \mathop{\tfrac{d}{dt}}\mathop{\langle\hat{x}^{2}\rangle}\nolimits_{\theta}\big{|}_{t{=}0} is sensitive to smoothness of $\tilde{f}_{k}(\bm{p})$. Specifically, the rhs of Eq. (15) is exploded by any highly oscillatory components of $\tilde{f}_{k}(\bm{p})$ and diverges if $\tilde{f}_{k}(\bm{p})$ is discontinuous. This entirely quantum effect is the origin of poor performance of the clipped solutions (14) seen in Fig. 2a. Equation (15) uncovers unavoidable errors in the potential energy. The optimal solutions (11) enforce error cancellation \mathop{\tfrac{d}{dt}}\mathop{\langle\tfrac{\hat{p}^{2}}{2m}\rangle}\nolimits_{\theta}\big{|}_{t{=}0}{=}{-}\mathop{\tfrac{d}{dt}}\mathop{\langle U(\hat{x})\rangle}\nolimits_{\theta}\big{|}_{t{=}0} between kinetic and potential energies leaving the total energy intact \mathop{\tfrac{d}{dt}}\mathop{\langle\hat{H}\rangle}\nolimits_{\theta}\big{|}_{t{=}0}{=}0. In fact, the error cancellation is achieved with a large class of physically feasible functions $\tilde{f}_{k}(\bm{p})$ that may substantially differ from the solutions (11) everywhere but the region of high probability density $\varpi(p){=}\mathop{\rm{Tr}}[\delta{(\hat{p}{-}p)\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}}]$ (however, note the remark in Appendix D). This is illustrated in Fig. 2a by dashed curves overlapping with solid curves.

The master equations (2) and (3) provide accurate non-perturbative description of collisions with a background gas of atoms or photons Poyatos et al. (1996); Vacchini, B. (2001); Vacchini (2009); Hornberger and Vacchini (2008); Barchielli and Vacchini (2015). Hence, the above theoretical conclusions can be directly tested in the laboratory using well-developed techniques, e.g., the setup shown in Fig. 2b. Here a two-level atom is subject to two orthogonally polarized counterpropagating monochromatic nonsaturating laser fields of the same amplitude $\cal E$ and frequency $\omega_{\mbox{\scriptsize l}}$. We show in Appendix D that the translational motion of the atom can be modeled using Eq. (2) with an isotropic friction term of form ${\cal L}_{\mbox{\scriptsize rel}}{=}{\cal B}_{\kappa,\tilde{g}^{\mbox{\scriptsize iso}}}$. Here

[TABLE]

and the parameters $\tilde{c}_{k}$ can be tuned by $\cal E$ and $\omega_{\mbox{\scriptsize l}}$.

Now we are ready to clarify why the deviations from canonical equilibrium increase with $|\kappa|$ in Fig. 2a. The parameters $\hbar|\kappa|$ and $\tilde{g}^{\mbox{\scriptsize iso}}(p)^{2}$ in Eq. (16) can be regarded as the change of atomic momentum after absorption of a photon and the absorption rate. The case of small $\hbar|\kappa|{\ll}\sqrt{\mathop{\langle\hat{p}^{2}\rangle}\nolimits}$ and large $\tilde{g}^{\mbox{\scriptsize iso}}(p)^{2}$ implies tiny and frequent momentum exchanges subject to the central limit theorem. The net result is a velocity-dependent radiation pressure with vanishing fluctuations. The opposite case of large $\hbar|\kappa|{\gg}\sqrt{\mathop{\langle\hat{p}^{2}\rangle}\nolimits}$ and small $\tilde{g}^{\mbox{\scriptsize iso}}(p)^{2}$ is the strong shot noise limit, where the stochastic character of light absorption is no longer averaged out, notably perturbing the thermal state. Note that a similar interpretation applies to quantum statistical forces in Ref. Vuglar et al. .

The dissipative model (2) and (3) with optimized parameters (11) is further analyzed in Fig. 1 using the same parameters as in Fig. 2a. Both Figs. 1 and 2a indicate that thermalization can be modeled for a wide range of recoil momenta \hbar\kappa\in\mathopen{}\mathclose{{}\left({-}(\hbar\sqrt{\beta})^{-1},(\hbar\sqrt{\beta})^{-1}}\right) and the higher the temperature, the better the accuracy. Thus, Eqs. (8) and (9) enable to simulate a variety of velocity dependences of friction and diffusion.

Finally, Fig. 1a benchmarks such simulations against the commonly used quantum optical master equation (QOME) Gardiner and Zoller (2004) defined by Eq. (2c) with $K{=}2$, $\hat{L}_{1}{=}\sqrt{2\Gamma\omega}({1{-}e^{{-}\frac{\hbar\omega}{\theta}}})^{{-}\frac{1}{2}}\hat{a}$, $\hat{L}_{2}{=}\sqrt{2\Gamma\omega}({e^{\frac{\hbar\omega}{\theta}}}{-}1)^{{-}\frac{1}{2}}\hat{a}^{\dagger}$, where $\hat{a}$ is the harmonic oscillator annihilation operator. For a correct comparison, the parameters of both models are adjusted to ensure identical decay rates in Eq. (12). Systematic errors in our model and QOME are comparable for the equilibrium displacements $\Delta x_{0}{\sim}\hbar\beta^{{-}\frac{1}{2}}$ at zero temperature and $\Delta x_{0}{\sim}0.1\hbar\beta^{{-}\frac{1}{2}}$ for $\theta{\sim}\hbar\omega$. For low-frequency molecular vibrational modes ($m{\sim}10^{4}$ atomic units, $\omega{\sim}200$ cm*-1*), these shifts are of order $0.4$ Å and $0.04$ Å, respectively, which are in the range of typical molecular geometry changes due to optical excitations or liquid environments. We found the displacement-independent errors in the model (2) and (3) to be very important for quantum control via reservoir engineering. Furthermore, the same feature can also be exploited for engineering the mechanical analogs of nonreciprocal optical couplings Metelmann and Clerk (2015) and energy-efficient molecular quantum heat machines Zhdanov et al. . These subjects will be explored in a forthcoming publication.

Acknowledgements.

We thank to Alexander Eisfeld for valuable discussions and drawing our attention to Refs. Holevo (1995, 1996); Petruccione and Vacchini (2005); Vacchini (2009); Vacchini and Hornberger (2009). We are grateful to the anonymous referees for important suggestions and insightful critic. T. S. and D. V. Zh. thank the National Science Foundation (Award number CHEM-1012207 to T. S.) for support. D. I. B. is supported by AFOSR Young Investigator Research Program (FA9550-16-1-0254).

SUPPLEMENTAL MATERIAL

Introduction

This supplemental material is organized as follows. In the Sections A, B and C we give the proofs of Lemma 1, first, and second no-go theorems, respectively. Finally, the supporting mathematical derivations for the Doppler cooling model briefly discussed in the main text are provided in Section D.

The Roman numbers in parentheses refer everywhere to the equations in the main text of the letter.

Appendix A The proof of lemma 1

The property of the translational invariance can be formulated as

[TABLE]

where

[TABLE]

is the superoperator of translational shift: $\forall g(\hat{\bm{x}}):{\cal R}^{\intercal}_{\bm{\delta x}}[g(\hat{\bm{x}})]{=}g(\hat{\bm{x}}{+}\bm{\delta x})$.

With the help of the canonical commutation relations, any operator $\hat{L}_{k}{=}L_{k}(\hat{\bm{p}},\hat{\bm{x}})$ can expanded in the series $\hat{L}_{k}{=}\sum_{l,m}c_{k,l,m}\hat{B}_{l,m}$, where $\hat{B}_{l,m}{=}e^{{-}i\bm{\kappa}_{l}\hat{\bm{x}}}g_{m}(\hat{\bm{p}})$ and the functions $g_{m}(\bm{p})$ constitute a set of (not necessarily orthogonal) basis functions. Using this expansion, any superoperator of form ${\cal L}_{\mbox{\scriptsize rel}}{=}\sum_{k}{\cal L}^{\mbox{\scriptsize lbd}}_{\hat{L}_{k}}$ can be rewritten as

[TABLE]

where

[TABLE]

It follows from Eq. (17) that if ${\cal L}_{\mbox{\scriptsize rel}}$ is translationally invariant then it should satisfy the identity

[TABLE]

where the Hermitian matrices $\tilde{c}^{(l)}$ are defined as

[TABLE]

Let us substitute in Eq. (21) the matrices $\tilde{c}^{(l)}$ with their Jordan decomposition $\tilde{c}^{(l)}{=}\tilde{u}^{(l)}\tilde{\gamma}^{(l)}{\tilde{u}^{(l)}}^{\dagger}$, where $\tilde{u}^{(l)}$ is unitary and $\tilde{\gamma}^{(l)}$ is diagonal. The result is

[TABLE]

where $A_{l,m}{=}\tilde{f}_{l,m}(\hat{\bm{p}})e^{{-}i\bm{\kappa}_{l}\hat{\bm{x}}}$ and $\tilde{f}_{l,m}(\hat{\bm{p}}){=}\sqrt{\tilde{\gamma}^{(l)}_{m,m}}\times\sum_{m^{\prime}}\tilde{u}^{(l)}_{m^{\prime},m}g_{m^{\prime}}(\hat{\bm{p}})$. Finally, note that Eq. (23) can be cast into the form (3) by replacing the compound index $\{l,m\}$ with the single consecutive index $k$. The lemma is proven.

• Remark 1.

In this work, the Gaussian (continuous) translationally invariant dissipators of form

[TABLE]

are treated as the limiting case of Eq. (3) with $\bm{\kappa}_{k}{=}\epsilon\bm{\mu}_{k}{\to}0$. Specifically, one can verify by direct calculation that

[TABLE]

• Remark 2.

The translation invariance criterion is generalized to non-Markovian dynamics in Ref. Gasbarri et al. .

Appendix B The proof of no-go theorem 1 (by contradiction)

Suppose that some eigenstate $\mathop{|{\Psi_{0}}\rangle}\nolimits$ of Hamiltonian $\hat{H}$ is also the fixed point of the quantum Liouvillian $\cal L$ defined by Eqs. (2) and (3). Since ${\cal L}_{\mbox{\scriptsize rel}}$ is assumed translation invariant, it should commute with spatial shift superoperator (18) for any $\bm{\delta x}$. Hence, ${\cal L}_{\mbox{\scriptsize rel}}[{\cal R}_{\bm{\delta x}}[\hat{\rho}_{0}]]={\cal R}_{\bm{\delta x}}[{\cal L}_{\mbox{\scriptsize rel}}[\hat{\rho}_{0}]]{=}0$, where $\hat{\rho}_{0}{=}\mathop{|{\Psi_{0}}\rangle\langle{\Psi_{0}}|}\nolimits$. Furthermore, the linearity of ${\cal L}_{\mbox{\scriptsize rel}}$ implies that

[TABLE]

Equation (27) can be further simplified using the identity

[TABLE]

where $\mathop{|{\bm{x}^{\prime}}\rangle}\nolimits$ is the eigenstate of position operator: $\hat{x}_{k}\mathop{|{\bm{x}^{\prime}}\rangle}\nolimits{=}x_{k}^{\prime}\mathop{|{\bm{x}^{\prime}}\rangle}\nolimits$, ${\mathop{\langle{\bm{x}^{\prime\prime}}|{\bm{x}^{\prime}}\rangle}\nolimits}{=}\delta(\bm{x}^{\prime\prime}{-}\bm{x}^{\prime})$. The validity of Eq. (28) can be verified by comparing the wavefunctions in momentum representation corresponding to its left and right sides. Identities (28) and $\int g(\bm{x}^{\prime})\mathop{|{\bm{x}^{\prime}}\rangle\langle{\bm{x}^{\prime}}|}\nolimits d^{N}\bm{x}^{\prime}{=}g(\hat{\bm{x}})$ allow to equivalently rewrite Eq. (27) as

[TABLE]

Consider the case $g(\bm{x}){=}g_{\bm{\lambda}}(\bm{x}){=}e^{-i\bm{\lambda}\bm{x}}$, where $\bm{\lambda}$ is some real $N$-dimensional vector. Note that the operator ${\cal L}_{\mbox{\scriptsize rel}}[\hat{w}_{g}]$ then includes the explicit dependence on coordinate operators $\hat{\bm{x}}$ only in forms of matrix exponentials $e^{-i\bm{\lambda}\hat{\bm{x}}}$, $e^{\pm i\bm{\kappa}_{k}\hat{\bm{x}}}$ and commutators $[\hat{x}_{k},\odot]$. Using the commutation relation $e^{-i\tilde{\bm{\lambda}}\hat{\bm{x}}}\hat{\bm{p}}{=}(\hat{\bm{p}}{+}\hbar\tilde{\lambda})e^{-i\tilde{\bm{\lambda}}\hat{\bm{x}}}$ with $\tilde{\bm{\lambda}}{=}\bm{\lambda},\pm\bm{\kappa}_{k}$, it is possible to group out the momentum and coordinate operators in ${\cal L}_{\mbox{\scriptsize rel}}[\hat{w}_{g}]$ and rewrite the condition (29) as:

[TABLE]

where

[TABLE]

and

[TABLE]

Note that the cancellation of all the matrix exponentials $e^{\pm i\bm{\kappa}_{k}\hat{\bm{x}}}$ at the rhs of Eq (30) is the consequence of the translation invariance of ${\cal L}_{\mbox{\scriptsize rel}}$. Condition (30) implies that

[TABLE]

and $\forall\bm{p},\bm{p}^{\prime}{\in}\mathbb{R}^{N}:G(\bm{p},\bm{p}^{\prime}){=}0$ (except a possible zero measure subset of points $\{\bm{p},\bm{p}^{\prime}\}$ where $\Psi_{0}(\hat{\bm{p}})\Psi_{0}(\hat{\bm{p}}^{\prime})^{\dagger}{=}0$). In particular, this means that

[TABLE]

Equality (35) can be satisfied only if $\forall k:F_{k}(\bm{p}){\propto}$const, i.e., if $\tilde{f}_{k}(\bm{p}){=}c_{k}\frac{\Psi_{0}(\bm{p}{-}\hbar\bm{\kappa}_{k})}{\Psi_{0}(\bm{p})}$, where $c_{k}$ is some real constant. Substitution of this expression and $\bm{\lambda}{=}\bm{0}$ into Eq. (31) gives \tilde{G}_{\bm{0}}(\bm{p}){=}\sum_{k}c_{k}^{2}\mathopen{}\mathclose{{}\left({|{\Psi_{0}(\bm{p})}|^{2}{-}|\Psi_{0}(\bm{p}{-}\hbar\bm{\kappa}_{k})}|^{2}}\right){+}\hbar\bm{\mu}_{\mbox{\scriptsize aux}}\mathop{\tfrac{\partial|{\Psi_{0}}(\bm{p})|^{2}}{\partial\bm{p}}} 444Using (25), it is straightforward to deduce that the corresponding summand for Gaussian dissipator (24) takes form $\frac{1}{2}\hbar^{2}\bm{\mu}\mathop{\frac{\partial}{\partial\bm{p}}}(\bm{\mu}\mathop{\frac{\partial}{\partial\bm{p}}}|{\Psi_{0}(\bm{p})}|^{2})$. The resulting contribution in the lhs of Eq. (B) is $|\bm{\mu}|^{2}\hbar^{2}$. . Multiplication of the both sides of Eq. (34) by $\bm{p}^{2}$ and subsequent integration over $\bm{p}$ gives:

[TABLE]

According to our assumption, ${\mathop{\langle{\Psi_{0}}|\hat{\bm{p}}\,|{\Psi_{0}}\rangle}\nolimits}{=}\bm{0}$555The equality ${\mathop{\langle{\Psi_{0}}|\hat{\bm{p}}\,|{\Psi_{0}}\rangle}\nolimits}{=}\bm{0}$ holds for any non-degenerate eigenstate of the time-reversal invariant Hamiltonian (2b). Hence, Eq. (B) implies that $\sum_{k}c_{k}^{2}|\bm{\kappa}_{k}|^{2}{=}0$. This equality holds only if $\forall k:\bm{\kappa}_{k}{=}\bm{0}$. However, in this case all functions $\tilde{f}_{k}(\bm{p}){=}c_{k}$ reduce to constants, so that ${\cal L}_{\mbox{\scriptsize rel}}{=}0$.

This result completes the proof.

Appendix C The proof of no-go theorem 2 (by contradiction)

Denote as $\Psi_{k}(\bm{p})$ and $E_{k}$ $(k{=}0,...,\infty)$ the momentum-space wavefunction and energy of the $k$-th eigenstate $\mathop{|{\Psi_{k}}\rangle}\nolimits$ of the Hamiltonian $\hat{H}$. The thermal state $\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}$ can be expressed in these notations as

[TABLE]

Suppose that there exists such relaxation superoperator of form (3) that ${\cal L}_{\mbox{\scriptsize rel}}[\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}]{=}0$. Owing to assumed linearity and translational invariance of ${\cal L}_{\mbox{\scriptsize rel}}$, the thermal state $\hat{\rho}^{\mbox{\scriptsize th}}_{\theta}$ should satisfy the relation similar to (27):

[TABLE]

where ${\cal R}_{\bm{\delta x}}$ is the spatial shift superoperator defined by Eq. (18). With the help of relation (28), one can apply to Eq. (38) the same procedure as was used to derive the equality (29) from Eq. (27). The result is

[TABLE]

where

[TABLE]

Consider the case $g(\bm{x}){=}g_{\bm{\lambda}}(\bm{x}){=}e^{{-}i\bm{\lambda}\bm{x}}$, where $\bm{\lambda}$ is some real $N$-dimensional vector. The result of application of ${\cal L}_{\mbox{\scriptsize rel}}$ to $\hat{w}_{\theta,g_{\bm{\lambda}}}$ can be represented after some algebra as

[TABLE]

where

[TABLE]

In derivation of (C) the identity

[TABLE]

was used which follows directly from the definition (5) of the Blokhintsev function.

Eqs. (39) and (41) require that

[TABLE]

and hence $\forall\bm{\lambda}:\bar{G}_{2}(\bm{\lambda}){=}\int_{{-}\infty}^{\infty}\ldots\int_{{-}\infty}^{\infty}\mathop{\mathrm{\mathstrut d}}\nolimits^{N}\!\bm{p}\,G_{2}(\bm{p},\bm{\lambda}){=}0$, where

[TABLE]

The last equality in (46) is obtained assuming that $B_{\theta}(\bm{p},{-}\bm{\lambda}){=}B_{\theta}(\bm{p},\bm{\lambda})$ (see Eq. (6a)). It is easy to check that the integrations over all terms in the last line of (46) cancel out, so that

[TABLE]

According to the assumption (6a), the integrand in (47) is nonnegative. Moreover, $\bar{G}_{2}(\bm{\lambda}){=}0$ iif $\forall k:\tilde{f}_{k}(\bm{p}){=}c_{k}{=}$const. Hence, the expression (C) for $G_{1}(\bm{p},\bm{\lambda})$ can be simplified as

[TABLE]

Note that the terms ${\cal L}^{\mbox{\scriptsize lbd}}_{\hat{A}_{k}}$ in Eq. (3) with $\tilde{f}_{k}(\bm{p}){=}$const will have non-trivial effect only if $\bm{\kappa}_{k}{\neq}0$666In the case of Gaussian dissipator (24) Eq. (48) reduces to

$\displaystyle G_{1}(\bm{p},\bm{\lambda})=\frac{1}{2}\hbar^{2}\sum_{k,m,n}\mu_{k,n}\mu_{k,m}\mathop{\frac{\partial^{2}}{\partial p_{n}p_{m}}}B(\bm{p},\bm{\lambda}).$

(48*)

By assumption (6b), the quadratic form $\mathop{\frac{\partial^{2}}{\partial p_{n}p_{m}}}B(\bm{p},\bm{\lambda})$ in (48*) is negative-definite at $\{\bm{p},\bm{\lambda}\}{=}\{\bm{0},\bm{0}\}$. Hence, $G_{1}(\bm{0},\bm{0}){<}0$, which contradicts Eq. (45) and completes the proof for this case. . However, it follows from (6b) that in this case $G_{1}(\bm{0},\bm{0}){<}0$ which contradicts Eq. (45). The theorem is proven.

Appendix D Testing the model (2) and (3) in the laboratory

In this section, we provide the detailed analysis of the Doppler cooling example introduced in the main text (see Fig. 2b in the main text) and prove that the cooling mechanism is the quantum friction of form (13).

In the proposed setup an atom is subject to two orthogonally polarized counterpropagating beams of the same field amplitude $\cal E$ and carrier frequency $\omega_{\mbox{\scriptsize l}}$ (hereafter in this section we will omit the subscript l for shortness since it will not cause any ambiguity). We assume that $\omega$ is close to the frequency $\omega_{\mbox{\scriptsize a}}$ of the transition ${\mathop{\mathtt{g}}}{{\mathop{\leftrightarrow}}}{\mathop{\mathtt{e}}}$ between the ground ${\mathop{\mathtt{g}}}$ and degenerate excited ${\mathop{\mathtt{e}}}$ electron states of $s$- and $p$-symmetries, respectively. Let $d$ be the absolute value of the transition dipole moment and $\gamma$ be the excited state spontaneous decay rate.

For the spatial arrangement depicted in Fig. 2b the translation motion of the atom along $x$-axis is coupled to the field-induced electron dynamics since each absorbed or coherently emitted photon changes the $x$-component of atomic momentum hereafter denoted as $p$. Furthermore, we will assume that the spontaneous decay does not affect the $x$-component of atomic momentum. The latter condition can be achieved using, e.g., an arrangement shown in Fig. 3.

The master equation which describes this coupled dynamics can be written within the rotating wave approximation in the form (2) with

[TABLE]

and

[TABLE]

Here $\xi_{k}(t){=}{-}\frac{1}{2}\vec{d}_{k}\vec{\cal E}_{k}(t)$, where $\vec{d}_{1}$ and $\vec{d}_{2}$ are the transition dipole moments associated with the $s{\to}p_{z}$ and $s{\to}p_{y}$ electronic transitions into degenerate electronically excited sublevels ${\mathop{\mathtt{e}}}_{1}$ and ${\mathop{\mathtt{e}}}_{2}$, respectively, and $\vec{\cal E}_{k}(t)$ is the slowly varying complex amplitude of the associated field component. The remaining notations are defined in the main text.

The mean value of any observable of form $\hat{O}{=}f(\hat{p},\hat{x})$ can be written in Heisenberg representation as:

[TABLE]

where we define:

[TABLE]

The symbol $\stackrel{{\scriptstyle\Rightarrow}}{{{\cal T}}}$ in (52) denotes the chronological ordering superoperator which arranges operators in direct (inverse) time order for $t{>}t_{0}$ ($t{<}t_{0}$). Let us also define the following notations for the interaction representation generated by arbitrary splitting ${\cal L}(t){=}{\cal L}_{0}+{\cal L}_{1}(t)$:

[TABLE]

where the interaction Liouvillian reads

[TABLE]

In the case ${\cal L}_{0}^{\prime}{=}\frac{-i}{\hbar}[\frac{\hat{p}^{2}}{2m}{-}\hbar\omega_{\mbox{\scriptsize a}}\mathop{|{{\mathop{\mathtt{g}}}}\rangle\langle{{\mathop{\mathtt{g}}}}|}\nolimits,\odot]$ the associated interaction liouvillian (54) in the rotating wave approximation takes the form:

[TABLE]

where

[TABLE]

and $\Delta{=}\omega{-}\omega_{\mbox{\scriptsize a}}$ is detuning of carrier frequency of radiation from atomic resonance in the case of system at rest. Repeated application of the transformation (53) to (55) with ${\cal L}_{0}^{\prime\prime}{=}{\cal L}_{\mbox{\scriptsize rel}}{=}\gamma\sum_{n{=}1}^{2}{\cal L}^{\mbox{\scriptsize lbd}}_{\mathop{|{{\mathop{\mathtt{g}}}}\rangle\langle{{\mathop{\mathtt{e}}}_{n}}|}\nolimits}$ leads to expression:

[TABLE]

so that

[TABLE]

where $\hat{P}_{{\mathop{\mathtt{g}}}}{=}\mathop{|{{\mathop{\mathtt{g}}}}\rangle\langle{{\mathop{\mathtt{g}}}}|}\nolimits$ and the last equality is due to the exponential damping of excited states populations induced by relaxation superoperator (50). Let us consider the evolution $\hat{O}(t)$ generated by the superoperator $\mathop{{{\cal U}^{{{\cal L}_{\mbox{\scriptsize I}}^{\prime\prime}}^{\intercal}}_{t+\delta t,t}}}$:

[TABLE]

Integrands in Eq. (62) include the terms oscillating at frequencies $|\Delta{\pm}\frac{k\mathop{\langle\hat{p}\rangle}\nolimits}{m}|$. In sequel we will consider the so-called weak-field limit when these oscillations are rapid relative to the characteristic timescales of the relevant processes, so that the contributions of the associated terms asymptotically vanish. In this limit, the second term in rhs of Eq. (62) disappears. The remaining terms constitute two decoupled evolution equations for the reduced density matrices $f_{{\mathop{\mathtt{x}}}}(\hat{p},\hat{x},t{+}\delta t){=}{\mathop{\langle{{\mathop{\mathtt{x}}}}|\hat{O}(t)\,|{{\mathop{\mathtt{x}}}}\rangle}\nolimits}$ (${\mathop{\mathtt{x}}}{=}{\mathop{\mathtt{g}}},{\mathop{\mathtt{e}}}$):

[TABLE]

The explicit form of $\cal G$ is irrelevant in view of Eq. (61). The first two terms in the curly brackets in Eq. (63) can be transformed as

[TABLE]

[TABLE]

[TABLE]

where $\Delta_{1}(p){=}\Delta{-}\frac{\kappa(p{+}\frac{\hbar\kappa}{2})}{m}.$ The extra displacements $\frac{\hbar\kappa}{m}\tau_{n}$ in the $x$-dependencies of $f_{{\mathop{\mathtt{g}}}}$ in Eqs. (65) account for the change of the velocity of atom after the photon absorption. These displacements are typically very small compared to the characteristic scales of spatial change of the function $f_{{\mathop{\mathtt{g}}}}$ and can be neglected. With this approximation, the exponentials and functions $f_{{\mathop{\mathtt{g}}}}$ in Eqs. (65) commute, which allows to write:

[TABLE]

[TABLE]

[TABLE]

where

[TABLE]

Substitution of approximations (66) into (63) gives:

[TABLE]

where

[TABLE]

[TABLE]

[TABLE]

Eq. (68) allows to calculate the averaging in (61) within the reduced Hilbert space which involves only the translational degree of freedom:

[TABLE]

Here $\hat{\rho}_{0}^{\mbox{\scriptsize red}}{=}\mathop{\rm{Tr}}[\hat{\rho}]_{\mbox{\scriptsize el}}$ whereas $\mathop{\rm{Tr}}[\odot]_{\mbox{\scriptsize el}}$ and $\mathop{\rm{Tr}}[\odot]_{\mbox{\scriptsize spatial}}$ denote the partial traces over the electronic and translational subsystems.

The dissipator (70) reduces to the isotropic friction of form (13) provided that

[TABLE]

It is easy to verify that this condition is realized in two important cases.

D.1 Weak coherent laser driving

In this regime, $\xi_{1}(t){=}\xi_{2}(t){=}\xi{=}$const, and there exists such $\delta t$ in the range of applicability of the second-order expansion (62) that $\delta t{\gg}\gamma^{-1}$. Thence, the integrals in (67) can be easily computed, which gives:

[TABLE]

Note what the Hamiltonian $\hat{H}_{\mbox{\scriptsize eff}}$ describes the effect of the optical quadratic Stark shift which also can induce the effective potential forces on the system in the case of spatially non-uniform fields $\xi{=}\xi(x)$.

D.2 Incoherent driving

Suppose that the the atom is illuminated by the two classical light sources with the equal spectral densities $I(\omega)$ at the atomic site and having coherence times in the range $\Delta_{1}^{-1}(p){\ll}t_{\mbox{\scriptsize coh}}{\ll}\gamma^{-1}$. In this case, $\xi_{1}(t)$ and $\xi_{2}(t)$ represent the uncorrelated stationary stochastic processes. This allows one to choose such $\delta t$, that $\gamma^{-1}{\gg}\delta t{\gg}t_{\mbox{\scriptsize coh}}$, and calculate the integrals in Eqs. (67) neglecting the terms $\frac{\gamma}{2}$ in the exponents, which gives

[TABLE]

where $I(\omega)$ is the spectral density of each beam. Also, here we assumed equal transition dipole momenta: $d{=}|\vec{d}_{1}|{=}|\vec{d}_{2}|$.

• Remark.

The setup sketched in Fig. 3 as well as in Fig. 2 of the main text in principle can be used to measure both the momenta and positions of the environmental photons by registering the scattered photons and the position of atom. This implies that there must exist the fundamental restrictions on the physically admissible shapes and smoothness of profiles $\tilde{f}^{\mbox{\scriptsize iso}}(p)$ and, more generally, on admissible forms of operators $\hat{L}_{k}$ in Eq. (2c), that would prevent these measurements from violating the Heisenberg uncertainty principle. The detailed analysis of implications of this important observation is way beyond the scope of this paper and will be the subject of future work.