# Objective trajectories in hybrid classical-quantum dynamics

Jonathan Oppenheim  Department of Physics and Astronomy, University College
London, Gower Street, London WC1E 6BT, United Kingdom  Carlo Sparaciari
Department of Physics and Astronomy, University College London, Gower Street,
London WC1E 6BT, United Kingdom  Barbara Šoda  Department of Physics and
Astronomy, University College London, Gower Street, London WC1E 6BT, United
Kingdom  Dept. of Physics, University of Waterloo, Waterloo, Ontario, Canada
Perimeter Institute for Theoretical Physics, Waterloo, Ontario, Canada
Zachary Weller-Davies  Department of Physics and Astronomy, University College
London, Gower Street, London WC1E 6BT, United Kingdom

###### Abstract

Consistent dynamics which couples classical and quantum degrees of freedom
exists, provided it is stochastic. This provides a way to study the back-
reaction of quantum fields on space-time which does not suffer from the
pathologies of the semi-classical equations. Here we introduce several toy
models in which to study hybrid classical-quantum evolution, including a qubit
coupled to a particle in a potential, and a quantum harmonic oscillator
coupled to a classical one. We present an unravelling approach to calculate
the dynamics. Unlike the purely quantum case, the trajectories (or histories)
of this unravelling can be unique, conditioned on the classical degrees of
freedom. This provides a potential solution to the “measurement problem” of
quantum theory; quantum systems become classical because they interact with a
genuinely classical field.

##  I Introduction

We are often interested in the dynamics of composite systems where one system
can be considered classical while the other must be treated quantum
mechanically. In quantum thermodynamics and quantum chemistry, we often have
small molecules interacting with a large thermal reservoir which can be
treated classically. In measurement theory the quantum system interacts with a
macroscopic device which can be considered classical. In gravity, macroscopic
objects such as evaporating black holes radiate thermally, and we imagine
there is a regime where we can treat the black hole space-time classically
even though the radiation must still be described by quantum field theory.

There is a lot of debate in the literature on whether one can consistently
couple quantum systems and classical ones, and many ways of doing so which
although useful in some regimes, are pathological Bohr and Rosenfeld (1933);
DeWitt and Rickles (2011); DeWitt (1962); Eppley and Hannah (1977); Caro and
Salcedo (1999); Salcedo (1996); Sahoo (2004); Terno (2006); Salcedo (2012);
Barceló et al. (2012); Marletto and Vedral (2017). Many proposals for such
dynamics Boucher and Traschen (1988); Diósi et al. (2000); Anderson (1995);
Kapral and Ciccotti (1999); Kapral (2006a); Prezhdo and Kisil (1997); Kapral
(2006b) are not completely positive, and since both the density matrix of a
quantum system and the phase-space density of a classical system are positive
functions or matrices, such maps lead to negative probabilities, and are at
best an approximation to the true dynamics. Other proposals, such as semi-
classical gravity, are non-linear and thus don’t respect the statistical
nature of the density matrix Page and Geilker (1981); Gisin (1989). However,
evolution laws introduced in Blanchard and Jadczyk (1993) are consistent and
have been used to study the collisionless Boltzmann equation Alicki and
Kryszewski (2003), discussed in the context of measurement theory Diósi
(2014), and have been proposed as a fundamental theory coupling quantum fields
to classical space-time Poulin and Preskill (2017); Oppenheim (2018). Because
we don’t have a quantum theory of gravity, the prospect of understanding how
classical space-time and quantum fields interact is exciting, whether the
dynamics are fundamental or merely effective. In addition to the potential to
better understand black hole physics, an effective theory opens up the
possibility of gaining new insights into cosmology, since the quantum nature
of vacuum fluctuations are important in structure formation in an expanding
universe.

We will not discuss gravity in the present work. Rather the aim here is to
introduce basic techniques, and study some simple systems in order to build up
an intuition for how classical systems couple to quantum ones. In particular,
we shall see that such coupling not only leads to decoherence of the quantum
system Hall and Reginatto (2005); Poulin (2017); Oppenheim (2018), but also a
”collapse” of the wavefunction, in the sense that the quantum system jumps to
a pure state which can be determined from the classical degrees of freedom. In
this sense, we will see that the state of the system can be considered to have
some objectivity, unlike the purely quantum case where a decomposition of the
system’s density matrix into pure states (or evolution into a history) is not
unique Kent (1998); Dowker and Kent (1995). We also study another feature of
hybrid dynamics, namely a trade-off between decoherence and diffusion. Quantum
systems which have long coherence times necessarily have high diffusion in the
classical degrees of freedom. Here, we shall exhibit systems which have this
feature. That there is necessarily a trade-off between decoherence of the
quantum system and diffusion of the classical one, will be shown in Oppenheim
et al. (2020).

We shall here consider a qubit and a quantum harmonic oscillator which are
coupled to classical degrees of freedom are described by a phase-space
manifold $\mathcal{M}$, and for simplicity in the following we take
$\mathcal{M}\equiv\mathbb{R}^{2n}$, the space of position and momentum $z=q,p$
for $n$ particles. The quantum degrees of freedom are described by a Hilbert
space $\mathcal{H}$, which might be finite or infinite dimensional. Given the
Hilbert space, we define the set of positive semi-definite operators whose
trace is at most unity as $S_{\leq 1}(\mathcal{H})$. The object describing the
state of this composite system at a given time is the map
$\varrho:\mathcal{M}\rightarrow S_{\leq 1}(\mathcal{H})$ such that, when the
quantum degrees of freedom are traced out, returns a valid probability
distribution over the phase-space. We will call this the hybrid-density.
Mathematically, we have that the distribution
$\text{prob}(z)=\mathrm{Tr}\left[{\varrho(z)}\right]$ is such that
$\text{prob}(z)\geq 0$ for all $z\in\mathcal{M}$, and it is normalized
$\int_{\mathcal{M}}\mathrm{d}z\,\text{prob}(z)=1$. From the above property it
is easy to show that, if the classical degrees of freedom are traced out, the
resulting state $\rho=\int_{\mathcal{M}}\mathrm{d}z\,\varrho(z)$ is a valid
quantum state (i.e., a positive semi-definite operator with unit trace). The
simplest such system is the hybrid qubit, whose state can be written as

$\varrho(q,p)=\begin{pmatrix}u_{0}(q,p)&c(q,p)\\\
c^{\star}(q,p)&u_{1}(q,p)\end{pmatrix}$

If we integrate over phase space, we obtain the density matrix of the qubit,
while if we take the trace of the density matrix, we are left with the phase-
space density of a classical particle.

A consistent and general form of dynamics for these composite classical-
quantum systems is given by a master equation, first introduced in Ref.
Blanchard and Jadczyk (1993) and here reported as Eq. (2), which couples the
quantum and classical degrees of freedom. This master equation posses several
necessary properties; it is linear in the state of the composite classical-
quantum system, and it is completely positive and trace preserving, so that
it’s action preserves the statistical interpretation of $\varrho(q,p)$. A
valid quantum-classical state is transformed to another valid quantum-
classical state, with the master equation preserving the separation between
classical and quantum degrees of freedom. The dynamics is necessarily
stochastic, and it is not generally the case that it arises from the dynamics
of a larger quantum system. This makes it attractive as a theory of gravity –
for a full discussion of this, along with a discussion on the consistency of
quantum-classical dynamics we refer the reader to Oppenheim (2018). Ref.
Oppenheim (2018) also introduced a particular form of Eq. (2) which has the
desired classical limit of Hamiltonian evolution, here presented as Eq. (4).

In this paper, we are interested in understanding the main features of the
dynamics generated by this master equation. In order to study these
properties, we generalise a well-known tools from open quantum systems, known
as the _unravelling_ of the master equation Belavkin (1990); Dalibard et al.
(1992); Gardiner et al. (1992). In open quantum systems Alicki and Lendi
(1987); Breuer and Petruccione (2002), one studies the evolution of a quantum
system coupled to an external environment. The most general Markovian
evolution of the quantum system given some reasonable assumptions is the GKSL
or Lindblad equation Lindblad (1976); Gorini et al. (1976), whose action over
the set of quantum states forms a semi-group Davies (1974). The dynamics of a
quantum state $\rho$ are given by

$\displaystyle\frac{\partial\,\rho}{\partial
t}=-i\left[H,\rho\right]+\sum_{\alpha,\beta}\lambda^{\alpha\beta}L_{\alpha}\rho
L^{\dagger}_{\beta}-\frac{1}{2}\sum_{\alpha,\beta}\lambda^{\alpha\beta}\left(z\right)\left\\{{L^{\dagger}_{\beta}L_{\alpha}},{\rho}\right\\}_{+}$
(1)

with $H$ the Hamiltonian of the system, $\left\\{{\cdot},{\cdot}\right\\}_{+}$
is the anti-commutator and $L_{\alpha},L_{\beta}\in\mathcal{B}(\mathcal{H})$
are Lindblad operators. Such dynamics arise, for example, when a quantum
system is weakly coupled to a very large environment.

The unravelling is a technique that allows us to numerically simulate the
evolution of a quantum system generated by the GKSL equation Belavkin (1990);
Dalibard et al. (1992); Gardiner et al. (1992); Carmichael (1993); Mølmer and
Castin (1996); Mølmer et al. (1993). Instead of considering the evolution of
the full density matrix describing the system, in the unravelling one consider
a single pure state, which is evolved according to a stochastic dynamics.
Multiple paths in the Hilbert space are thus generated by the stochastic
dynamics, and by averaging over the paths we recover the evolution of the
quantum system. This procedure is computationally advantageous compared to
evolving the full density matrix, in particular when the dimension of the
system is large. Moreover, the unravelling technique provides a different
perspective on the evolution of open quantum systems; indeed, the evolution of
the system can be understood as generated by a continuous, deterministic
dynamics (when the system does not interact with the environment), scattered
by discrete jumps of the wave-function (when system and environment interact),
which occur stochastically in time Gardiner et al. (1992).

In the quantum case, there is generally not a unique unravelling for the
master equation, so one cannot regard a particular trajectory as actually
happening. Indeed, the decomposition of the dynamics in terms of Lindblad
operators $L_{\alpha}$ and a Hamiltonian $H$ in Eq. (1), is not in general
unique. One can rewrite the $L_{\alpha}$ in terms of some other basis of
operators. This presents a problem Kent (1998); Dowker and Kent (1995) for
interpretations of quantum theory which take the trajectories as giving a
theory of microscopic state reduction Wiseman (1996), and likewise for the
decohering histories approachOmnes (1992); Griffiths (1993); Gell-Mann and
Hartle (1993); Diósi et al. (1995). In contrast, by extending the tools of
unravelling to the master equation which generates the dynamics of classical-
quantum systems, we shall find that trajectories can be unique conditioned on
the classical degrees of freedom (c.f. Paz and Zurek (1993) where the
environment plays this role). This imbues the unravelled trajectories with an
ontological significance and allows one to interpret the state of the system
as ”collapsing” to a particular state. This collapse of the wave-function is
here generated by the interaction between classical and quantum degrees of
freedom rather than the introduction of an ad-hoc field as in spontaneous
collapse models Ghirardi et al. (1986); Pearle (1989); Ghirardi et al. (1990);
Bassi et al. (2013). Classical systems can be in definite states, and we find
that a classical-quantum interaction causes the quantum system to inherit this
property. This provides a new resolution of the measurement problem of quantum
theory – quantum systems become classical because they interact with another
system which is itself intrinsically classical (e.g. space-time as in
Oppenheim (2018)).

However, one needn’t make such an interpretational commitment. The unravelling
can merely be considered to be a calculational tool which allows us to better
simulate the dynamics of some simple hybrid systems and we present an
algorithm and code for this purpose. We study both the hybrid qubit and hybrid
oscillator numerically and analytically, and explore a relation between the
diffusion of the classical degrees of freedom in phase space and the rate of
decoherence of the quantum degrees of freedom. Interestingly, this relation
provides a way to experimentally test whether the master equation can be used
to construct a ”post-” quantum theory of gravity Oppenheim (2018). We also
study the energy of the composite system, defined in terms of the sum of the
classical and quantum Hamiltonian operators. Because this composite energy is
not the generator of time-translations, there is no reason to expect it to be
conserved Banks et al. (1984) – Noether’s theorem doesn’t apply if the
dynamics is not unitary. Nonetheless, insofar as one can have a classical-
quantum Hamiltonian, one might want to keep violations of its conservation
small. There are various mechanisms one can employ for this purpose which we
will briefly mention in the discussion. In the gravitational setting, the more
natural question to consider is whether the constraints are preserved in time,
a discussion which we reserve for the future Oppenheim and Weller-Davies
(2020).

The paper is structured as follow. In Sec. II we briefly review the hybrid
master equation describing the evolution of classical-quantum systems, and we
extend the unravelling technique to this setting. In Sec. III we study two
different toy models where the classical degrees of freedom interact with a
finite-dimensional quantum system (a qubit). We use the unravelling technique
to simulate the hybrid dynamics, and to identify its main features. In Sec. IV
we consider the case in which the quantum system is infinite-dimensional (a
harmonic oscillator), and we repeat our analysis in this setting, with a
special focus on the decoherence of this system. We conclude in Sec. V, and
suggest various ways these models could naturally be extended, for example, by
having the quantum dynamics always dependent on phase space degrees of freedom
or by adding in friction terms which would control the diffusion of the
classical degrees of freedom in phase space.

##  II The classical-quantum master equation and its unravelling

In this section, we show that the most general master equation governing the
dynamics of a classical-quantum (CQ) system can be unravelled, allowing us to
efficiently simulate the evolution of different finite-dimensional hybrid
systems. Recasting the dynamics under the unravelling formalism provides a
useful perspective for understanding the evolution of CQ systems. Indeed, from
this point of view the hybrid system evolves continuously in its classical and
quantum degrees of freedom, and its evolution is interrupted at random times
by jumps in both the classical phase-space and the Hilbert space. It is worth
noting that the jumps are here due to the interactions between the classical
and the quantum system, and they are not, as in the case of open quantum
systems, due to the presence of an external environment. Although one can
always consider the state of a hybrid classical-quantum system to be the
restriction of a pure quantum state on part of an enlarged Hilbert-space, one
is not guaranteed that there will be a unitary dynamics on the entire system
whose restricted dynamics is of the required form.

Another difference between the unravelling for hybrid dynamics and the one
used in open system dynamics concerns the objectiveness of the jumps. In the
case of open quantum systems, a given dynamics can be obtained from an entire
family of jumps Breuer and Petruccione (2002). There are usually many ways in
which one can decompose the dynamics in terms of jumps in Hilbert space, and
so one cannot think of the unravelling as representing actual trajectories
that the system makes. In contrast, for some classes of hybrid dynamics, each
classical jumps in phase-space is uniquely defined and one can imagine
continuously monitoring the classical system through time without disturbing
the system. We can thus think of the classical trajectories as being ones
which objectively occur.

Before we can extend the unravelling technique to case of classical-quantum
dynamics, we need to formally introduce the class of hybrid systems we
consider, and the form of the master equation describing their evolution. The
most general Markovian CQ dynamics in the case of bounded operators
$L_{\alpha}$ is given by Blanchard and Jadczyk (1993); Diósi (2014); Alicki
and Kryszewski (2003)

$\displaystyle\frac{\partial\,\varrho(z,t)}{\partial t}=$
$\displaystyle-i\left[H(z),\varrho(z,t)\right]$
$\displaystyle+\int\mathrm{d}\Delta\,\sum_{\alpha,\beta}W^{\alpha\beta}\left(z|z-\Delta\right)\,L_{\alpha}\,\varrho(z-\Delta,t)\,L^{\dagger}_{\beta}-\frac{1}{2}\sum_{\alpha,\beta}W^{\alpha\beta}\left(z\right)\left\\{{L^{\dagger}_{\beta}L_{\alpha}},{\varrho(z,t)}\right\\}_{+},$
(2)

where $z=\left(q_{1},p_{1},q_{2},p_{2},\ldots\right)\in\mathbb{R}^{2n}$ is the
vector of phase space coordinates for $n$ systems, and $\varrho(z,t)$ is a
classical-quantum state at time $t$. The Hamiltonian $H(z)$, which appears in
the commutator with the CQ state, controls the unitary evolution of the
quantum system, and in principle depends on the classical degrees of freedom.
For each $\alpha$ and $\beta$, the rate
$W^{\alpha\beta}\left(z|z-\Delta\right)$ governs the transition of the
classical degrees of freedom from $z-\Delta$ to $z$, as well as the jump of
the quantum state due to the map $L_{\alpha}\cdot L^{\dagger}_{\beta}$. For
the master equation to preserve the norm of the quantum state (once the
classical degrees of freedom have been traced out), we need to require that

$W^{\alpha\beta}\left(z\right)=\int\mathrm{d}\Delta\,W^{\alpha\beta}\left(z+\Delta|z\right)\quad\forall\,\alpha,\beta.$
(3)

The hybrid master equation (2) is completely positive (CP) over the classical
and quantum degrees of freedom, and trace-preserving (TP) in the sense that
the normalisation condition $\int dz\mathrm{Tr}\left[{\varrho(z)}\right]=1$ is
preserved. Furthermore, its structure is such that the separation between
classical and quantum degrees of freedom is preserved i.e. $\varrho(z)$ is
always a subnormalised density matrix. Finally, the master equation is clearly
linear in the CQ state. All these properties ensure that a CQ state is mapped
into another CQ state by the dynamics generated by Eq. (2).

Although Eq. (2) gives the general form of the dynamics, we are interested in
evolution which in some limit reproduces classical Hamiltonin dynamics. This
was done in Oppenheim (2018) by expanding $W^{\alpha\beta}$ in terms of
moments of $\delta$ and demanding that the 1st moment reduce to the Poisson
bracket of the state with some Hamiltonian. The particular form we will
explore here, is given by taking matrix
$W^{\alpha\beta}\left(z|z-\Delta\right)$ such that

$\int\mathrm{d}\Delta\,\sum_{\alpha}W^{\alpha\beta}\left(z|z-\Delta\right)\,L_{\alpha}\,\varrho(z-\Delta,t)\,L^{\dagger}_{\beta}=\sum_{\alpha\beta}\frac{1}{\tau_{\alpha\beta}}e^{\tau_{\alpha\beta}\,\left\\{{h^{\alpha\beta}(z)},{\,\cdot\,}\right\\}}\,L_{\alpha}\,\varrho(z,t)\,L^{\dagger}_{\beta},$
(4)

where $\tau_{\alpha\beta}>0$ can be understood as a rate,
$\left\\{{\cdot},{\cdot}\right\\}$ is the Poisson bracket, and the functions
$h^{\alpha}(z)$ are associated with the CQ interaction Hamiltonian

$\displaystyle H(z)=\sum h^{\alpha\beta}(z)L_{\beta}^{\dagger}L_{\alpha}$ (5)

The resulting master equation is

$\frac{\partial\,\varrho(z,t)}{\partial
t}=-i\left[H(z),\varrho(z,t)\right]+\sum_{\alpha\beta}\frac{1}{\tau_{\alpha\beta}}\left(e^{\tau_{\alpha\beta}\,\left\\{{h^{\alpha}(z)},{\,\cdot\,}\right\\}}\,L_{\alpha}\,\varrho(z,t)\,L^{\dagger}_{\beta}-\frac{1}{2}\left\\{{L^{\dagger}_{\beta}L_{\alpha}},{\varrho(z,t)}\right\\}_{+}\right).$
(6)

and we see that to 0’th order in the $\tau_{\alpha\beta}$, this master
equation reproduces Hamiltonian dynamics in the appropriate limit. Namely, if
we expand the exponential in the right hand side of the equation for small
values of the parameters $\tau_{\alpha\beta}$’s, we obtain

$\displaystyle\frac{\partial\,\varrho(z,t)}{\partial t}=$
$\displaystyle-i\left[H(z),\varrho(z,t)\right]+\sum_{\alpha}\frac{1}{\tau_{\alpha\beta}}\left(L_{\alpha}\,\varrho(z,t)\,L^{\dagger}_{\beta}-\frac{1}{2}\left\\{{L^{\dagger}_{\beta}L_{\alpha}},{\varrho(z,t)}\right\\}_{+}\right)$
$\displaystyle+\sum\\{h^{\alpha\beta}(z),L_{\alpha}\varrho(z)L^{\dagger}_{\beta}\\}+\cdots$
(7)

where the parameters $\tau_{\alpha\beta}$’s play the role of the relaxation
rates The first term of Eq. (7) is the dynamics of the quantum system, which
can in principle depend on the classical degrees of freedom. In this way, it
represents the evolution of the quantum system as controlled by the classical
one. The next two terms are purely Lindbladian dynamics, and in the limit that
$\tau_{\alpha\beta}\rightarrow 0$, it can lead to rapid decoherence of the
quantum system. The final term is the interaction Hamiltonian dynamics we
desire, namely under the trace, it is
$\mathrm{Tr}\left[{\\{H,\varrho\\}}\right]$. We can think of it as the
counter-part to the commutator in the sense that it gives the back-action of
the quantum system on the classical one, while the commutator reflects the
influence of the classical system on the quantum one. Indeed, in the classical
limit, we expect $-i[H(z),\varrho]$ to act as the Poisson bracket with respect
to the degrees of freedom of the quantum system, while in this limit one can
verify that ${\\{h^{\alpha\beta},L_{\alpha}\varrho L^{\dagger}_{\beta}\\}}$
are the additional terms of the Poisson bracket with respect to the degrees of
freedom of the classical system. We can take purely classical dynamics to
correspond to identity Lindblad operators which we take to be
$L_{\alpha=0}=\mathbb{I}$. If we wish this classical dynamics to be purely
deterministic, then we can take $\tau_{00}\rightarrow 0$, so that we have a
term in the dynamics given by $\\{H_{C},\varrho\\}$ with the classical
Hamiltonian is defined as

$H_{C}(z):=h^{00}(z)\,\mathbb{I}.$ (8)

Higher order terms can lead to diffusion in phase space. For example, at first
order in $\tau$, Eq. (7) will acquire terms of the form

$\displaystyle{\cal D}(\varrho)=\tau
D^{\alpha\beta}_{p,ij}(z)L_{\alpha}\frac{\partial^{2}\varrho}{\partial
p_{i}\partial p_{j}}L_{\beta}^{\dagger}+\tau
D^{\alpha\beta}_{qp,ij}(z)L_{\alpha}\frac{\partial^{2}\varrho}{\partial
q_{i}\partial p_{j}}L_{\beta}^{\dagger}+\tau
D^{\alpha\beta}_{q,ij}(z)L_{\alpha}\frac{\partial^{2}\varrho}{\partial
q_{i}\partial q_{j}}L_{\beta}^{\dagger}$ (9)

where for ease of presentation we have taken all $\tau_{\alpha\beta}=\tau$.
Such terms lead to diffusion in phase-space, with $D^{\alpha\beta}_{p,ij}(z)$
being the diffusion term usually found in the Fokker-Planck equation, or
arising from Langevin dynamics. When $W^{\alpha\beta}(z|z-\Delta)$ is
deterministic in $q$, but diffusive in $p$, we shall refer to it as dynamics
of the Langevin type. It is possible to also have terms appearing in the
expansion which give rise to friction. We discuss such master equations in
greater detail in Sec. V. If the dynamics is of the Langevin type, then one
can verify that the equation of motion for $\dot{q}$ is unchanged by the
diffusive terms appearing in Eq. (7). Inverting this expression gives
$p(q,\dot{q})$. On the other hand, master equations which contain terms
$D^{\alpha\beta}_{qp,ij}$ and $D^{\alpha\beta}_{q,ij}(z)$ in their expansion,
we shall call q-diffusive and we will discuss such master equations in more
detail elsewhere. Care should be taken in any expansion of
$W^{\alpha\beta}(z|z-\Delta)$ – although the full dynamics is completely
positive, if we truncate our expansion to any finite order, it is not.

###  II.1 Unravelling of the CQ dynamics

In the quantum open system formalism, the unravelling of a master equation
takes the master equation governing the evolution of the density matrix, and
divides it up into trajectories – i.e. probabilities assigned to an ensemble
of histories, where each history is the evolution of a pure quantum states
which either undergoes a continuous evolution or jumps from the one state to
another as a result of a Lindblad operator being applied to it. Due to the
stochastic nature of the unravelling, the initial state is mapped into
different final states by the same evolution. If we record each evolution in
time (or trajectory) and we average over them, we re-obtain the evolution law
for the density matrix, as given via the master equation (in the limit of
infinite trajectories and infinitesimal time steps). The unravelling of the CQ
master equation proceeds in analogy with the quantum one, and in this section
we provide an algorithm to unravel the dynamics given by Eq. (2), where in
addition we consider a fully classical term, given by the Poisson bracket
between a classical Hamiltonian $H_{C}(z)$ and the hybrid state,

$\displaystyle\frac{\partial\,\varrho(z,t)}{\partial t}=$
$\displaystyle-i\left[H(z),\varrho(z,t)\right]+\left\\{{H_{C}(z)},{\varrho(z,t)}\right\\}$
$\displaystyle+\int\mathrm{d}\Delta\,\sum_{\alpha,\beta}W^{\alpha\beta}\left(z|z-\Delta\right)\,L_{\alpha}\,\varrho(z-\Delta,t)\,L^{\dagger}_{\beta}-\frac{1}{2}\sum_{\alpha,\beta}W^{\alpha\beta}\left(z\right)\left\\{{L^{\dagger}_{\beta}L_{\alpha}},{\varrho(z,t)}\right\\}_{+}.$

Specifically, we will provide the tools to generate different classical-
quantum trajectories, and we will show that the average of these trajectories
coincides with the evolution of the CQ state under the above master equation.

In the following, we consider the case in which the matrix of rates
$W^{\alpha\beta}\left(z|z-\Delta\right)$ are diagonal in $\alpha,\beta$. If
this matrix is not diagonal in the same basis at each point, we can introduce
unitary operators $u(z)$ that diagonalize it,

$W^{\prime}_{\gamma}\left(z|z-\Delta\right)=u_{\gamma,\alpha}(z)\,W^{\alpha\beta}\left(z|z-\Delta\right)\,u^{\dagger}_{\beta,\gamma}(z),$
(10)

and use the new Lindblad operators in the unravelling procedure,
$L^{\prime}_{\gamma}(z)=u_{\gamma,\alpha}(z)\,L_{\beta}$. In this way, the
same Lindblad operator acts on the left and the right of the CQ state in Eq.
(2), and we can use the same procedure explained in the following, with the
difference that now the Lindblad operators additionally depends on the phase
space coordinate $z$. The central object of this study is a pure CQ state,
that is,

$\ket{\varrho(z,t)}=\delta(z-\bar{z}(t))\,\ket{\phi(t)},$ (11)

where $\ket{\phi(t)}$ is a quantum state independent of the classical degrees
of freedom. This object represent a classical-quantum system, centred in the
point $\bar{z}(t)$ of the phase space, and described by the quantum state
$\ket{\phi(t)}$.

The first step in the unravelling procedure is the discretisation of the
dynamics of the CQ state; we fix a time interval $\delta t$, and the evolution
of the CQ state is given for multiples of this interval. It is worth noting
that the unravelling evolution is exact at first order in this parameter, and
coincides with the one given by the master equation in the limit of $\delta
t\rightarrow 0$. We now outline the procedure for generating each trajectory,
and we will then show that when we average over trajectories, we obtain the CQ
master equation. For a given point $z$ in phase space, at each time step we
either update the state with a continuous evolution, or with a jump in the
classical and quantum degrees of freedom. When the evolution is continuous, we
apply to the quantum part of the CQ state a unitary evolution given by the
(non-Hermitian) _effective Hamiltonian_

$H_{\text{eff}}(z)=H(z)-\frac{i}{2}\sum_{\alpha}W^{\alpha}\left(z\right)L_{\alpha}^{\dagger}L_{\alpha},$
(12)

which is composed by the Hamiltonian $H(z)$, which appears in the commutator
of Eq. (2), and by the operator appearing in the anti-commutator of the same
equation. The classical part of the CQ state is instead subjected to a
translation in phase space generated by the classical Hamiltonian $H_{C}(z)$,
which appears in the Poisson bracket in the master equation. In order to
evolve the CQ state with a jump, we first need to choose at random a value of
the tuple $\left(\alpha,\Delta\right)$ from the (un-normalised) distribution
$W^{\alpha}\left(z|z-\Delta\right)$. The $\alpha$ parameter specifies which
Lindblad operator $L_{\alpha}$ is performed over the quantum state, while the
$\Delta$ parameter quantifies the shift performed over the classical degrees
of freedom.

The unravelling algorithm is then defined by the following set of
instructions, that specify how to update a pure CQ state at time $t$ under
continuous or jumping evolution,

$\ket{\varrho(z,t+\delta
t)}=\begin{cases}\frac{1}{\sqrt{N_{0}(z,t)}}\,\left(\mathbb{I}-i\,\delta
t\,H_{\text{eff}}(z)\right)\ket{\varrho(z-\delta t\,\Omega\,\nabla
H_{C}(z),t)}&\text{with}\ p_{0}(z,t),\\\
\frac{1}{\sqrt{N_{\alpha,\Delta}(z,t)}}\,L_{\alpha}\ket{\varrho(z-\Delta,t)}&\text{with}\
p_{\alpha,\Delta}(z,t).\end{cases}$ (13)

At each time step, the hybrid pure state can be updated either continuously,
first line in the right hand side of the above equation, or with a jump,
second line. The type of evolution is chosen at random; the probability of
updating the state with a continuous evolution is given by $p_{0}(z,t)$, while
the probability of evolving it with a jump is given by
$p_{\alpha,\Delta}(z,t)$, where $\alpha$ specifies the quantum jump and
$\Delta$ specifies the classical one. In the continuous evolution, the quantum
degrees of freedom are evolved by the unitary operator $U(z)=e^{-i\delta
t\,H_{\text{eff}}(z)}$, that we have expanded at first order in $\delta t$,
while the continuous shift in the classical degrees of freedom is generated by
the operator $e^{\delta t\,\left\\{{H_{C}(z)},{\cdot}\right\\}}$. The gradient
vector $\nabla
H_{C}(z)=\left(\partial_{q_{1}}H_{C}(z),\partial_{p_{1}}H_{C}(z),\ldots\right)^{T}\in\mathbb{R}^{2n}$,
while

$\Omega=\bigoplus_{k=1}^{n}\,\begin{pmatrix}0&1\\\ -1&0\end{pmatrix}\in
M_{2n,2n}(\mathbb{R}),$ (14)

is the symplectic matrix. By iterating the procedure in Eq. (13),we can
produce a single trajectory. Different iterations provide different
trajectories, since at every time step the CQ state is evolved either
continuously or with a jump.

The pure CQ states produced by the unravelling technique are normalized. In
the following, we provide the normalization coefficient for both the
continuous and jumping evolution, and the probability that either evolution
occurs. For simplicity, in the rest of the section we suppress the time $t$
dependence. The notmalization coefficient are given as

$\displaystyle N_{0}(z)$ $\displaystyle=\delta(z-\bar{z})\left(1-\delta
t\sum_{\alpha}W^{\alpha}\left(z\right)\bra{\phi}L_{\alpha}^{\dagger}L_{\alpha}\ket{\phi}\right),$
(15a) $\displaystyle N_{\alpha,\Delta}(z)$
$\displaystyle=\delta(z-(\bar{z}+\Delta))\bra{\phi}L_{\alpha}^{\dagger}L_{\alpha}\ket{\phi},\quad\forall\,\alpha,\Delta,$
(15b)

where we have used the form of the CQ pure state, given in Eq. (11), and we
are only considering terms in the first order of the time step $\delta t$. In
the next section we explicitly derive these coefficients. At a given time
step, the probability of evolving the pure CQ state continuously or with a
jump is given by,

$\displaystyle p_{0}(z)$ $\displaystyle=N_{0}(z),$ (16a) $\displaystyle
p_{\alpha,\Delta}(z)$ $\displaystyle=\delta
t\,W^{\alpha}\left(z|z-\Delta\right)\,N_{\alpha,\Delta}(z),\quad\forall\,\alpha,\Delta.$
(16b)

It is easy to see that the above distribution is a valid probability
distribution over the phase space. Clearly, each probability distribution is
non-negative, and the distribution is normalised,

$\displaystyle\int\mathrm{d}z\,\left(p_{0}(z)+\sum_{\alpha}\int\mathrm{d}\Delta\,p_{\alpha,\Delta}(z)\right)$
$\displaystyle=1-\delta
t\sum_{\alpha}\left(\int\mathrm{d}z\,W^{\alpha}\left(z\right)\delta(z-\bar{z})\right)\bra{\phi}L_{\alpha}^{\dagger}L_{\alpha}\ket{\phi}$
$\displaystyle+\delta
t\sum_{\alpha}\int\mathrm{d}\Delta\,\left(\int\mathrm{d}z\,W^{\alpha}\left(z|z-\Delta\right)\delta(z-(\bar{z}+\Delta))\right)\bra{\phi}L_{\alpha}^{\dagger}L_{\alpha}\ket{\phi}$
$\displaystyle=1-\delta
t\sum_{\alpha}\left(W^{\alpha}\left(\bar{z}\right)-\int\mathrm{d}\Delta\,W^{\alpha}\left(\bar{z}+\Delta|\bar{z}\right)\right)\bra{\phi}L_{\alpha}^{\dagger}L_{\alpha}\ket{\phi}=1,$
(17)

where the last equality follows from Eq. (3).

Before concluding the section, let us comment on an interesting feature of the
unravelling of the hybrid master equation, which is not shared by its quantum
counterpart. If a unique shift $\Delta$ in the classical degrees of freedom is
associated with each Lindblad operator $L_{\alpha}$, then by monitoring the
classical degrees of freedom we have complete knowledge on the jumps that
occur in the quantum part of the hybrid state. Furthermore, the quantum
evolution conditioned on the classical degrees of freedom is unique, since
only a specific sequence of jumps could produce the sequence of shifts in the
classical part. The hybrid state can be said to actually make the transitions
given in Eq. (13). This is in contrast with the unravelling for open quantum
systems, where the lack of a classical system which can be monitored makes it
impossible to single out different trajectories in the Hilbert space.

###  II.2 From unravelling to the CQ master equation

We can now show that the update rule presented in Eq. (13) reproduces, at
first order in $\delta t$, the same dynamics of the hybrid master equation.
First, let us derive the normalization coefficients shown in Eqs. (15). To
compute the normalization of the continuously evolved state, we need to Taylor
expand it and truncate the expansion at first order in $\delta t$. It is
straightforward to show that

$\ket{\varrho(z-\delta t\,\Omega\,\nabla H_{C}(z))}=\left(1+\delta
t\,\left\\{{H_{C}(z)},{\,\cdot\,}\right\\}\right)\ket{\varrho(z)}+O(\delta
t^{2}).$ (18)

The normalization coefficient is then given by

$\displaystyle N_{0}(z,t)$ $\displaystyle=\bra{\varrho(z-\delta
t\,\Omega\,\nabla H_{C}(z),t)}\left(\mathbb{I}+i\,\delta
t\,H_{\text{eff}}(z)^{\dagger}\right)\left(\mathbb{I}-i\,\delta
t\,H_{\text{eff}}(z)\right)\ket{\varrho(z-\delta t\,\Omega\,\nabla
H_{C}(z),t)}$ $\displaystyle=1+\delta
t\left\\{{H_{C}(z)},{\braket{\varrho(z,t)}{\varrho(z,t)}}\right\\}-\delta
t\sum_{\alpha}W^{\alpha}(z)\bra{\varrho(z,t)}L_{\alpha}^{\dagger}L_{\alpha}\ket{\varrho(z,t)}+O(\delta
t^{2})$ $\displaystyle=1-\delta
t\sum_{\alpha}W^{\alpha}\left(z\right)\bra{\varrho(z,t)}L_{\alpha}^{\dagger}L_{\alpha}\ket{\varrho(z,t)}+O(\delta
t^{2}),$

where the second line follows from the Taylor expansion in Eq. (18), from the
form of the effective Hamiltonian $H_{\text{eff}}(z)$ and from Leibniz’s rule
of the Poisson bracket, while the last line follows from the fact that
$\ket{\varrho(z,t)}$ is normalized. The normalization of the CQ pure state
after a jump is simply

$N_{\alpha,\Delta}(z)=\bra{\varrho(z-\Delta,t)}L_{\alpha}^{\dagger}L_{\alpha}\ket{\varrho(z-\Delta,t)}.$
(19)

We can now derive the master equation (2) (equipped with the additional fully-
classical term) from the update rule given in Eq. (13). To so so, we express
the state $\ket{\varrho(z,t+\delta t)}$ as a density operator, and we re-write
it as a mixture over the possible evolutions the system undergoes. Let us
first consider the continuous update for the CQ density operator, which we
refer to as $\rho_{0}(z,t+\delta t)$, where we truncate at first order in
$\delta t$,

$\displaystyle\rho_{0}(z,t+\delta t)=\frac{1}{N_{0}(z,t)}\Bigg{(}$
$\displaystyle\ket{\varrho(z,t)}\bra{\varrho(z,t)}+\delta
t\,\left\\{{H_{C}(z)},{\ket{\varrho(z,t)}\bra{\varrho(z,t)}}\right\\}$
$\displaystyle-i\,\delta
t\left[H(z),\ket{\varrho(z,t)}\bra{\varrho(z,t)}\right]-\frac{\delta
t}{2}\sum_{\alpha}W^{\alpha}\left(z\right)\left\\{{L_{\alpha}^{\dagger}L_{\alpha}},{\ket{\varrho(z,t)}\bra{\varrho(z,t)}}\right\\}_{+}\Bigg{)},$

where we have used the definition of $H_{\text{eff}}(z)$, see Eq. (12), and
$[\cdot,\cdot]$ is the commutator while $\left\\{{\cdot},{\cdot}\right\\}_{+}$
is the anti-commutator. When the evolution is given by an $\alpha$-jump in the
quantum degrees of freedom, and by a $\Delta$-jump in the classical ones, we
find the following density operator

$\rho_{\alpha,\Delta}(z,t+\delta
t)=\frac{1}{N_{\alpha,\Delta}(z,t)}\,L_{\alpha}\ket{\varrho(z-\Delta,t)}\bra{\varrho(z-\Delta,t)}L^{\dagger}_{\alpha}.$
(20)

We can now express the overall evolution of the state by weighting the
different updates with the correct probabilities, given in Eq. (16). Thus, we
find that

$\displaystyle\ket{\varrho(z,t+\delta t)}\bra{\varrho(z,t+\delta t)}$
$\displaystyle=p_{0}(z,t)\,\rho_{0}(z,t+\delta
t)+\sum_{\alpha}\int\mathrm{d}\Delta\,p_{\alpha,\Delta}(z,t)\,\rho_{\alpha,\Delta}(z,t+\delta
t)$ $\displaystyle=\ket{\varrho(z,t)}\bra{\varrho(z,t)}+\delta
t\,\left\\{{H_{C}(z)},{\ket{\varrho(z,t)}\bra{\varrho(z,t)}}\right\\}$
$\displaystyle-i\,\delta
t\left[H(z),\ket{\varrho(z,t)}\bra{\varrho(z,t)}\right]-\frac{\delta
t}{2}\sum_{\alpha}W^{\alpha}\left(z\right)\left\\{{L_{\alpha}^{\dagger}L_{\alpha}},{\ket{\varrho(z,t)}\bra{\varrho(z,t)}}\right\\}_{+}$
$\displaystyle+\delta
t\sum_{\alpha}\int\mathrm{d}\Delta\,W^{\alpha}\left(z|z-\Delta\right)L_{\alpha}\ket{\varrho(z-\Delta,t)}\bra{\varrho(z-\Delta,t)}L^{\dagger}_{\alpha},$
(21)

and if we rearrange the above equation by moving
$\ket{\varrho(z,t)}\bra{\varrho(z,t)}$ from the right to the left hand side,
we divide by $\delta t$, and we send $\delta t\rightarrow 0$, we obtain

$\displaystyle\frac{\partial\ }{\partial
t}\ket{\varrho(z,t)}\bra{\varrho(z,t)}=$
$\displaystyle-i\,\left[H(z),\ket{\varrho(z,t)}\bra{\varrho(z,t)}\right]+\left\\{{H_{C}(z)},{\ket{\varrho(z,t)}\bra{\varrho(z,t)}}\right\\}$
$\displaystyle+\int\mathrm{d}\Delta\,\sum_{\alpha}W^{\alpha}\left(z|z-\Delta\right)L_{\alpha}\ket{\varrho(z-\Delta,t)}\bra{\varrho(z-\Delta,t)}L^{\dagger}_{\alpha}$
$\displaystyle-\frac{1}{2}\sum_{\alpha}W^{\alpha}\left(z\right)\left\\{{L_{\alpha}^{\dagger}L_{\alpha}},{\ket{\varrho(z,t)}\bra{\varrho(z,t)}}\right\\}_{+}.$
(22)

As a last step, we need to move from a CQ pure state $\ket{\varrho(z,t)}$ to a
(possibly mixed) CQ state $\varrho(z,t)$. To do so, we perform an average over
the classical-quantum configurations of the initial state, that is, we average
over a phase space distribution $P_{\text{ps}}$ of initial points
$\bar{z}(t=0)=\bar{z}_{0}$, as well as a distribution $P_{\text{q}}$ over the
initial quantum states $\ket{\phi(t=0)}=\ket{\phi_{0}}$. Thus, by applying
this average over the CQ pure state $\ket{\varrho(z,t)}$, we obtain

$\varrho(z,t)=\mathbb{E}[\ket{\varrho(z,t)}\bra{\varrho(z,t)}]=\int\mathrm{d}\bar{z}_{0}\,\mathrm{d}\phi_{0}\,P_{\text{ps}}(\bar{z}_{0})P_{\text{q}}(\phi_{0})\,\ket{\varrho(z,t)}\bra{\varrho(z,t)}.$
(23)

If we now perform the average in Eq. (II.2), we obtain the CQ master equation.

##  III Main features of hybrid dynamics in qubit toy models

In this section we present a few toy models to illustrate the main features of
the hybrid dynamics generated by the master equation (2). The simplest hybrid
dynamics we can consider is that of a spin half particle in a classical
potential, which in our examples is taken to be linear. For each of the toy
models we derive, both analytically and numerically, the dynamics of the
corresponding hybrid system. To numerically evolve the CQ state we use the
unravelling procedure presented in the previous section. We tailor the
unravelling to the specific master equation used for the toy models, see Eq.
(III). This unravelling procedure is shown in Appendix B.

The main features of the hybrid dynamics, that are highlighted in our toy
models, are i) the presence of stochastic collapses of the quantum degrees of
freedom, ii) a trade-off between the classical diffusion in phase space and
the quantum decoherence, and iii) the fact that energy is not conserved by the
dynamics. More specifically, we find that the interaction between classical
and quantum degrees of freedom generates a sudden collapse in the latter. This
collapse occurs when the CQ state is subjected to a jump, that is, when the
classical degrees of freedom are shifted in phase space and a Lindblad
operator acts over the quantum state. If the Lindblad operators are associated
with unique shifts in the classical degrees of freedom, we find that the
quantum dynamics can be unambiguously recorded by monitoring the classical
one, see for example Fig. 1. Furthermore, we find that the rate of dechoerence
of the quantum degrees of freedom is linked to the rate of diffusion of the
classical ones; the faster the quantum state decoheres, the slower the
classical degrees of freedom spread in phase space. This is not a specific
feature of the toy models we study, but rather a general property of the
hybrid dynamics, as shown in Ref. Oppenheim et al. (2020). Finally, we find
that the average energy of the toy models studied is not conserved by the
hybrid dynamics. This should not come as a surprise, since the Hamiltonian
operator appearing in Eq. (2) is not the generator of the dynamics. A detailed
study of the conservation laws and symmetries present in the hybrid dynamics
is performed in Ref. Oppenheim and Weller-Davies (2020).

Let us introduce the toy models we study in this section, namely, a single
spin half particle in a classical potential. The position and momentum of the
particle are taken to be the classical degrees of freedom
$z=\left(q,p\right)$, while the spin of the particle ($s=\frac{1}{2}$) is the
quantum one. To simplify the notation, we take the transition matrix in the
hybrid master equation (6) to be diagonal in $\alpha,\beta$ (if the transition
matrix is not diagonal, we have shown in the previous section that we can
always diagonalise it).

For a spin-$\frac{1}{2}$ particle, the most general quantum Hamiltonian can be
written as $\boldsymbol{\lambda}\cdot\boldsymbol{\sigma}$ where
$\boldsymbol{\lambda}$ is a vector in $\mathbb{R}^{4}$ and
$\boldsymbol{\sigma}=\left\\{\mathbb{I},\sigma_{x},\sigma_{y},\sigma_{z}\right\\}$
is the vector of Pauli matrices. By having this depend on the classical
degrees of freedom, we cause the classical and quantum system to interact. Let
us take the interaction term of the CQ-Hamiltonian to be

$H_{I}(q,p)=B(q)\left(\omega_{0}\ket{0}\bra{0}+\omega_{1}\ket{1}\bra{1}\right),$
(24)

where the coupling constant $B(q)$ depends on the position of the particle
$q$, and the operator is diagonal in the basis
$\left\\{\ket{0},\ket{1}\right\\}$. By comparing the above interaction
Hamiltonian with the general form of Eq. (5), we can easily see that we are
left some freedom in the choice of the Lindblad operators
$\left\\{L_{\alpha}\right\\}_{\alpha\neq 0}$ and of the functions
$h^{\alpha}(z)$. By exploiting this freedom, we will study two different CQ
models describing a particle in a linear potential. In both cases, we choose a
classical Hamiltonian $H_{C}(q,p)=\frac{p^{2}}{2m}$, that is, the one for a
free particle of mass $m$. One could imagine a purely quantum evolution
determined by a Hamiltonian which is not dependent on classical degrees of
freedom, but in analogy with the gravitational case we will not consider it
here. Finally, we demand the rates $\tau_{\alpha\neq 0}$’s to be all equal to
$\tau>0$, so that no non-trivial Lindblad operator acts more often than the
others.

The CQ master equation (6) for the case we are considering is then given by

$\displaystyle\frac{\partial\,\varrho(z,t)}{\partial t}$
$\displaystyle=-i\left[H_{I}(z),\varrho(z,t)\right]+\frac{1}{\tau_{0}}\left(e^{\tau_{0}\,\left\\{{H_{C}(z)},{\,\cdot\,}\right\\}}-1\right)\varrho(z,t)$
$\displaystyle+\frac{1}{\tau}\sum_{\alpha\neq
0}\left(e^{\tau\,\left\\{{h^{\alpha}(z)},{\,\cdot\,}\right\\}}\,L_{\alpha}\,\varrho(z,t)\,L^{\dagger}_{\alpha}-\frac{1}{2}\left\\{{L^{\dagger}_{\alpha}L_{\alpha}},{\varrho(z,t)}\right\\}_{+}\right).$
(25)

The first term on the right is the evolution of the quantum degrees of freedom
and it depends on the position of the classical particle. The The next term is
the purely classical evolution. If we take the rate $\tau_{0}$ to be a finite
positive value, this classical evolution is stochastic, while if it tends to
$0$ it is determinstic. In the former case, we will see that the equation can
be solved analytically. When $\tau_{0}\rightarrow 0$, the purely classical
evolution is given by the Poisson brackets between the classical Hamiltonian
and the CQ state, $\left\\{{H_{C}(z)},{\,\varrho(z,t)\,}\right\\}$. In this
second case, we analyse numerically the dynamics generated by the master
equation, using the unravelling method shown in App. II.1. The final set of
terms, gives the backreaction of the qubit on the classical degrees of
freedom.

###  III.1 Qubit evolution with diagonal Lindblad operators in a linear
potential

The first toy model we consider is a spin half particle interacting with a
linear potential through Lindblad operators that are diagonal in the
interaction Hamiltonian eigenbasis $\left\\{\ket{0},\ket{1}\right\\}$. We make
the following choice for the operators,

$\displaystyle L_{\alpha=1}$ $\displaystyle=\ket{0}\bra{0},$ (26a)
$\displaystyle L_{\alpha=2}$ $\displaystyle=\ket{1}\bra{1}.$ (26b)

Once the Lindblad operators are defined, they fix uniquely the functions
$h^{\alpha}(q,p)$. Indeed, in order to identify the interaction Hamiltonian
$H_{I}(q,p)$ in Eq. (24) with the expression in Eq. (5), we need to define
these functions as

$\displaystyle h^{\alpha=1}(q,p)$ $\displaystyle=\omega_{0}\,B(q),$ (27a)
$\displaystyle h^{\alpha=2}(q,p)$ $\displaystyle=\omega_{1}\,B(q).$ (27b)

It is worth noting that, for a given interaction Hamiltonian, the above choice
of Lindblad operators and $h^{\alpha}(q,p)$ functions is not unique. In the
next section, for instance, we consider non-diagonal operators that
nonetheless return the same interaction Hamiltonian, but generate a completely
different dynamics. Furthermore, one could chose the same Lindblad operators
as we did above, but with a different normalization. To re-obtain the correct
Hamiltonian $H_{I}(q,p)$ we then need to scale the functions $h^{\alpha}(q,p)$
consistently. Notice that the normalization of the Lindblad operators
influences the rate at which the operator is applied, while the function
$h^{\alpha}(q,p)$ regulates the amplitude of the shifts in the classical
degrees of freedom. This freedom in the Lindblad operators and the functions
$h^{\alpha}(q,p)$ is linked to the fact that the interaction Hamiltonian is
not the generator of the hybrid dynamics of Eq. (2), and this is highlighted
by the fact that this operator is in fact not conserved by the evolution.

In this example, we consider a potential linearly dependent in position, so
that $B(q)=qB$, where $B$ has a constant value. We can express the CQ state
$\varrho(q,p,t)$ as

$\varrho(q,p,t)=\begin{pmatrix}u_{0}(q,p,t)&c(q,p,t)\\\
c^{\star}(q,p,t)&u_{1}(q,p,t)\end{pmatrix},$ (28)

where $u_{i}(q,p,t)$ is the population of the quantum state $\ket{i}$ (for
$i=0,1$), and $c(q,p,t)$ is the coherence.

We now re-write the CQ master equation given in Eq. (III), in terms of the
populations and the coherence of the state $\varrho(q,p,t)$. In the equations
below, we take $\tau_{0}\rightarrow 0$, so as to obtain a system of non-local
differential equations. When we solve the equations analytically, however, we
will require $\tau_{0}$ to be a positive (albeit small) constant. The dynamics
of the populations and coherence is thus given by

$\displaystyle\frac{\partial\,u_{i}(q,p,t)}{\partial t}$
$\displaystyle=-\frac{p}{m}\frac{\partial\,u_{i}(q,p,t)}{\partial
q}+\frac{1}{\tau}\big{(}u_{i}(q,p+\omega_{i}B\tau,t)-u_{i}(q,p,t)\big{)},\quad
i\in\left\\{0,1\right\\},$ (29a)
$\displaystyle\frac{\partial\,c(q,p,t)}{\partial t}$
$\displaystyle=-iqB\left(\omega_{0}-\omega_{1}\right)c(q,p,t)-\frac{p}{m}\frac{\partial\,c(q,p,t)}{\partial
q}-\frac{1}{\tau}\,c(q,p,t).$ (29b)

It is worth noting that the equation for the coherence can be solve
analytically, and the solution has the form

$c(q,p,t)=\tilde{c}(q-\frac{p}{m}t)\,e^{-iqB(\omega_{0}-\omega_{1})t-\frac{t}{\tau}},$
(30)

for any function $\tilde{c}(q-\frac{p}{m}t)$ in $C^{1}$, the set of all
continuously differentiable functions. We thus see that the coherence term
decays exponentially fast. We will see that this is a result of the quantum
system collapsing to the $0$ or $1$ state while making a momentum jump.

Figure 1: The evolution of the population $u_{0}(q,p,t)$ and $u_{1}(q,p,t)$ in
phase space, under the CQ qubit dynamics with diagonal Lindblad operators. The
initial CQ state is $\delta(q)\delta(p)\ket{+}$, where
$\ket{+}=\frac{1}{\sqrt{2}}\left(\ket{0}+\ket{1}\right)$. The interaction with
the field makes the quantum state stochastically collapse into either
$\ket{0}$ or $\ket{1}$ at a rate given by $\frac{\omega}{\tau}$. The same
collapse shifts the classical momentum of the particle by $\pm B\tau$,
depending on the quantum state it projects onto. By monitoring the classical
degrees of freedom, we can establish the state of the quantum system. Indeed,
a particle in the state $\ket{0}$ has negative momentum, whereas a particle in
$\ket{1}$ has positive momentum. The evolution of the populations is shown for
$t=0.02\,\mathrm{s}$, $0.06\,\mathrm{s}$ and $0.1\,\mathrm{s}$. The state is
evolved for $t=0.1\,\mathrm{s}$, with time steps of $\delta
t=10^{-4}\,\mathrm{s}$, and the jump rate is $\tau=10^{-2}\,\mathrm{s}$. The
constant $B=1\,\mathrm{J}\cdot\mathrm{s}\cdot\mathrm{m}^{-1}$, the mass
$m=1\,\mathrm{kg}$, and the frequencies are
$\omega_{0}=-\omega_{1}=\omega=1\,\mathrm{s}^{-1}$.

####  III.1.1 Analytical and numerical evolution of populations

In order to study the evolution of the populations analytically, we can re-
express Eq. (29a) as a stochastic equation, see App. A. To do so, we do not
send to $0$ the rate $\tau_{0}$ in Eq. (III), and we obtain the following
evolution for the populations,

$\frac{\partial\,u_{i}(q,p,t)}{\partial
t}=\frac{1}{\tau_{0}}\big{(}u_{i}(q-\frac{p}{m}\tau_{0},p,t)-u_{i}(q,p,t)\big{)}+\frac{1}{\tau}\big{(}u_{i}(q,p+\omega_{i}B\tau,t)-u_{i}(q,p,t)\big{)},\quad
i\in\left\\{0,1\right\\},$ (31)

The solution of this equation is given in the appendix for general
Hamiltonians, see Eq. (91). For the specific case we have are considering, the
solution is

$u_{i}(q,p,t)=\sum_{k,n=0}^{\infty}P_{0}(k)\,P_{1}(n)\,u_{i}(n,k)$ (32)

where $k$ is the number of position shifts, $n$ is the number of momentum
shifts, and both $P_{0}(k)$ and $P_{1}(n)$ are Poisson distributions with mean
value $\lambda_{0}=\frac{t}{\tau_{0}}$ and $\lambda_{1}=\frac{t}{\tau}$,
respectively. The function $u_{i}(n,k)$ is given by

$u_{i}(k,n)=\frac{1}{\binom{n+k}{k}}\sum_{i}\pi_{i}\left(\underbrace{e^{-\frac{p}{m}\tau_{0}\partial_{q}}\ldots
e^{-\frac{p}{m}\tau_{0}\partial_{q}}}_{k}\,\underbrace{e^{B\omega_{i}\tau\partial_{p}}\ldots
e^{B\omega_{i}\tau\partial_{p}}}_{n}\right)u_{i}(q,p,0),$ (33)

where each $\pi_{i}$ is a permutation that creates a different combination of
the $n+k$ shifts.

The model under consideration can additionally be solved numerically, using
the unravelling method presented in the previous section. In particular, we
can use Eq. (III) to better understand the evolution of the CQ state. In this
equation, both classical and quantum degrees of freedom evolve either
continuously, or with a sudden jump. The latter evolution is the most
interesting, since during a jump one of the Lindblad operators is applied to
the quantum state, and this state is projected either in $\ket{0}$ or in
$\ket{1}$, so that the coherence of the state is destroyed, and the classical
momentum receives a kick proportional to $B\tau$, whose direction depends on
the quantum state. From this point of view, this hybrid model mimics some of
the features of a standard Stern-Gerlach experiment, namely the fact that when
the magnetic field measures the quantum state (and collapses it into $\ket{0}$
or $\ket{1}$), it kicks the particle (that is, increase the classical
momentum) either upward or downward, depending on the state the spin is
collapsed into. In Fig. 1 we provide the numerical solution for the evolution
of the populations, obtained using the unravelling technique of Sec. (II.1).

Figure 2: The mean value and standard deviation of the population
$u_{0}(q,p,t)$ in phase space. Left. The mean value $\langle q\rangle$ and
standard deviation $\sigma_{q}$ of the marginal distribution $u_{0}(q,t)$
(where we traced out over the momentum degree of freedom). The solid lines are
the values obtained from the numerical simulation, whereas the dashed line are
obtained from the analytical solution. Centre. The mean value $\langle
p\rangle$ and standard deviation $\sigma_{p}$ of the marginal distribution
$u_{0}(p,t)$. Right. The coherence $\left|c\right|^{2}$ of the quantum state
located in $q=0\,\mathrm{m}$ and
$p=0\,\mathrm{kg}\cdot\mathrm{m}\cdot\mathrm{s}^{-1}$ , as a function of time.
It is worth noting that the numeric solution is obtained with the unravelling
method, in which the position is updated during the continuous evolution step,
see Eq. (115). Instead, the analytic solution is obtained for the case in
which the position evolves through jumps, see Eq. (31). In order to match the
two evolution, we have to ask the stochastic rate $\tau_{0}$ to be of the
order of the infinitesimal time steps $\delta t$ used in the simulation. The
mean value and variance of position and momentum are therefore computed using
Eqs. (34) and (35) for a value of $\tau_{0}=\delta t$.

####  III.1.2 Diffusion in phase-space and trade-off with decoherence

We can also study the diffusion in phase-space of the populations of the CQ
state, both numerically and analytically. Let us take the initial state of the
CQ system to be at the origin of phase-space (position $q(t=0)=0$ and momentum
$p(t=0)=0$), and in the quantum state
$\ket{+}=\frac{1}{\sqrt{2}}\left(\ket{0}+\ket{1}\right)$. In the unravelling
picture when a quantum jump occurs, the state jumps to either $\ket{0}$ or
$\ket{1}$ due to the application of the Lindblad operators of Eq. (26) to the
state. After the jump, the quantum state cannot change anymore, but additional
jumps do increase the momentum of the particle. Since the initial condition is
symmetric in $\ket{0}$ and $\ket{1}$, as the evolution of the populations is,
we can focus on the population $u_{0}(q,p,t)$. The kind of decoherence
affecting the quantum degrees of freedom in this model can be understood in
terms of a leakage of information about the quantum degrees of freedom into
the classical degrees of freedom. Indeed, each Lindblad operator $L_{\alpha}$
is here associated to a distinct classical jump in phase-space, which allows
an observer monitoring the classical degrees of freedom to know exactly in
which quantum state the system is in and when the transition occurred. Thus,
the possibility of monitoring the quantum system using its position in phase
space removes any coherence in the basis $\left\\{\ket{0},\ket{1}\right\\}$.
Furthermore, by conditioning over the classical degrees of freedom, the
quantum state remains pure and does not become mixed.

This is different to the standard decoherence found in the usual Lindblad
equation. There, the decomposition of the dynamics into Lindblad operators is
not unique, and there is no physical meaning to the jumps – rather, the
density matrix of the system evolves continuously in time with the off-
diagonal matrix elements slowly decaying with time. In contrast, here, the
quantum jump is accompanied by a discontinuous jump in $p$. The quantum jump
thus has physical meaning, since an observer who is monitoring the classical
degree of freedom can verify the quantum jump. If they observe a sudden
increase (decrease) in momentum, they will expect that the quantum state is
now in the $\ket{1}$ ($\ket{0}$) state, and can verify this by measuring the
quantum state. The crucial ingredient here is that each classical jump
corresponds to a single Lindblad operator being applied to the quantum state.

We will now see that the quantum system’s coherence time is related to the
diffusion in phase space. In order to obtain the diffusion constant from the
numeric results, we compute the marginal distribution for position and
momentum as a function of time, and we fit them with a Gaussian distribution
to extrapolate mean value and variance. These quantities can also be computed
using the solution in Eq. (32), as we show in App. A.2.1. From the analytic
solution we find that the mean value of position and momentum is given by

$\displaystyle\langle q\rangle(t)$
$\displaystyle=\frac{1}{2}\frac{B\omega_{0}}{m}t^{2}+\langle q\rangle(0),$
(34a) $\displaystyle\langle p\rangle(t)$ $\displaystyle=B\omega_{0}t+\langle
p\rangle(0),$ (34b)

which is the expected solution for the position and momentum of a particle in
a linear potential given by $V(q)=qB\omega_{0}$. We notice that the
expectation values of position and momentum are independent of the parameter
$\tau_{0}$, which can therefore be sent to $0$. The variance of position and
momentum can be obtained too, as we show in the Appendix A, and are given by

$\displaystyle\sigma_{q}^{2}(t)$
$\displaystyle=\frac{1}{3}\left(\frac{B\omega_{0}t}{m}\right)^{2}\left((\tau_{0}+\tau)t+5\tau_{0}\tau\right)\xrightarrow{\tau_{0}\rightarrow
0}\frac{1}{3}\left(\frac{B\omega_{0}t}{m}\right)^{2}\tau t,$ (35a)
$\displaystyle\sigma_{p}^{2}(t)$ $\displaystyle=(B\omega_{0})^{2}\tau t.$
(35b)

The variation in momentum is due to the fact that the number of momentum jumps
will be normally distributed and as with a random walk, increase like
$\sqrt{t}$, while the variation in position is a consequence of the momentum
having a variation. We can additionally link the variance in momentum to the
diffusion term appearing in the expansion of the hybrid master equation, see
Eq. (9). To do so, we need to expand the exponential operator in Eq. (III) to
second order, obtaining the following diffusion term,

$\mathcal{D}\left(u_{0}(q,p,t)\right)=\frac{1}{\tau}\left\\{{\tau\,h^{\alpha=1}(q,p)},{\left\\{{\tau\,h^{\alpha=1}(q,p)},{u_{0}(q,p,t)}\right\\}}\right\\}=(B\omega_{0})^{2}\tau\,\frac{\partial^{2}}{\partial\,p^{2}}u_{0}(q,p,t).$
(36)

It is then easy to see that the variance in momentum arises from the diffusion
coefficient $D^{\alpha=1}_{p}(q,p)=(B\omega_{0})^{2}\tau$. This is consistent
with the fact that the expansion of the master equation at first order in
$\tau$ has the form of a Fokker-Planck equation with the above diffusion term.

In Fig. 2 we compare the numerical results with the ones obtained
analytically, and we show that indeed the numeric simulation accurately
describes the evolution of the CQ state. It is worth noting that, while the
decoherence of the quantum system is inversely proportional to $\tau$ as seen
from Eq. (30), the diffusion in phase space is directly proportional to it.
Thus, we can have two opposite situations; for $\tau\ll 1$, the state of the
quantum system quickly collapses and does not diffuses much in phase space,
following an almost Liouvillian deterministic evolution. When $\tau\gg 1$,
instead, the quantum system slowly decoheres, but the classical degrees of
freedom quickly diffuse in phase space. We thus have a trade-off between the
decoherence rate and the diffusion rate. This turns out to be a standard
feature of hybrid dynamics and can be understood in terms of the moments of
the Kramers-Moyal expansion of the CQ master equation, as proven in Oppenheim
et al. (2020).

####  III.1.3 Energy conservation

We now turn to the question of energy conservation in this model of hybrid
dynamics. We will see that there is a violation of energy conservation at a
rate proportional to $\tau$, although it can be made arbitrarily small. The
average energy of the system can be computed using the total Hamiltonian
$H(q,p)=H_{C}(q,p)+H_{I}(q,p)$, and the solution obtained in Eq. (32) for the
populations. Indeed, the average energy of the CQ system is given by

$\langle
H\rangle(t)=\int\mathrm{d}q\mathrm{d}p\,\mathrm{Tr}\left[{H(q,p)\varrho(q,p,t)}\right]=\sum_{i=0,1}\int\mathrm{d}q\mathrm{d}p\left(\frac{p^{2}}{2m}+q\omega_{i}B\right)u_{i}(q,p,t).$
(37)

Let us consider, for simplicity, that the initial state of the CQ system is
$\varrho(q,p,0)=\delta(q)\delta(p)\ket{0}\bra{0}$, so that we only need to
consider the energy contribution of the population $u_{0}$. Furthermore, given
the specific initial state of the system, we can express the ensemble
$u_{0}(n,k)$, shown in Eq. (33), as
$u_{0}(n,k)=\int\mathrm{d}\bar{q}\mathrm{d}\bar{p}\
P(\bar{q},\bar{p}\,|\,n,k)\ u_{0}(\bar{q},\bar{p})$, where
$u_{0}(\bar{q},\bar{p})=\delta(q-\bar{q})\delta(p-\bar{p})$. The evolution of
the population $u_{0}$ can therefore be expressed as

$u_{0}(q,p,t)=\int\mathrm{d}\bar{q}\mathrm{d}\bar{p}\left(\sum_{k,n=0}^{\infty}P_{0}(k)\,P_{1}(n)\,P(\bar{q},\bar{p}\,|\,n,k)\right)u_{0}(\bar{q},\bar{p})=\int\mathrm{d}\bar{q}\mathrm{d}\bar{p}\,P(\bar{q},\bar{p})\,u_{0}(\bar{q},\bar{p}),$
(38)

and by replacing it in Eq. (37) we find that the average energy is

$\langle
H\rangle(t)=\int\mathrm{d}\bar{q}\mathrm{d}\bar{p}\,P(\bar{q},\bar{p})\left(\frac{\bar{p}^{2}}{2m}+\bar{q}\omega_{0}B\right)=\int\mathrm{d}\bar{p}\,P(\bar{p})\,\frac{\bar{p}^{2}}{2m}+\int\mathrm{d}\bar{q}\,P(\bar{q})\,\bar{q}\omega_{0}B,$
(39)

where the mean value and variance of the marginal probabilities $P(\bar{q})$
and $P(\bar{p})$ have been computed in the previous section, in Eqs. (34) and
(35) respectively.

Using these results, we can show that the energy associated with the classical
Hamiltonian is given by $\langle
H_{C}\rangle(t)=\frac{(B\omega_{0}\tau)^{2}}{2m}\left(\frac{t}{\tau}+\left(\frac{t}{\tau}\right)^{2}\right)$,
while the quantum coupling contributes to the energy as $\langle
H_{I}\rangle(t)=-\frac{(B\omega_{0})^{2}}{2m}t^{2}$ (notice that in the
previous section we have considered the absolute value of the position for
convenience). As a result, we find that the average energy of the system is
not conserved, and instead grows linearly in time from the initial value of
$\langle H\rangle(0)=0$,

$\langle H\rangle(t)=\frac{(B\omega_{0})^{2}\tau}{2m}t,$ (40)

see Fig. 3. In the limit when $\tau$ tends to $0$, one recovers the classical
evolution given by the Liouville equation, and energy is conserved, but
otherwise the energy increases at a rate proportional to the jump distance
$\tau$.

As we stressed at the beginning of the section, the fact that the average
energy is not conserved does not come as a surprise. Indeed, we are taking the
energy to be given by the operator $H(z)$ even though it is not the generator
of the symmetry of time-translation. In this model however, the failure of the
system to conserve energy is not due to changes in the quantum system – the
interaction of the qubit with the classical system merely causes it to
decohere in a basis which commutes with the total CQ Hamiltonian. Instead, the
term appearing in Eq. (40) is due to dispersion in the momentum of the
classical system. This can be prevented by adding in a friction term, as
discussed in Sec. V.

Figure 3: The different contributions to the total energy of the CQ dynamics,
as a function of time. The solid lines are the numerical results, while the
dashed ones are the analytical ones. The energy contribution given by the
kinetic energy, $E_{C}=\langle H_{C}\rangle$, grows faster than the one given
by the potential energy, $E_{I}=\langle H_{I}\rangle$. As a result, the total
energy of the system, $E_{\mathrm{tot}}=\langle H\rangle$, increases linearly
as a function of time.

###  III.2 Qubit evolution with non-diagonal Lindblad operators in linear
potential

We now consider the same hybrid system, with the difference that now the
Lindblad operators are non-diagonal in the interaction Hamiltonian eigenbasis.
We will see that the dynamics of this CQ system is qualitatively and
quantitatively different from the one in the previous section. The form of the
Lindblad operators is

$\displaystyle L_{\alpha=1}$ $\displaystyle=\ket{1}\bra{0},$ (41a)
$\displaystyle L_{\alpha=2}$ $\displaystyle=\ket{0}\bra{1}.$ (41b)

In order to obtain the interaction Hamiltonian $H_{I}(q,p)$ of Eq. (24) from
these operators, the functions $h^{\alpha}(q,p)$ are

$\displaystyle h^{\alpha=1}(q,p)$ $\displaystyle=\omega_{0}\,qB,$ (42a)
$\displaystyle h^{\alpha=2}(q,p)$ $\displaystyle=\omega_{1}\,qB,$ (42b)

where again the particle is in a linear potential, i.e., $B(q)=qB$, with
$B=\text{const}$. If we express the CQ state $\varrho(q,p,t)$ as in Eq. (28),
and we take the rate $\tau_{0}\rightarrow 0$, then the master equation given
in Eq. (III) can be expressed as a system of non-local differential equations.
In particular, the evolution of the populations is given by of the following
differential equations,

$\frac{\partial\,u_{i}(q,p,t)}{\partial
t}=-\frac{p}{m}\frac{\partial\,u_{i}(q,p,t)}{\partial
q}+\frac{1}{\tau}\big{(}u_{i\oplus
1}(q,p+\omega_{i}B\tau,t)-u_{i}(q,p,t)\big{)},\quad i\in\left\\{0,1\right\\},$
(43)

where $\oplus$ is addition modulo $2$. Notice that the above equations have
the same form as the ones obtained in the previous section, see Eq. (29a), but
now each population depends non-locally on the other. The evolution for the
coherence is unchanged from the one in Eq. (29b).

Figure 4: The evolution of the population $u_{0}$ (blue) and $u_{1}$ (green)
under the CQ dynamics with non-diagonal Lindblad operators. The evolution of
the populations is shown for $t=0\,\mathrm{s}$, $0.04\,\mathrm{s}$ and
$0.1\,\mathrm{s}$. The initial CQ state is $\delta(q)\delta(p)\ket{+}$, where
$\ket{+}=\frac{1}{\sqrt{2}}\left(\ket{0}+\ket{1}\right)$ is a coherent
superposition of the eigenstates of the interaction Hamiltonian $H_{I}$. This
interaction makes the quantum degrees of freedom stochasitcally collapse from
$\ket{0}$ to $\ket{1}$ and vice versa, at a rate given by $\tau$. When the
quantum state is mapped as $\ket{0}\rightarrow\ket{1}$, the classical momentum
is increased by $B\tau$, while the transformation $\ket{1}\rightarrow\ket{0}$
is accompanied by a change in momentum equal to $-B\tau$. The simulation
depicted here uses time steps of $\delta t=10^{-4}\,\mathrm{s}$, and a jump
rate $\tau=10^{-2}\,\mathrm{s}$. The interaction constant is
$B=1\,\mathrm{J}\cdot\mathrm{s}\cdot\mathrm{m}^{-1}$, the mass is
$m=1\,\mathrm{kg}$, and the frequencies are
$\omega_{0}=-\omega_{1}=\omega=1\,\mathrm{s}^{-1}$.

####  III.2.1 Analytical and numerical evolution of populations

To solve the model analytically we can express Eq. (43), which describes the
evolution of the populations, as a stochastic equation, see App. A. To do so,
we fix the rate $\tau_{0}$ in Eq. (III) to be a (small) positive constant, so
that the equation for the populations can be expressed as

$\frac{\partial\,u_{i}(q,p,t)}{\partial
t}=\frac{1}{\tau_{0}}\big{(}u_{i}(q-\frac{p}{m}\tau_{0},p,t)-u_{i}(q,p,t)\big{)}+\frac{1}{\tau}\big{(}u_{i\oplus
1}(q,p+\omega_{i}B\tau,t)-u_{i}(q,p,t)\big{)},\quad i\in\left\\{0,1\right\\},$
(44)

If $\tau_{0}\rightarrow 0$, we come back to the original equation for the
populations shown in Eq. (43). In order to solve this equation we can make use
of the tools described in appendix, see Sec. A.3, and the solution we obtain
is,

$u_{i}(q,p,t)=\sum_{\ell,j=0}^{\infty}P_{0}(\ell)\,P_{1}(2j)\,u_{i}(\ell,j)+\sum_{\ell^{\prime},j^{\prime}=0}^{\infty}P_{0}(\ell^{\prime})\,P_{1}(2j^{\prime}+1)\,u_{i\oplus
1}(\ell^{\prime},j^{\prime}),\quad i\in\left\\{0,1\right\\},$ (45)

where $P_{0}$ and $P_{1}$ are Poisson distributions with mean value
$\frac{t}{\tau_{0}}$ and $\frac{t}{\tau}$, respectively. For the specific case
under consideration, we can express the state in terms of shift operators

$\displaystyle u_{i}(\ell,j)$
$\displaystyle=\frac{1}{\binom{\ell+2j}{\ell}}\sum_{n}\pi_{n}\left(\left(e^{-\frac{p}{m}\tau_{0}\partial_{q}}\right)^{\ell}\left(e^{B\omega_{i}\tau\partial_{p}}\,e^{B\omega_{i\oplus
1}\tau\partial_{p}}\right)^{j}\right)u_{i}(q,p,t=0),$ (46a) $\displaystyle
u_{i\oplus 1}(\ell,j)$
$\displaystyle=\frac{1}{\binom{\ell+2j+1}{\ell}}\sum_{n}\pi_{n}\left(\left(e^{-\frac{p}{m}\tau_{0}\partial_{q}}\right)^{\ell}\left(e^{B\omega_{i}\tau\partial_{p}}\,e^{B\omega_{i\oplus
1}\tau\partial_{p}}\right)^{j}e^{B\omega_{i}\tau\partial_{p}}\right)u_{i\oplus
1}(q,p,t=0),$ (46b)

where the permutations $\pi_{n}$ create different combinations of the shift
operators for position and momentum, while preserving the relative order of
the shift operators for the momentum.

To better understand the above solution, let us consider the case in which the
initial CQ state is $\delta(q)\delta(p)\ket{0}\bra{0}$, so that at $t=0$ only
the level $\ket{0}$ is populated. Additionally, we fix
$\omega_{0}=-\omega_{1}=\omega>0$, so that opposite shifts in momentum cancel
each others, see Fig. 4. At time $t$ we find that both quantum levels are
populated, and in partiular we have

$\displaystyle u_{0}(q,p,t)$
$\displaystyle=\sum_{\ell,j=0}^{\infty}\frac{1}{\binom{\ell+2j}{\ell}}P_{0}(\ell)\,P_{1}(2j)\sum_{n}\pi_{n}\left(\left(e^{-\frac{p}{m}\tau_{0}\partial_{q}}\right)^{\ell}\left(e^{B\omega\tau\partial_{p}}\,e^{-B\omega\tau\partial_{p}}\right)^{j}\right)\delta(q)\delta(p),$
(47a) $\displaystyle u_{1}(q,p,t)$
$\displaystyle=\sum_{\ell,j=0}^{\infty}\frac{1}{\binom{\ell+2j+1}{\ell}}P_{0}(\ell)\,P_{1}(2j+1)\sum_{n}\pi_{n}\left(\left(e^{-\frac{p}{m}\tau_{0}\partial_{q}}\right)^{\ell}\left(e^{B\omega\tau\partial_{p}}\,e^{-B\omega\tau\partial_{p}}\right)^{j}e^{B\omega\tau\partial_{p}}\right)\delta(q)\delta(p).$
(47b)

By counting the number of shift operators in momentum, we can see that the
population $u_{0}$ spreads in position $q$ while keeping the momentum $p$
fixed at $0$. On the other hand, the population $u_{1}$ has an odd number of
momentum shifts, and therefore its momentum is fixed at $-B\omega\tau$.

We can also use the unravelling technique to obtain a numerical solution for
the CQ dynamics of Eq. (43). In this equation, we see that the CQ state can
evolve either continuously, or through a jump in both the quantum and
classical degrees of freedom. When a jump occur, the quantum state is
projected in either $\ket{0}$ or $\ket{1}$, and this change is accompanied by
a positive or negative shift in momentum, respectively. Furthermore, due to
the form of the Lindblad operators, we have that a CQ state whose quantum
degree of freedom is described by $\ket{0}$ can only jump when it is hit by
the operator $L_{0}$, and vice versa for $\ket{1}$. In our simulation, we
require these jumps to be associated with an opposite change in momentum. In
this case, an initial CQ state with well-defined momentum (for instance the
one used in Fig. 4) can spread to at most 3 different values of momentum.

####  III.2.2 Diffusion in phase-space

We can additionally study the spreading, due to the interaction between
classical and quantum degrees of freedom, of the position of the particle as
time passes by. The numerical simulation shows that at later times, the
populations $u_{0}(q,p,t)$ and $u_{1}(q,p,t)$ divide into two Gaussian
distributions, one with zero momentum, the other with non-zero momentum, see
Fig. 4. Here, we consider the same scenario as in the previous section, where
the initial CQ state is $\varrho(q,p)=\delta(q)\delta(p)\ket{0}\bra{0}$. As we
noticed before, this state evolves into two different ensembles, one
associated with the quantum level $\ket{0}$, and the other with the quantum
level $\ket{1}$. In particular, the population $u_{0}$ spreads in position
while keeping the momentum fixed at $p=0$, while the population $u_{1}$
spreads in position with a fixed momentum $p=-\omega B\tau$. In Fig. 5, we
show the average position and standard deviation of the two populations as a
function of time. We obtain these values by fitting the numerical data, and by
the following analytical considerations.

Figure 5: The mean value $\langle q\rangle$ (or better, its absolute value)
and standard deviation $\sigma_{q}$ of populations $u_{0}$ and $u_{1}$. For
the numerical simulation we choose the frequencies
$\omega_{0}=-\omega_{1}=\omega=1\,\mathrm{s}^{-1}$, the interaction constant
$B=1\,\mathrm{J}\cdot\mathrm{s}\cdot\mathrm{m}^{-1}$, the mass
$m=1\,\mathrm{kg}$, and a jump rate $\tau=10^{-2}\,\mathrm{s}$. The time step
used in the simulation is $\delta t=10^{-4}\,\mathrm{s}$. The solid lines are
obtained by fitting the numerical data with a normal distribution, while the
dashed lines are obtained using analytical considerations discussed in the
main text. It is interesting to notice that, while $\sigma_{q}$ is roughly the
same for the two populations, the mean position is shifted. This is because
the distribution of $u_{1}$ is the first to be populated, since only one jump
in momentum is needed when starting from the initial quantum state
$\ket{0}\bra{0}$. On the other side, the distribution of $u_{0}$ is populated
(for values of $q\neq 0$) only after two jumps in momentum are preformed.
Since the rate of jump is given by $\tau$, we find that the mean position of
the two populations is actually shifted (in time) by this parameter.

Let us consider the evolution in phase space of the initial state under
consideration. At each time step, this state has a probability
$P_{\text{jump}}=\frac{\delta t}{\tau}$ of jumping in momentum, and a
probability $1-P_{\text{jump}}$ of jumping in position. As we noticed before,
given the initial state under consideration, the momentum of the CQ state at
time $t$ can be either $0$ or $p_{\min}=-B\omega\tau$. When the momentum is
non-zero, a position jump modifies the position of the state by $\Delta
q=p_{\min}\delta t=-B\omega\tau\delta t$. We can then compute the average
change in position between two jumps, when the momentum is $p_{\min}$. The
(normalised) probability that $n$ time steps separate two jumps is given by

$\text{prob}\left(\\#\,\text{steps}=n\right)=P_{\text{jump}}(1-P_{\text{jump}})^{n}.$
(48)

Using the above distribution, we can compute the average distance travelled by
the state in between two momentum jumps, which is given by

$\Delta Q=\sum_{n=0}^{\infty}n\,\Delta
q\,\text{prob}\left(\\#\,\text{steps}=n\right)=\frac{1-P_{\text{jump}}}{P_{\text{jump}}}\Delta
q=-B\omega\tau\left(\tau-\delta t\right).$ (49)

After $n$ time steps, the average number of jumps is
$N_{\text{jump}}=nP_{\text{jump}}$. Furthermore, the position of the state can
only change after an odd number of jumps, since only in that case the momentum
is non-zero. As a result, we find that the average position of the state at
time $t$ is

$\langle q\rangle(t)=\frac{1}{2}N_{\text{jump}}\Delta Q+\langle
q\rangle(0)=\frac{1}{2}\left(\tau-\delta t\right)t+\langle q\rangle(0),$ (50)

where we define $t=n\,\delta t$. It is worth noting that the above result only
applies to the average position of the population $u_{1}$. In this case, since
the initial quantum state is $\ket{0}$, a single jump in momentum is required
for the CQ state to populate $\ket{1}$. However, the CQ state starts to
populate $\ket{0}$ (with a position $q\neq 0$), only after two jumps in
momentum are performed. Since the rate of jumping is given by the parameter
$\tau$, we have that the average position of the population $u_{0}$ is roughly
delayed by this amount of time, and it is therefore equal to $\langle
q\rangle(t-\tau)$, see Fig. 5.

To compute the variance of the position for the populations, we notice that
the number of jumps $N_{\text{jump}}$ is distributed according to a binomial
distribution, and its variance is given by
$\sigma_{\text{jump}}^{2}=nP_{\text{jump}}\left(1-P_{\text{jump}}\right)$.
Then, keeping the other terms constant in Eq. (50), we have

$\sigma_{q}(t)=\frac{1}{2}\sigma_{\text{jump}}\left|\Delta
Q\right|=\frac{B\omega}{2}(\tau-\delta t)^{\frac{3}{2}}\sqrt{t}.$ (51)

As we noticed before, the CQ state starts in position $q=0$ with momentum
$p=0$, and roughly starts moving only after $t=\tau$, that is, the time
interval before a momentum jumps occurs (on average). As a result, the
standard deviation of the populations needs to be modified so as to account
for this delay, and we find that for both $u_{0}$ and $u_{1}$ the standard
deviation is given by $\sigma_{q}(t-\tau)$, at least for $t\gg\tau$. As a
final remark we notice that also in this model, like in the previous one, the
average energy, as computed by averaging the operator $H(z)=H_{C}(z)+H_{I}(z)$
over the classical and quantum degrees of freedom, is not conserved. Perhaps
surprisingly, though, the energy increases from the initial time until it
saturates at a positive value, see Fig. 6.

Figure 6: The energy of the CQ state as a function of time, when the initial
state is $\delta(q)\delta(p)\ket{0}\bra{0}$. The different contributions to
the energy are shown, such as the kinetic energy $E_{C}=\langle H_{C}\rangle$
and the potential energies $E_{I}^{(i)}=\langle H_{I}^{(i)}\rangle$. As in the
first model in the previous section, the total energy of the system is not
conserved. However, in this case the total energy $E_{\text{tot}}=\langle
H_{C}+H_{I}\rangle$ increases during the initial part of the evolution, until
it reaches the maximum value $\frac{1}{2}\tau^{2}$.

##  IV Quantum harmonic oscillator in a classical potential

In this section, we will consider the dynamics of a quantum harmonic
oscillator whose internal degrees of freedom are coupled to its classical
position and momentum. The model we focus on here demonstrates some basic
features of the classical-quantum dynamics, though it describes a somewhat
artificial situation where raising the energy level of the harmonic oscillator
kicks the momentum towards the left, and lowering the energy kicks the
momentum towards right. One quickly sees that energy will not be conserved.

The simpler model we explore in this section has one important and universal
feature: starting from a state that is a delta function in momentum and
position, and a superposition of two states in the harmonic oscillator, as
expected, the numerical simulations show the initially pure state going
through a process of decoherence, with the parts of the state associated to
different energy levels becoming more and more distinguishable in the phase
space as time goes on. We numerically simulated the dynamics and analytically
found the decoherence rate, for the large $n$ approximation (where $n$ marks
the energy levels of the harmonic oscillator).

We will begin with a simple model with interaction Hamiltonian,
$H_{I}(q,p)=B(q)\left(\omega_{0}aa^{\dagger}-\omega_{1}a^{\dagger}a\right)$.
This simple model is analogous to the one used for the spin-$\frac{1}{2}$
particle. We consider a model with two Lindblad operators, one of which is
proportional to the creation operator and the other is proportional to the
annihilation operator:

$\displaystyle L_{\alpha=1}$ $\displaystyle=\sqrt{|\omega_{0}|}\,a,$ (52a)
$\displaystyle L_{\alpha=2}$ $\displaystyle=\sqrt{|\omega_{1}|}\,a^{\dagger}.$
(52b)

As a result, the functions $h^{\alpha}(q,p)$ are defined as

$\displaystyle h^{\alpha=1}(q,p)$
$\displaystyle=\text{sign}(\omega_{0})\,B(q),$ (53a) $\displaystyle
h^{\alpha=2}(q,p)$ $\displaystyle=\text{sign}(\omega_{1})\,B(q).$ (53b)

In the following, we focus on the case in which $B(q)=qB$, where $B$ is a
constant, and for simplicity we set the coupling constant
$\omega_{0}=-\omega_{1}=1\,\mathrm{s}^{-1}$. Notice that, for this choice of
constants, the interaction Hamiltonian reduces to $H_{I}(q,p)=qB\,\mathbb{I}$,
and the commutator in Eq. (III) drops out. We will discuss more general
models, including two coupled oscillators in Sec. V.

The resulting CQ master equation, where we additionally account for the free-
particle classical Hamiltonian $H_{C}=\frac{p^{2}}{2m}$, is given by

$\frac{\partial\varrho(q,p,t)}{\partial
t}=-\frac{p}{m}\frac{\partial\,\varrho(q,p,t)}{\partial
q}+\frac{1}{\tau}\left(\gamma_{\uparrow}a^{\dagger}\,\varrho(q,p-B\tau,t)\,a+\gamma_{\downarrow}a\,\varrho(q,p+B\tau,t)\,a^{\dagger}-\frac{1}{2}\left\\{{\gamma_{\downarrow}a^{\dagger}a+\gamma_{\uparrow}aa^{\dagger}},{\varrho(q,p,t)}\right\\}_{+}\right),$
(54)

where $\tau$ is the rate of jump in both the classical and quantum degrees of
freedom, and the constants $\gamma_{\uparrow}$, $\gamma_{\downarrow}$
determine the rate of damping of the quantum oscillator compared with the rate
of pumping. For simplicity we henceforth take
$\gamma_{\uparrow}=\gamma_{\downarrow}=1$ although from a physical point of
view, one might want to have the damping term be larger, to drive the quantum
oscillator to its ground state111In a field theory, taking the damping term to
be larger than the anti-damping term, can result in the theory violating
causalityPoulin and Preskill (2017).

In this model, when the energy of the particle is increased with the raising
operator $a^{\dagger}$, its momentum is lowered, and vice versa, when the
energy is lowered with the annihilation operator $a$, its momentum is
increased. It is worth noting that the momentum is here one dimensional; thus,
increasing the momentum does not necessarily mean increasing also its absolute
value, and therefore, for the free-particle Hamiltonian $H_{C}$, we cannot
conserve the energy in this model. An example of the dynamics generated by
this master equation is shown in Fig. 7.

Figure 7: The evolution of the populations $u_{n}(q,p,t)$ of the harmonic
oscillator in phase space,. The system is initially in a superposition over
the state $\ket{n_{1}=12}$ and $\ket{n_{2}=16}$, centered in the origin of the
phase space. The evolution of the populations is shown for
$t=0.01\,\mathrm{s}$, $0.03\,\mathrm{s}$ and $0.5\,\mathrm{s}$. The state is
updated every $\delta t=2.5\cdot 10^{-5}\,\mathrm{s}$, and the jump rate of
the model is $\tau=0.1\,\mathrm{s}$. The constant
$B=1\,\mathrm{J}\cdot\mathrm{s}\cdot\mathrm{m}^{-1}$, and the mass of the
harmonic oscillator is $m=1\,\mathrm{kg}$. Top panels. The evolution of the
population for $n=n_{1}$ is shown. Initially, this population is peaked around
the origin of phase space, and it starts to diffuse in position due to the
action of the rising and lowering operators. At later time, a second packet
can be seen, with momentum
$p=0.4\,\mathrm{kg}\cdot\mathrm{m}\cdot\mathrm{s}^{-1}$. This packet
originated from the population at $n=n_{2}$, and reaches $n=n_{1}$ thanks to
the repeated action of lowering operators. Central panels. The evolution of
the population for $n=\frac{n_{1}+n_{2}}{2}$ is shown. The symmetric packets
with negative and positive momentum ($p=\pm
0.2\,\mathrm{kg}\cdot\mathrm{m}\cdot\mathrm{s}^{-1}$) originated from
$n=n_{1}$ and $n=n_{2}$, respectively. Bottom panels. The evolution of the
population for $n=n_{2}$ is shown. The general behaviour in phase space is
analogous to that of the population for $n=n_{1}$, but with a second packet
with negative momentum.

Because the CQ master equation introduces a new contribution to decoherence,
coming from the interaction between the classical and quantum degrees of
freedom, we now aim to study the decoherence rate of the harmonic oscillator.
In order to analytically find this rate, we have to compute the time evolution
of the CQ state; to simplify the analysis, however, we will neglect the
contribution of the Poisson bracket with the classical Hamiltonian $H_{C}$,
and only focus on the jumping dynamics, which is here responsible for the
decoherence of the system. First, we write out the CQ state $\varrho(q,p,t)$
in terms of the Fock basis,

$\varrho(q,p,t)=\sum_{n,m}u_{n,m}(q,p,t)\ket{n}\bra{m}.$ (55)

We can additionally perform the Fourier transform of the above state with
respect to the momentum $p$. Thus, we get the following equation

$\displaystyle\frac{\partial\tilde{\varrho}(q,x,t)}{\partial
t}=\frac{1}{\tau}\left(a^{\dagger}\,\tilde{\varrho}(q,x,t)\,e^{iB\tau
x}\,a+a\,\tilde{\varrho}(q,x,t)\,e^{-iB\tau
x}\,a^{\dagger}-\frac{1}{2}\left\\{{a^{\dagger}a+aa^{\dagger}},{\tilde{\varrho}(q,x,t)}\right\\}_{+}\right).$
(56)

This step simplifies our analysis, since in this way we remove from Eq. (54)
the non-locality in momentum. However, at the end of our analysis we will have
to perform an inverse-Fourier transform to recover the state in the momentum
domain.

If we sandwich the above equation with two energy eigenstates
$\bra{n}\,\cdot\,\ket{m}$, we obtain the following equation for the elements
$u_{n,m}(q,x,t)$ of the CQ state,

$\displaystyle\frac{\partial u_{n,m}(q,x,t)}{\partial
t}=\frac{\sqrt{nm}}{\tau}\,u_{n-1,m-1}(q,x,t)\,e^{iB\tau
x}+\frac{\sqrt{(n+1)(m+1)}}{\tau}\,u_{n+1,m+1}(q,x,t)\,e^{-iB\tau
x}-\frac{n+m+1}{\tau}\,u_{n,m}(q,x,t).$ (57)

The fact that the matrix element $(n,m)$ interacts only with itself, and its
nearest neighbours (on a line parallel to the matrix diagonal) $(n-1,m-1)$ and
$(n+1,m+1)$ will prove to be convenient. Indeed, this implies that if we are
interested in the problem where the initial state of the CQ density matrix is
a superposition of two energy eigenstates,

$\varrho(q,p,t=0)=\delta(p)\,\delta(q)\,\frac{1}{2}(\ket{n_{1}}+\ket{n_{2}})(\bra{n_{1}}+\bra{n_{2}}),$
(58)

then we only need to solve three differential vector equations, where the
vector takes values along the three diagonal lines involved; the main
diagonal, among whose elements we find $(n_{1},n_{1})$ and $(n_{2},n_{2})$,
the second line parallel to the diagonal marked by the matrix element
$(n_{1},n_{2})$, and the third one marked by $(n_{2},n_{1})$.

###  IV.1 Solving the differential equations for the CQ harmonic oscillator

We now attempt to solve the differential equation with the previously
mentioned initial conditions, where we know that the three diagonals along
which we have to solve are the ones marked with the matrix elements
$(n_{2},n_{1})$, $(n_{1},n_{1})$ and $(n_{1},n_{2})$. For example, we can
focus on the differential equation of the third diagonal,

$\displaystyle\frac{\partial u_{n_{1}+k,n_{2}+k}(q,x,t)}{\partial t}=$
$\displaystyle\frac{\sqrt{(n_{1}+k)(n_{2}+k)}}{\tau}\,u_{n_{1}+k-1,n_{2}+k-1}(q,x,t)\,e^{iB\tau
x}$ $\displaystyle-$
$\displaystyle\frac{n_{1}+n_{2}+1}{\tau}\,u_{n_{1}+k,n_{2}+k}(q,x,t)$
$\displaystyle+$
$\displaystyle\frac{\sqrt{(n_{1}+k+1)(n_{2}+k+1)}}{\tau}\,u_{n_{1}+k+1,n_{2}+k+1}(q,x,t)\,e^{-iB\tau
x}.$ (59)

Assuming that $n_{1}<n_{2}$, we can now rearrange the functions
$u_{n_{1}+k,n_{2}+k}(q,x,t)$, for $k\in\\{-n_{1},-n_{1}+1,\ldots\\}$, into a
vector

$\vec{v}=(u_{0,n_{2}-n_{1}},u_{1,n_{2}-n_{1}+1},...)^{T}=\sum_{k}u_{n_{1}+k,n_{2}+k}\,\vec{e}_{k},$
(60)

where the $\\{\vec{e}_{k}\\}_{k}$ are unit vectors, and we suppressed the time
and phase space dependence for simplicity. The goal is then to find the
solution of the equation

$\frac{\partial\vec{v}}{\partial t}=\mathcal{M}\,\vec{v},$ (61)

where $\mathcal{M}$ is a tridiagonal matrix, whose entries are determined by
Eq. (59). The size of this matrix is infinite, but we can obtain an
approximate solution to the above equation by considering only part of the
elements in the vector $\vec{v}$. In particular, we require the index $k$ to
take values in an interval $\left(-\frac{N}{2},\frac{N}{2}\right]$, for a
given integer $N$. This approximation is acceptable for short enough times,
when the coherence $u_{n_{1},n_{2}}$ has not yet spread too much. By imposing
these lower and upper bounds to the $k$ index, we make the system effectively
finite dimensional, and therefore the finite-dimensional matrix $\mathcal{M}$
can be diagonalised. In the following, we refer to the $i$-th eigenvector of
this matrix as $\vec{w}_{i}$, with eigenvalue $\lambda_{i}$.

The equation regulating the dynamics of the eigenvectors $\vec{w}_{i}$ is then
given by

$\frac{\partial\vec{w}_{i}}{\partial t}=\lambda_{i}\vec{w}_{i},$ (62)

and the solution is given by
$\vec{w}_{i}(t)=e^{\lambda_{i}\,t}\vec{w}_{i}(0)$. The same procedure can be
followed for the main diagonal, marked with the matrix elements
$(n_{1},n_{1})$ and $(n_{2},n_{2})$, so as to obtain the evolution for the
populations. The diagonalisation of the matrix $\mathcal{M}$ can be performed
numerically, so as to obtain the decay rates for the coherences of the
(Fourier transformed) CQ state. However, obtaining the evolution of the actual
CQ state $\rho(q,p,t)$ would require us to Fourier transform back to the
momentum domain, which can be computational expensive, given that such
transformation should be performed on every element $u_{n,m}(q,x,t)$ of the CQ
state. For this reason, we preferred to make use of the unravelling technique
to obtain the evolution of the CQ harmonic oscillator, and to compute the
dechoerence rate.

###  IV.2 Decoherence rate in the large $n$ approximation

We now compute an analytical expression for the decoherence rate in the case
in which the Fock numbers $n_{1}$ and $n_{2}$ are large compared to the range
of the $k$ index. In this situation, we can approximate the matrix elements
along the diagonals of $\mathcal{M}$ to be equal, so as to obtain a Toeplitz
matrix. Concretely, if the $k$ index takes value between
$\left(-\frac{N}{2},\frac{N}{2}\right]$, where $N\ll n_{1},n_{2}$ is an even
number, then the matrix $\mathcal{M}$ can be approximated as

$\mathcal{M}\approx\begin{pmatrix}\dots&\dots&0&\dots&\dots\\\
\dots&-\frac{n_{1}+n_{2}}{\tau}&\frac{e^{i\phi}}{\tau}\sqrt{n_{1}n_{2}}&0&\dots\\\
0&\frac{e^{-i\phi}}{\tau}\sqrt{n_{1}n_{2}}&-\frac{n_{1}+n_{2}}{\tau}&\frac{e^{i\phi}}{\tau}\sqrt{n_{1}n_{2}}&\dots\\\
0&0&\frac{e^{-i\phi}}{\tau}\sqrt{n_{1}n_{2}}&-\frac{n_{1}+n_{2}}{\tau}&\dots\\\
\dots&0&0&\dots&\dots\end{pmatrix},$ (63)

where we define $e^{i\phi}=e^{iB\tau x}$ for convenience. We can now find the
eigenvalues $\lambda_{m}$ of the above Toeplitz matrix, which are give by
Noschese et al. (2013)

$\lambda_{m}=-\frac{n_{1}+n_{2}}{\tau}+\frac{2}{\tau}\sqrt{n_{1}n_{2}}\cos{\frac{m\,\pi}{N+1}},\quad
m=1,\ldots,N.$ (64)

The eigenvectors of the matrix $\left\\{\vec{\omega}_{m}\right\\}^{N}_{m=1}$
are instead given by

$\vec{\omega}_{m}=\sum_{r=1}^{N}B_{m,r}\,\vec{e}_{r-\frac{N}{2}},$ (65)

where $B_{m,r}=\sqrt{\frac{2}{N+1}}\,e^{-ir\phi}\,\sin{\frac{m\,r\,\pi}{N+1}}$
is the unitary matrix diagonalizing $\mathcal{M}$.

Under the dynamics described by Eq. (61), the eigenvectors evolve as
$\vec{\omega}_{m}(t)=e^{\lambda_{m}\,t}\,\vec{\omega}_{m}(0)$. Thus, by re-
writing the vector $\vec{v}$, defined in Eq. (60), in terms of the
eigenvectors of $\mathcal{M}$, we can compute the approximate time evolution
of the coherence of the state, which is given by

$u_{n_{1}-\frac{N}{2}+\ell,n_{2}-\frac{N}{2}+\ell}(q,x,t)=\sum_{r,m=1}^{N}u_{n_{1}-\frac{N}{2}+r,n_{2}-\frac{N}{2}+r}(q,x,0)\,B^{-1}_{r,m}\,e^{\lambda_{m}\,t}\,B_{m,\ell},\quad\ell=1,\ldots,N.$
(66)

There is a lot of debate in the literature on whether one can consistently
couple quantum systems and classical ones, and many ways of doing so which
although useful in some regimes, are pathological Bohr and Rosenfeld (1933);
DeWitt and Rickles (2011); DeWitt (1962); Eppley and Hannah (1977); Caro and
Salcedo (1999); Salcedo (1996); Sahoo (2004); Terno (2006); Salcedo (2012);
Barceló et al. (2012)

We can now specialise the above approximate solution to the initial condition
we introduced in Eq. (58). In this case, the sole contribution to the
coherence at time $t=0$ is
$u_{n_{1},n_{2}}(q,x,0)=\frac{1}{2}\,\delta(q)\delta(p)=\frac{1}{2}\,\delta(q)\int
dxe^{ipx}$, and its time evolution is given by

$u_{n_{1},n_{2}}(q,x,t)=\frac{e^{-\frac{n_{1}+n_{2}}{\tau}t}}{N+1}\sum_{m=1}^{N}\sin^{2}\left(\frac{N}{2(N+1)}\,m\,\pi\right)e^{\frac{2}{\tau}\sqrt{n_{1}n_{2}}\cos\left({\frac{m\pi}{N+1}}\right)t}\,\delta(q)\delta(p),$
(67)

where the dominant rate of decoherence is given by
$\Gamma=\frac{n_{1}+n_{2}}{\tau}$, with corrections of the order
$\sqrt{n_{1}n_{2}}$. The above analytical result can be compared to the
numerical one obtained with the unravelling technique, see Fig. 8.

Figure 8: The evolution of the coherence for two different choices of initial
state for the CQ system. In both cases, the system is evolved up to the time
$t=2.5\cdot 10^{-2}\,\mathrm{s}$, with time steps of $\delta t=2.5\cdot
10^{-6}\,\mathrm{s}$, and a jumping rate $\tau=2.5\cdot 10^{-1}\,\mathrm{s}$.
The interaction constant is
$B=1\,\mathrm{J}\cdot\mathrm{s}\cdot\mathrm{m}^{-1}$, and the mass is
$m=1\,\mathrm{kg}$. Left. The initial quantum state is an even superposition
between $\ket{n_{1}=15}$ and $\ket{n_{1}=30}$. The evolution of the quantum
coherence $\left|u_{m_{1},m_{2}}\right|$ is shown for different values of
$m_{1}$ and $m_{2}$, and the classical degrees of freedom have been traced
out. The numerical evolution is represented with the solid line, and the
dashed line is the analytical solution given in Eq. (66), for $N=10$. Centre.
The initial quantum state is a even superposition between $\ket{n_{1}=22}$ and
$\ket{n_{1}=23}$, and the evolution of the coherence is shown for those
$m_{1}$ and $m_{2}$ close to $n_{1}$ and $n_{2}$. Again, the solid lines are
obtained through the numerical simulation, using the unravelling code, and the
dashed ones are the (approximate) analytical solutions. Right The evolution of
the coherence $\left|u_{n_{1},n_{2}}\right|$ for the two choices of initial
state. When the difference between $n_{1}$ and $n_{2}$ is big, the coherence
reaches lower values than in the case in which $n_{1}$ and $n_{2}$ are close
together; this behaviour of the coherence is analogous to the case in which
the dynamics is fully quantum, and the system decoheres due to the coupling
with an external bath, rather than to the interaction with the classical
degrees of freedom.

Notice that, in order to obtain the decoherence as a function of phase space
variables, we have to perform the inverse Fourier transform on the coherence
given in Eq. (67). However, since the time evolution of $u_{n_{1},n_{2}}$ is
constant in the $x$-space, we find that in $p$-space the solution picks up a
delta function $\delta(p)$, and the decoherence rate is not modified. The
presence of a delta function for the momentum is consistent with the dynamics
we are considering. Indeed, from Eq. (54) it is easy to see that the coherence
$u_{n_{1},n_{2}}$ is non-zero for $p=0$ only, while for any $k\neq 0$ the
coherence $u_{n_{1}+k,n_{2}+k}$ is non-zero for $p=-kB\tau$. Likewise, the
position dependence is described by $\delta(q)$, since the master equation we
are considering does not include the classical Hamiltonian
$H_{C}=\frac{p^{2}}{2m}$.

It is worth noting that the evolution of the coherence induced by the CQ
dynamics described in Eq. (54) is analogous to the evolution induced by a
fully quantum master equation where the phase space degrees of freedom are
neglected. The decoherence can however be interpreted in two different ways;
for the CQ dynamics, this is due to the interaction with the classical degrees
of freedom, and the fact that (quantum) information is leaking in the
classical phase space. For the fully quantum dynamics, instead, the
decoherence is due to the interaction with an external environment. The
coherence evolution for these two settings can still be distinguished if we
introduce the classical Hamiltonian $H_{C}$ in the CQ dynamics, since this
term provides a non-trivial dependence on the position $q$ in phase space. As
a final remark, we notice that coherence in this model is “moving” from the
initial superposition between $\ket{n_{1}}$ and $\ket{n_{2}}$ to superposition
of nearby Fock states. This is shown in the left and centre plots of Fig. 8;
the coherence of the initial state $\ket{n_{1}}+\ket{n_{2}}$ monotonically
decreases, while the coherence of the superposition $\ket{n_{1}\pm
1}+\ket{n_{2}\pm 1}$ increases at early times and then start decreasing as
well.

##  V Discussion

We have here considered general dynamics which consistently couples a quantum
system with classical degrees of freedom. The combined system is described by
a hybrid-density which gives the probability of the system to be at any point
in phase space and if so, what the density matrix of the quantum system is. We
introduced an unravelling approach to solve the equations of motion. In this
approach, both the classical and quantum systems evolve along a particular
trajectory with some probability, and when we average over these trajectories,
we recover the dynamics for the hybrid-density. Unlike the purely quantum
case, the trajectories considered here can lead to objectively certain state
transformations, in which the quantum system transitions to a particular pure
state, conditioned on the classical degrees of freedom. This happens when the
Lindblad operators $L_{\alpha}$ are uniquely determined by the classical
transition $W^{\alpha}(z|z-\Delta)$ rate, so that the observation of a
particular classical transition unambiguously informs us about the quantum
transition. Since the classical system has a definite trajectory, the
trajectory of the quantum system is also definite. This is contrast to the
purely quantum case, where the decomposition of a density matrix into an
ensemble of pure states is not generally unique. A particularly simple and
illustrative example is the master equation

$\displaystyle\frac{\partial\varrho}{\partial
t}=L_{\alpha}\rho(q,p-\Delta_{\alpha})L^{\dagger}_{\alpha}-\frac{1}{2}\left\\{{L^{\dagger}_{\alpha}L_{\alpha}},{\varrho(q,p)}\right\\}_{+}$
(68)

where, if the $\Delta_{\alpha}$ are different for each $\alpha$, observing a
sudden jump in momentum by an amount $\Delta_{\alpha}$ implies that the
quantum state has changed from $\ket{\psi}$ to a suitably normalized
$L_{\alpha}\ket{\psi}$. The classical system encodes the trajectory of the
quantum system.

Furthermore, there is a very real sense in which the quantum system jumps from
being in a superposition of several states, to being in one of those states.
Ironically, this is in contrast to spontaneous collapse models Ghirardi et al.
(1986); Pearle (1989); Ghirardi et al. (1990) where collapse is meant to be
caused by an unobservable field. In these theories, as in ordinary quantum
theory, the density matrix of the quantum system fully describes the system,
and so, one cannot distinguish between decoherence (which we think of as the
gradual decay of off-diagonal elements of the density matrix in a chosen
basis), and the case where the wavefunction has a probability per unit time of
being projected to a particular states, such that on average, off-diagonal
matrix elements of the density matrix decay with time. In both cases, the
density matrix of the quantum system is the same. In contrast, in the theory
studied here, by monitoring the classical degree of freedom one can determine
that a sudden change has occurred.

The system which is treated classically could be merely a large quantum system
that is effectively a classical system such as an environment or a measuring
apparatus. But it could also be fundamentally classical, for example, space-
time, as in Refs. Hall and Reginatto (2005); Poulin and Preskill (2017);
Oppenheim (2018), or it could be a system whose classical nature is somehow
assumed, such as the experimenter herself. This provides a potential solution
to the measurement or ”reality” problem of quantum theory. Namely, quantum
systems in superposition collapse to a particular state in that superposition
because they interact with a system which is classical. Which set of possible
states they can collapse to (sometimes called the pointer basis) is determined
by the interaction.

Turning towards more pedestrian phenomena, we have found that the toy models
we studied exhibited a trade-off between the amount of classical diffusion
that occurred in phase space and the amount of quantum decoherence or noise
that occurred in the quantum system. To minimise the disturbance on the
quantum state, there needs to be a lot of diffusion in the classical system.
In Ref. Oppenheim et al. (2020) we will formally prove that this is indeed the
case.

We also found that because the Hamiltonian was not the generator of time-
translations, it typically is not conserved Banks et al. (1984). We saw that
this effect does not have to be large, and there are a number of ways to keep
it under control so that we can have stable solutions. Violations of energy
conservation present a possible experimental signature of fundamental hybrid
dynamics, but also feature in spontaneous collapse models, and gravitationally
induced decoherence Diosi (1987); Penrose (1996)222Although both Refs. Diosi
(1987); Penrose (1996) and Hall and Reginatto (2005); Poulin (2017); Oppenheim
(2018) consider gravity’s role in decoherence or state collapse, the reasons
put forward for the two effects are very different. The former is usually
attributed to non-linearities in fully quantum evolution law, the latter due
to the gravitational field being treated classically.. However, the way these
violations appear in the theories discussed here is different; indeed, the
dynamics is not Markovian if one considers just the quantum system, because
the classical system can act as a memory (likewise the quantum system acts as
a memory for dynamics of the classical system). Furthermore, here the effect
is renormalised in terms of the interaction itself Oppenheim et al. (2020). A
full study on the issue of energy conservation is beyond the scope of the
current work, but in the following we present a brief discussion of the
problem, at leading order.

It is instructive to recall the analogous issue in classical mechanics. If one
has diffusion as in the Fokker-Planck equation, then the magnitude of the
momentum will gradually drift and increase, resulting in energy increase. The
introduction of a friction term prevents this, by causing the system to settle
down to an equilibrium state. The Fokker-Planck equation for the phase space
density $\rho(q,p)$ is of the form

$\displaystyle\frac{\partial\,\rho(q,p)}{\partial t}=$
$\displaystyle\left\\{{H_{C}(q,p)},{\rho}\right\\}+\frac{\partial^{2}D(q,p)\rho(q,p)}{\partial
p^{2}}+\frac{\partial}{\partial p}\big{[}\eta(q,p)\rho(q,p)\big{]},$ (69)

with the friction coefficient $\eta:=\gamma\frac{\partial H_{C}}{\partial p}$.
This is just the theory of Brownian motion as developed by Einstein and
Smoluchowski, with a thermal stationary distribution given by

$\displaystyle\rho_{\beta}=\frac{e^{-\beta H_{C}}}{Z},$ (70)

where $Z$ is the partition function and the inverse temperature is

$\displaystyle\beta=\gamma/D.$ (71)

By adding a friction term, we can stop the diffusion in the classical degrees
of freedom.

For the type of hybrid model presented in Sec. III.1, with Lindblad operators
that are a basis of orthogonal projectors
$L_{\alpha}=|\alpha\rangle\\!\langle\alpha|$ onto energy eigenstates of
$H_{Q}$, and $h^{\alpha}=\omega Bq$ (to take a simple example) we may add a
friction term such as

$\displaystyle\cal{F}(\rho)$ $\displaystyle=\frac{\partial}{\partial
p}\big{[}\gamma(q)\frac{\partial H_{C}}{\partial
p}\rho\big{]}=\frac{\gamma(q)}{m}\Big{[}p\frac{\partial\rho}{\partial
p}+\rho\Big{]}$ (72)

to the right hand side of the evolution law of Eq. (III). The friction
coefficient $\gamma$ can be chosen to depend on $q$. To first order in $\tau$,
the diffusion term is given by Eq. (36), and so from Eq. (71), we expect that
to first order, the classical system will equilibriate to a thermal state at
inverse temperature $\beta=\gamma(q)/(B\omega\langle\alpha\rangle)^{2}\tau$.
Because the Lindblad operators project onto eigenstates of $H_{Q}$, the
quantum system also equilibrates by decohering into one of it’s energy
eigenstates.

The case where the Lindblad operations cause noise is more interesting, since
then one has diffusion not only in phase space but also in the Hilbert space.
What’s more, the amount of diffusion in phase space can depend on the quantum
state; for example, for a quantum harmonic oscillator the diffusion can grow
with $n$, the energy level of the oscillator. This means that the diffusion
term can be arbitrarily large, and overcome any classical friction term. We
would thus like the friction term to also increase with $n$ so that the
classical system can reach a steady state. To explore this, let us consider a
more sophisticated model than Eq. (69), where a classical harmonic oscillator
with Hamiltonian $H_{C}=\omega_{c}(p^{2}+q^{2})/2$ is coupled to a quantum
harmonic oscillators with free Hamiltonian $H_{Q}=\omega_{Q}a^{\dagger}a$; the
Lindblad operators are now $L_{\alpha}=a$, the annihilation operator. A
natural coupling between the two oscillators would be to take the interaction
Hamiltonian to be $H_{I}=-2\omega_{c}\omega_{Q}qQ$ with
$Q=\frac{1}{\sqrt{2}}(a+a^{\dagger})$, but it is hard to make the hybrid
dynamics of such a coupling completely positive. An alternative, is to take
the quantum part to be positive definite (while leaving the classical part
arbitrary), e.g. $H_{I}=B(q)Q^{2}$, which is still a local coupling. Let us be
slightly more general, and consider coupling to $Q^{2}+\kappa P^{2}$ with
$P=\frac{i}{\sqrt{2}}(a^{\dagger}-a)$ and $\kappa\geq 0$. It will prove
convenient to use $[a,a^{\dagger}]=1$ to rewrite this as

$\displaystyle
H_{I}=B(q)\frac{1}{2}\Big{(}(1-\kappa)(aa+a^{\dagger}a^{\dagger})+(1+\kappa)(a^{\dagger}a+1)\Big{)}$
(73)

where we have dropped the purely classical terms from the interaction. A
master equation for such a theory is then

$\displaystyle\frac{\partial\varrho(q,p)}{\partial t}$
$\displaystyle=\\{H_{c}(q,p),\varrho(q,p)\\}+\mathcal{F}(\varrho(q,p))-i\omega_{Q}\big{[}a^{\dagger}a,\varrho(q,p)\big{]}-iB(q)\big{[}Q^{2}+kP^{2},\varrho(q,p,t)\big{]}$
$\displaystyle+\int d\Delta
W_{\uparrow}(q,p|q,p-\Delta)a^{\dagger}\,\varrho(q,p-\Delta)\,a+\int d\Delta
W_{\downarrow}(q,p|q,p-\Delta)a\,\varrho(q,p-\Delta)\,a^{\dagger}$
$\displaystyle-\frac{1}{2}\left\\{{W_{\downarrow}(q,p)\,a^{\dagger}a+W_{\uparrow}(q,p)\,aa^{\dagger}},{\varrho(q,p)}\right\\}_{+}$
(74)

where we are still to specify the transition rates $W_{\uparrow}$ and
$W_{\downarrow}$, associated to the pumping and dumping of the quantum
oscillator, respectively, and the classical friction term $\mathcal{F}$. The
last term in the first line is the influence of the classical system on the
harmonic oscillator, while the second and third lines is the influence of the
quantum system on the classical one. In the following, we will require the
expansion of $W_{\uparrow}$ and $W_{\downarrow}$ to contain a term
proportional to $\tau B^{\prime}(q)\frac{\partial\varrho(q,p)}{\partial p}$,
so as to balance the energy increase of the classical system due to the
interaction with the quantum one.

On the other hand, at leading order in $\Delta$ the second and third lines
give rise to a Lindblad equation with a damping and pumping term,

$\displaystyle\mathcal{D}(\varrho):=W_{\downarrow}(z)\,a^{\dagger}\,\varrho\,a+W_{\downarrow}(z)\,a\,\varrho\,a^{\dagger}-\frac{1}{2}\left\\{{W_{\downarrow}(z)a^{\dagger}a+W_{\uparrow}(z)aa^{\dagger}},{\varrho}\right\\}_{+}$
(75)

and so when $W_{\uparrow}(z)\geq W_{\downarrow}(z)$, the state of the harmonic
oscillator will continue to increase in energy. The quantum analogue of adding
a friction term is to take the damping term to be strictly larger than the
pumping term, and if the classical system were to remain at rest, the quantum
oscillator will thermalise to an inverse temperature proportional to
$\log{W_{\downarrow}(z)/W_{\uparrow}(z)}$, at which point one can verify that
$\mathcal{D}(\varrho)=0$. We can thus have the system reach equilibrium by
taking $\log{W_{\downarrow}(z)/W_{\uparrow}(z)}$ to be strictly positive and
adding a purely classical friction term $\mathcal{F}(\varrho)$.

If on the other hand, we take $B(q)$ and the rates $W_{\downarrow}(z)$ and
$W_{\uparrow}(z)$ such that the classical system equilibriates in a region of
phase space with sufficient probability of having $W_{\downarrow}(z)\leq
W_{\uparrow}(z)$, then energy can be pumped into the harmonic oscillator
indefinitely. In such a case, we would like to control this, so that it
happens slowly. To explore this, let us consider the specific realisations,

$\displaystyle\int d\Delta W_{\uparrow}(p|p-\Delta)\,\varrho(q,p-\Delta)$
$\displaystyle=\frac{1}{\tau}h_{\uparrow}(q)e^{\tau\partial_{p}(X_{\uparrow}\cdot)}\varrho(q,p),$
(76a) $\displaystyle\int d\Delta
W_{\downarrow}(p|p-\Delta)\,\varrho(q,p-\Delta)$
$\displaystyle=\frac{1}{\tau}h_{\downarrow}(q)e^{\tau\partial_{p}(X_{\downarrow}\cdot)}\varrho(q,p),$
(76b)

where the differential operator $\partial_{p}(X\cdot)$ acts over an hybrid
operator $\rho(q,p)$ as $\partial_{p}(X\,\rho(q,p))$. The two realisations are
positive as long as $h_{\uparrow}$ and $h_{\downarrow}$ is everywhere
positive, and the master equation will preserve normalisation with
$W_{\uparrow}(q,p)=h_{\uparrow}(q)/\tau$ and
$W_{\downarrow}(q,p)=h_{\downarrow}(q)/\tau$. Let us now choose
$X_{\downarrow}=\Delta_{\downarrow}(q)+\eta_{\downarrow}$,
$X_{\uparrow}=\Delta_{\uparrow}(q)+\eta_{\uparrow}$, where $\eta_{\uparrow}$
and $\eta_{\downarrow}$ are friction terms as in Eq. (69). Then, at zeroth
order in $\tau$ we have the Lindblad equation (75) with temperature
proportional to $\log{h_{\downarrow}(q)/h_{\uparrow}(q)}$. However, even if we
don’t have $h_{\downarrow}(q)\geq h_{\uparrow}(q)$, the rate at which energy
is pumped into the oscillator is

$\delta H_{Q}=\omega\frac{d{\langle
a^{\dagger}a\rangle}}{dt}=\frac{2}{\tau}\omega^{2}(h_{\uparrow}(q)-h_{\downarrow}(q))n-\frac{\omega}{\tau}h_{\uparrow}(q),$
(77)

which is easily calculated from the master equation. This term can potentially
be made small in the phase space region where the classical system
equilibriates. The next order terms give,

$\displaystyle\ldots$
$\displaystyle+h_{\downarrow}(q)\Delta_{\downarrow}(q)\,a\frac{\partial\varrho}{\partial
p}a^{\dagger}+h_{\uparrow}(q)\Delta_{\uparrow}(q)\,a^{\dagger}\frac{\partial\varrho}{\partial
p}a+\frac{\tau}{2}h_{\downarrow}(q)\Delta_{\downarrow}(q)^{2}\,a\frac{\partial^{2}\varrho}{\partial
p^{2}}a^{\dagger}+\frac{\tau}{2}h_{\uparrow}(q)\Delta_{\uparrow}(q)^{2}\,a^{\dagger}\frac{\partial^{2}\varrho}{\partial
p^{2}}a$
$\displaystyle+h_{\downarrow}(q)\,a\frac{\partial\gamma_{\downarrow}\varrho}{\partial
p}a^{\dagger}+h_{\uparrow}(q)\,a^{\dagger}\frac{\partial\gamma_{\uparrow}\varrho}{\partial
p}a+\ldots$ (78)

where we take the last two terms to be of roughly the same order as the
previous diffusion terms.

The $\kappa=1$ case is a natural one to consider, since it is reminiscent of
the master equation of general relativity coupled to a scalar field Oppenheim
(2018). The classical oscillator couples to the energy of the quantum one, as
the classical gravitational field couples to the scalar field expanded in
terms of momentum modes. Likewise, analogy with the gravitational case
suggests taking the pure commutator term of Eq. (V) since gravity is always
coupling to the energy. In this case, we can take

$\displaystyle
h_{\downarrow}(q)\Delta_{\downarrow}(q)=h_{\uparrow}(q)\Delta_{\uparrow}(q)=\frac{\omega}{2}B^{\prime}(q)$
(79)

so that the first two terms in Eq. (78) will give the Poisson bracket
$\\{B(q)\langle\frac{1}{2}Q^{2}+\frac{1}{2}P^{2}\rangle,\varrho\\}$ in the
classical limit. When we put this all together, the master equation (V)
expands as

$\displaystyle\frac{\partial\varrho(q,p)}{\partial t}$
$\displaystyle=\\{H_{c}(q,p),\varrho\\}+\mathcal{F}(\varrho)-iB(q)\,\omega\,\big{[}a^{\dagger}a,\varrho(q,p,t)\big{]}$
$\displaystyle+\frac{\omega}{\tau}\Big{(}h_{\uparrow}(q)a^{\dagger}\,\varrho\,a+h_{\downarrow}(q)a\,\varrho\,a^{\dagger}-\frac{1}{2}\left\\{{\left(h_{\uparrow}(q)+h_{\downarrow}(q)\right)a^{\dagger}a+h_{\downarrow}(q)},{\varrho}\right\\}_{+}\Big{)}+\frac{B^{\prime}(q)\,\omega}{2}\Big{(}a^{\dagger}\,\frac{\partial\varrho}{\partial
p}\,a+a\,\frac{\partial\varrho}{\partial p}\,a^{\dagger}\Big{)}$
$\displaystyle+\frac{\tau}{4}B^{\prime}(q)\omega\Big{(}\Delta_{\uparrow}(q)a^{\dagger}\,\frac{\partial^{2}\varrho}{\partial
p^{2}}\,a+\Delta_{\downarrow}(q)a\,\frac{\partial^{2}\varrho}{\partial
p^{2}}\,a^{\dagger}\Big{)}+\omega\Big{(}h_{\uparrow}(q)a^{\dagger}\,\frac{\partial\eta_{\uparrow}\varrho}{\partial
p}\,a+h_{\downarrow}(q)a\,\frac{\partial\eta_{\downarrow}\varrho}{\partial
p}\,a^{\dagger}\Big{)}$ (80)

In analogy with the field theory case, we would like to take
$h_{\downarrow}(q)=h_{\uparrow}(q)=h(q)$, since in that context, the theory is
local with the consequence that the oscillator will not equilibrate and will
instead increase in energy Poulin and Preskill (2017). This can be slowed down
by any desired amount by taking $h(q)/\tau$ to be small. However, anticipating
the discussion in Ref. Oppenheim et al. (2020), slowing the energy increase
down results in a greater amount of classical diffusion. In particular, the
requirement that Eq. (79) be satisfied, which in this case is

$\displaystyle h(q)\Delta(q)=\frac{\omega}{2}B^{\prime}(q)$ (81)

encodes this trade-off. The right hand side governs the rate of evolution of
$p$ and is set by the dynamics, and the left hand side then requires that a
realisation of this form have a trade-off between $h(q)$ and $\Delta(q)$.

Let us consider the simplest example of $B(q)=Gq$, where $G$ is a positive
constant, so that the classical system acts as if it is subject to a constant
force $F\approx\omega Gn$ when the quantum system is in the Fock state
$\ket{n}$, with $n$ large. The rate at which energy is pumped into the quantum
harmonic oscillator is $\delta H_{Q}=\omega h(q)/\tau$, as shown in Eq. (77).
Recalling the Fokker-Planck equation, the effective diffusion coefficient $D$
is the term in front of $\frac{1}{2}\frac{\partial^{2}\varrho}{\partial
p^{2}}$ in Eq. (V), going as $D\approx F\tau\Delta(q)/2$. Then, Eq. (81) gives
$\delta H_{Q}D\approx\frac{1}{4n}\omega F^{2}$ so that a small amount of
energy being pumped into the oscillator requires large diffusion in the
classical system in relation to the total back-reaction exerted by the quantum
system on the classical one. As $n$ increases the relative trade-off, when the
force $F$ is held fixed, becomes less pronounced.

There appears to be significant freedom in the choice of $h(q)$ and
$\Delta(q)$, the main requirement other than Eq. (81) being that $h(q)$ be
positive. One can even take $h(q)=|B^{\prime}(q)|$ and $\Delta(q)$
proportional to $\mathrm{sign}(B^{\prime}(q))$. If on the other hand,
$\Delta(q)$ is large, we can then get energy being pumped into the classical
system, since large amounts of diffusion in momentum will increase the energy
of the classical system. This can be controlled by increasing the friction
term, to keep the momentum low. However, since the diffusion term goes like
$D\propto n\omega G\tau\Delta(q)$ and thus increases with $n$, the friction
term would also need to scale with $n$ or risk getting overwhelmed by the
diffusion. For this reason, the purely classical friction term
$\mathcal{F}(\varrho)$ will not be sufficient, while the hybrid friction terms
in Eq. (V),

$\mathcal{F}_{H}(\varrho)=\omega\Big{(}h(q)a^{\dagger}\,\frac{\partial\eta_{\uparrow}\varrho}{\partial
p}\,a+h(q)a\,\frac{\partial\eta_{\downarrow}\varrho}{\partial
p}\,a^{\dagger}\Big{)},$ (82)

should allow us to compensate the increase. This is seen explicitly by tracing
out the oscillator, resulting in dynamics for the classical phase space
density $\rho_{n}$ conditioned on the the oscillator being in state $n$

$\displaystyle\frac{\partial\rho_{n}(q,p)}{\partial t}$
$\displaystyle=\\{H_{c},\rho_{n}\\}+\mathcal{F}(\rho_{n})+B^{\prime}(q)\omega(n+\frac{1}{2})\frac{\partial\rho_{n}}{\partial
p}+\frac{\tau\Delta}{2}B^{\prime}(q)\omega(n+\frac{1}{2})\frac{\partial^{2}\rho_{n}}{\partial
p^{2}}+2h(q)\omega(n+\frac{1}{2})\frac{\partial\eta\rho_{n}}{\partial p}+...$
(83)

where we have taken $\eta_{\uparrow}=\eta_{\downarrow}=\eta$ for simplicity.
To this order, this is just the Fokker-Planck equation, and if
$\eta=\tilde{\gamma}\frac{\partial H_{c}}{\partial p}$, then we expect the
effective diffusion coefficient to be
$D\approx\tau\Delta(q)B^{\prime}(q)\omega n$, and the effective friction
coefficient to be $\gamma\approx nh(q)\tilde{\gamma}(q)$. Eq. (71) would then
give the effective inverse “temperature” that the classical system is held at,
although here, since it could depend on $q$, it is not a temperature but still
defines the equilibrium state
$\beta=h(g)\tilde{\gamma}(q)/\tau\Delta(q)B^{\prime}(q)\omega$, which is now
independent of $n$. Note that this system could act as if the classical system
is being held at one temperature while the oscillator is held at a different
temperature. Were we to have access to the environment that each system
appears immersed in, we would presumably see heat flow from one environment to
the other. Understanding how this would work in detail, either analytically or
numerically would be an interesting research direction to pursue. One would
like to better understand the various trade-offs involved in systems like
this, and what realisations lead to equilibrium or steady states, or states
which only change slowly.

## Acknowledgements

We would like to thank Robert Alicki and Joan Camps for valuable discussions.
JO is supported by an EPSRC Established Career Fellowship, a Royal Society
Wolfson Merit Award, C.S. and Z.W.D. acknowledges financial support from
EPSRC. This research was supported by the National Science Foundation under
Grant No. NSF PHY11-25915 and by the Simons Foundation It from Qubit Network.
Research at Perimeter Institute is supported in part by the Government of
Canada through the Department of Innovation, Science and Economic Development
Canada and by the Province of Ontario through the Ministry of Economic
Development, Job Creation and Trade.

## References

  * Bohr and Rosenfeld (1933) N. Bohr and L. Rosenfeld, Quantum theory and measurement pp. 478–522 (1933).
  * DeWitt and Rickles (2011) C. M. DeWitt and D. Rickles, _The role of gravitation in physics: Report from the 1957 Chapel Hill Conference_ , vol. 5 (epubli, 2011).
  * DeWitt (1962) B. S. DeWitt, Journal of Mathematical Physics 3, 619 (1962).
  * Eppley and Hannah (1977) K. Eppley and E. Hannah, Foundations of Physics 7, 51 (1977).
  * Caro and Salcedo (1999) J. Caro and L. Salcedo, Physical Review A 60, 842 (1999).
  * Salcedo (1996) L. Salcedo, Physical Review A 54, 3657 (1996).
  * Sahoo (2004) D. Sahoo, Journal of Physics A: Mathematical and General 37, 997 (2004).
  * Terno (2006) D. R. Terno, Foundations of Physics 36, 102 (2006).
  * Salcedo (2012) L. Salcedo, Physical Review A 85, 022127 (2012).
  * Barceló et al. (2012) C. Barceló, R. Carballo-Rubio, L. J. Garay, and R. Gómez-Escalante, Physical Review A 86, 042120 (2012).
  * Marletto and Vedral (2017) C. Marletto and V. Vedral, npj Quantum Information 3, 29 (2017).
  * Boucher and Traschen (1988) W. Boucher and J. Traschen, Physical Review D 37, 3522 (1988).
  * Diósi et al. (2000) L. Diósi, N. Gisin, and W. T. Strunz, Physical Review A 61, 022108 (2000).
  * Anderson (1995) A. Anderson, Physical Review Letters 74, 621 (1995).
  * Kapral and Ciccotti (1999) R. Kapral and G. Ciccotti, The Journal of chemical physics 110, 8919 (1999).
  * Kapral (2006a) R. Kapral, Annu. Rev. Phys. Chem. 57, 129 (2006a).
  * Prezhdo and Kisil (1997) O. V. Prezhdo and V. V. Kisil, Physical Review A 56, 162 (1997).
  * Kapral (2006b) R. Kapral, Annu. Rev. Phys. Chem. 57 (2006b).
  * Page and Geilker (1981) D. N. Page and C. Geilker, Physical Review Letters 47, 979 (1981).
  * Gisin (1989) N. Gisin, Helv. Phys. Acta 62, 363 (1989).
  * Blanchard and Jadczyk (1993) P. Blanchard and A. Jadczyk, Physics Letters A 175, 157 (1993).
  * Alicki and Kryszewski (2003) R. Alicki and S. Kryszewski, Physical Review A 68, 013809 (2003).
  * Diósi (2014) L. Diósi, Physica Scripta 2014, 014004 (2014).
  * Poulin and Preskill (2017) D. Poulin and J. Preskill (2017), Frontiers of Quantum Information Physics, KITP, URL http://online.kitp.ucsb.edu/online/qinfo-c17/poulin/.
  * Oppenheim (2018) J. Oppenheim, arXiv:1811.03116 [hep-th, physics:quant-ph] (2018).
  * Hall and Reginatto (2005) M. J. Hall and M. Reginatto, Physical Review A 72, 062109 (2005).
  * Poulin (2017) D. Poulin (2017), private communication (result announced in Poulin and Preskill (2017)).
  * Kent (1998) A. Kent, Physica Scripta T76, 78 (1998), ISSN 0031-8949, URL http://dx.doi.org/10.1238/Physica.Topical.076a00078.
  * Dowker and Kent (1995) F. Dowker and A. Kent, Physical Review Letters 75, 3038 (1995).
  * Oppenheim et al. (2020) J. Oppenheim, B. Soda, C. Sparaciari, and Z. Weller-Davies (2020), manuscript in preparation.
  * Belavkin (1990) V. P. Belavkin, Letters in Mathematical Physics p. 85–89 (1990).
  * Dalibard et al. (1992) J. Dalibard, Y. Castin, and K. Mølmer, Phys. Rev. Lett. 68 (1992).
  * Gardiner et al. (1992) C. W. Gardiner, A. S. Parkins, and P. Zoller, Phys. Rev. A 46 (1992).
  * Alicki and Lendi (1987) R. Alicki and K. Lendi, _Quantum Dynamical Semigroups and Applications_ (Springer, Berlin, 1987).
  * Breuer and Petruccione (2002) H. P. Breuer and F. Petruccione, _The Theory of Open Quantum Systems_ (Oxford University Press, 2002).
  * Lindblad (1976) G. Lindblad, Commun. Math. Phys. p. 48 (1976).
  * Gorini et al. (1976) V. Gorini, A. Kossakowski, and E. Sudarsahan, J. Math. Phys. p. 821 (1976).
  * Davies (1974) E. B. Davies, Commun. Math. Phys. p. 91–110 (1974).
  * Carmichael (1993) H. Carmichael, _An Open Systems Approach to Quantum Optics_ (Springer, Berlin, 1993).
  * Mølmer and Castin (1996) K. Mølmer and Y. Castin, Quantum and Semiclassical Optics: Journal of the European Optical Society Part B 8 (1996).
  * Mølmer et al. (1993) K. Mølmer, Y. Castin, and J. Dalibard, Journal of the Optical Society of America B 10, 524 (1993).
  * Wiseman (1996) H. Wiseman, Quantum and Semiclassical Optics: Journal of the European Optical Society Part B 8, 205 (1996).
  * Omnes (1992) R. Omnes, Reviews of Modern Physics 64, 339 (1992).
  * Griffiths (1993) R. B. Griffiths, Physical review letters 70, 2201 (1993).
  * Gell-Mann and Hartle (1993) M. Gell-Mann and J. B. Hartle, Physical Review D 47, 3345 (1993).
  * Diósi et al. (1995) L. Diósi, N. Gisin, J. Halliwell, and I. C. Percival, Physical review letters 74, 203 (1995).
  * Paz and Zurek (1993) J. P. Paz and W. H. Zurek, Physical Review D 48, 2728 (1993).
  * Ghirardi et al. (1986) G. C. Ghirardi, A. Rimini, and T. Weber, Physical Review D 34, 470 (1986), publisher: American Physical Society, URL https://link.aps.org/doi/10.1103/PhysRevD.34.470.
  * Pearle (1989) P. Pearle, Physical Review A 39, 2277 (1989), publisher: American Physical Society, URL https://link.aps.org/doi/10.1103/PhysRevA.39.2277.
  * Ghirardi et al. (1990) G. C. Ghirardi, P. Pearle, and A. Rimini, Physical Review A 42, 78 (1990), publisher: American Physical Society, URL https://link.aps.org/doi/10.1103/PhysRevA.42.78.
  * Bassi et al. (2013) A. Bassi, K. Lochan, S. Satin, T. P. Singh, and H. Ulbricht, Reviews of Modern Physics 85, 471 (2013).
  * Banks et al. (1984) T. Banks, M. E. Peskin, and L. Susskind, Nuclear Physics B 244, 125 (1984).
  * Oppenheim and Weller-Davies (2020) J. Oppenheim and Z. Weller-Davies (2020), manuscript in preperation.
  * Noschese et al. (2013) S. Noschese, L. Pasquini, and L. Reichel, Numer. Linear Algebra Appl. 20, 302 (2013).
  * Diosi (1987) L. Diosi, Physics letters A 120, 377 (1987).
  * Penrose (1996) R. Penrose, General relativity and gravitation 28, 581 (1996).

##  Appendix A The CQ master equation as a stochastic equation

In this appendix, we provide an analytical solution for the CQ master
equations considered in Sec. III. These equations are special cases of the
more general one shown in Eq. (III), which we reproduce here for simplifying
the exposition,

$\displaystyle\frac{\partial\,\varrho(z,t)}{\partial t}$
$\displaystyle=-i\left[H_{I}(z),\varrho(z,t)\right]+\frac{1}{\tau_{0}}\left(e^{\tau_{0}\,\left\\{{H_{C}(z)},{\,\cdot\,}\right\\}}-1\right)\varrho(z,t)$
$\displaystyle+\frac{1}{\tau}\sum_{\alpha\neq
0}\left(e^{\tau\,\left\\{{h^{\alpha}(z)},{\,\cdot\,}\right\\}}\,L_{\alpha}\,\varrho(z,t)\,L^{\dagger}_{\alpha}-\frac{1}{2}\left\\{{L^{\dagger}_{\alpha}L_{\alpha}},{\varrho(z,t)}\right\\}_{+}\right),$

where we recall that the interaction Hamiltonian $H_{I}(z)=\sum_{\alpha\neq
0}h^{\alpha}(z)L^{\dagger}_{\alpha}L_{\alpha}$ drives a stochastic evolution
of the classical and quantum degrees of freedom with a rate $\tau>0$, while
the classical Hamiltonian $H_{C}(z)$ governs the stochastic evolution of the
sole classical degrees of freedom with a different rate $\tau_{0}$.

###  A.1 Solution of a stochastic equation

We first present the solution of Eq. (III) in the special case in which no
classical-quantum coupling is present, and the whole master equation only
describes the evolution of the classical degrees of freedom. The classical
system is here described by the ensemble $\rho(q,p,t)$, and its evolution is
given by

$\frac{\partial\rho(q,p,t)}{\partial
t}=\frac{1}{\tau_{0}}\left(e^{\tau_{0}\left\\{{H_{C}(q,p)},{\cdot}\right\\}}-1\right)\rho(q,p,t),$
(84)

where the finite parameter $\tau_{0}>0$ can be understood as the rate of
jumping in phase space, due to the effect of the shift operator
$e^{\tau\left\\{{H_{C}(q,p)},{\cdot}\right\\}}$.

It is easy to show that the solution of Eq. (84) is given by

$\rho(q,p,t)=\sum_{k=0}^{\infty}P(k)\,\rho(k\tau_{0}),$ (85)

where $P(k)$ is the Poisson distribution which carries the time dependence of
the solution,

$P(k)=\left(\frac{t}{\tau_{0}}\right)^{k}\frac{e^{-\frac{t}{\tau_{0}}}}{k!}$
(86)

and the ensemble $\rho(k\tau_{0})$ is obtained from the initial ensemble
$\rho(q,p,t=0)$ by applying the operator $k$ times,

$\rho(k\tau_{0})=e^{k\tau_{0}\left\\{{H_{0}(q,p)},{\cdot}\right\\}}\rho(q,p,0).$
(87)

It is worth noting that each ensemble $\rho(k\tau_{0})$ is the solution of the
Liouville equation

$\frac{\partial\rho(q,p,t)}{\partial
t}=\left\\{{H_{C}(q,p)},{\rho(q,p,t)}\right\\},$ (88)

at time $t=k\tau_{0}$, so that we can interpret it as the state of the system
after $k$ jumps in phase space. Thus, the solution of Eq. (84) can be
understood as a mixture of different ensembles, where each of these ensemble
are the solution of Eq. (88) at discrete times $k\tau_{0}$, for
$k\in\mathbb{N}$, and the probability of making $k$ jumps is given by the
Poisson distribution whose mean value is $\frac{t}{\tau_{0}}$.