Decoherence, discord and the quantum master equation for cosmological perturbations
Abstract
We examine environmental decoherence of cosmological perturbations in order to study the quantumtoclassical transition and the impact of noise on entanglement during inflation. Given an explicit interaction between the system and environment, we derive a quantum master equation for the reduced density matrix of perturbations, drawing parallels with quantum Brownian motion, where we see the emergence of fluctuation and dissipation terms. Although the master equation is not in Lindblad form, we see how typical solutions exhibit positivity on superhorizon scales, leading to a physically meaningful density matrix. This allows us to write down a Langevin equation with stochastic noise for the classical trajectories which emerge from the quantum system on superhorizon scales. In particular, we find that environmental decoherence increases in strength as modes exit the horizon, with the growth driven essentially by white noise coming from local contribution to environmental correlations. Finally, we use our master equation to quantify the strength of quantum correlations as captured by discord. We show that environmental interactions have a tendency to decrease the size of the discord and that these effects are determined by the relative strength of the expansion rate and interaction rate of the environment. We interpret this in terms of the competing effects of particle creation versus environmental fluctuations, which tend to increase and decrease the discord respectively.
I Introduction
The inflationary paradigm G1 (); G2 () successfully describes the presentday uniformity of the Universe on large scales, as well as accounting for small fluctuations in the cosmic microwave background (CMB) P1 (); P2 (); P3 (); P4 (), leading to a nearlyscaleinvariant power spectrum of temperature anisotropies. According to the inflationary hypothesis, the origin of these inhomogeneities can be traced to primordial quantum fluctuations produced during the inflationary epoch, which become stretched and amplified across the sky by the quaside Sitter expansion. In this sense, primordial fluctuations sow the seeds for the largescale structure observed today. However, both largescale structure and even the CMB are treated in an essentially classical way, in that they are described by classical probabilities. A question of fundamental importance is therefore to understand how, and in what sense these fluctuations are rendered classical as they journey outside the horizon.
For free fields in (quasi)de Sitter there is a standard answer to this question. Since the fluctuations are isolated, they can be described by a pure state whose corresponding wavefunction
\Psi[\chi_{\bm{k}}]=\left(\frac{2\text{Re}\,\Omega_{{\bm{k}}}}{\pi}\right)^{1/% 2}\exp\left(\Omega_{{\bm{k}}}\chi_{\bm{k}}^{2}\right),  (1) 
becomes highly squeezed (in the momentum direction) as modes exit the horizon. As explained in KP (); Kiefer:2008ku (); Martin:2012pea () expectation values calculated in this squeezed state can be described as averages over a set of classical stochastic trajectories in analogy with a Gaussian wavepacket of a freeparticle. Nevertheless, the wave function is still a pure state with long range coherence in field space and the interpretation in terms of classical stochastic trajectories is seen to be special to a free theory. It is interesting to go beyond the free theory and study the quantumtoclassical transition when there are non trivial interactions.
Indeed, in reality, fluctuations are not truly isolated and will interact with other fields (at the very least they must interact with metric fluctuations) playing the role of an “environment”. Therefore, one necessarily has to understand the emergence of classicality in open systems, and in what sense a classical stochastic description emerges. This is the subject of open quantum systems BP (), which provides a formalism to analyse the reduced density matrix of the particular system (and its evolution via a quantum master equation) by tracing out the environmental degrees of freedom. In particular, this formalism allows one to study the process of environmental decoherence, which at the simplest level, leads to the decay of quantum interference in the system. It also allows one to understand how, starting from an initially quantum state, classical stochastic dynamics emerges due to interaction with an environment – this is of course more challenging than the situation in the free theory where the classical description is rather trivial.
Steps towards understanding the evolution of the density matrix during inflation in the presence of an external environment have already been made in refs. KP (); KPS2000 (); Kiefer:2006je () where a phenomenologically motivated ansatz of the form
\displaystyle\bra{\chi_{{\bm{k}}}}\hat{\rho}\ket{\tilde{\chi}_{{\bm{k}}}}% \equiv\rho(\chi_{{\bm{k}}},\tilde{\chi}_{{\bm{k}}})=\left(\frac{2\text{Re}\,% \Omega_{{\bm{k}}}}{\pi}\right)  
\displaystyle\times\exp\big{[}\Omega_{{\bm{k}}}(\eta)\chi_{\bm{k}}^{2}% \Omega^{*}_{{\bm{k}}}(\eta)\tilde{\chi}_{\bm{k}}^{2}\frac{\xi}{2}\left\chi% \tilde{\chi}\right^{2}\big{]}\ ,  (2) 
was considered, where \Omega_{\bm{k}} are simply the freesolutions associated to eq. (1). This ansatz consists essentially of adding by hand an additional factor \exp[\xi/2\left\chi\tilde{\chi}^{\prime}\right^{2}] to the free density matrix associated to (1), where \xi is assumed to be a constant determined by the details of the environment. Plainly eq. (I) is heuristic and warrants further investigation in a concrete model. For instance, \Omega_{\bm{k}}(\eta) will obviously receive corrections of order the coupling strength with the environment, thereby incorporating new features into the power spectrum (which is proportional to [\text{Re}\,\Omega_{\bm{k}}]^{1}). Similarly, one would expect, and as we shall see, that \xi has some time and momentum dependence \xi_{\bm{k}}(\eta). What is needed is to consider particular interactions in order to formulate a more comprehensive analysis of the evolution of the density matrix allowing one to build on the work of, e.g., KPS1998 (); BLM () which sketched the implications of particular environmentsystem couplings. In particular one would like to examine the timeevolution of \Omega_{{\bm{k}}}(\eta) and particularly \xi_{\bm{k}}(\eta) in a specific model coupling fluctuations to another field.
Under favourable conditions, one can then obtain a quantum master equation for the timeevolution of \rho into which the ansatz (I) can be substituted to determine the evolution equations for \Omega_{{\bm{k}}}(\eta) and \xi_{\bm{k}}(\eta).
In order to have a model of environmental interactions which is amenable to explicit calculation, it is very useful to have one that maintains the description in terms of individual momentum modes of the perturbations \chi_{\bm{k}}. This leads us to consider the model developed by Boyanovksy Boy1 (); Boyanovsky:2015jen () where the interaction is between the “system” field \phi and another “environment field” \varphi coupled via
S_{\text{int}}=\int d^{4}x\,\sqrt{g}\lambda\phi\varphi^{2}=\int d^{4}x\,a(% \eta)\lambda\,\chi\psi^{2},  (3) 
where \chi=a\phi and \psi=a\varphi are the usual conformally rescaled fields. Note that the special form of the coupling here, linear in \chi (i.e. \phi), means that the interaction does not mix the different momentum modes \chi_{\bm{k}} up to and including \mathcal{O}(\lambda^{2}), and so at lowest order we can still consider the problem modebymode providing a great computational simplification. In Boyanovsky:2015jen (), Boyanovksy was focussed more on using the master equation to derive corrections to the power spectrum rather than the details of decoherence and the quantumtoclassical transition. In this paper, we shall develop his analysis and use the interaction (3) as a toymodel in which to address the evolution of decoherence and the onset of classical stochastic behaviour as modes exit the horizon.
There are different formalisms for modelling the effects of interactions with an environment involving apparently rather different approximations. On the one hand, via a master equation for the reduced density matrix, and on the other, via the influence action. The latter approach can be shown, under suitable circumstances, to be equivalent to the master equation approach in flat space problems Boy1 (), but it is by no means clear whether this will be true in the inflationary context. Although it seems the key observable of the fluctuations – the power spectrum – can be calculated in either approach with agreement Boyanovsky:2015jen (). Some previous works on decoherence in the inflationary context have concentrated on the influence action formalism Sakagami:1987mp (); Lombardo:2005iz (); however, the master equation route offers the advantage of allowing one to follow the quantum state explicitly as expansion and decoherence proceeds and provides a way to calculate the growth of entropy and assess the degree to which the perturbations are classical. Other important work on cosmological perturbations, decoherence and the quantumtoclassical transition includes Albrecht:1992kf (); Polarski:1995jg (); Lesgourgues:1996jc (); Kiefer:1998qe (); Kiefer:1998pb (); Campo:2008ij (); Kiefer:2006je (); Martineau:2006ki (); Burgess:2006jn (); Kiefer:2008ku (); Martin:2012pea (); Burgess:2014eoa (); Burgess:2015ajz (); Matacz:1996gk (); EN ().
Key results
1. By deriving and solving evolution equations for \xi_{{\bm{k}}}(\eta) and \Omega_{{\bm{k}}}(\eta) we are able to show that \Omega_{{\bm{k}}}(\eta) does indeed receive additional \lambdadependent corrections. For instance we find that
\displaystyle\text{Re}\left[\Omega_{{\bm{k}}}(\eta)\right]  
\displaystyle=k^{3}\eta^{2}\exp\left[\frac{\lambda^{2}}{12\pi^{2}H^{2}}\Big{(}% \log^{2}\frac{\eta}{\eta_{0}}\frac{4}{3}\log\frac{\eta}{\eta_{0}}\Big{)}% \right]+O\left(\eta^{3}\right).  (4) 
Therefore by direct computation of the density matrix we are able to reproduce, via an alternate calculational route, the resummed corrections to the power spectrum (which is proportional to 1/\text{Re}[\Omega_{{\bm{k}}}]) found by Boyanovksy Boyanovsky:2015jen (). On superhorizon scales these corrections are dominated by nonlocal correlations of the “environment field” \varphi, corresponding to longtime correlations, i.e. memory effects, as noted in Boyanovsky:2015jen ().
2. We find that decoherence, as characterised by \xi_{\bm{k}}(\eta), exhibits timedependence, oscillating on subhorizon scales (k\eta\gg 1) and growing monotonically outside the horizon (k\eta\ll 1) – a manifestation of the fact that by virtue of the operator being relevant the factor a(\eta) in the RHS of (3) leads essentially to an interaction strength which grows in time. By contrast to \Omega_{\bm{k}}(\eta), we find that the superhorizon evolution of the decoherence parameter, \xi_{\bm{k}}(\eta) is driven by the local environment correlators, in other words, by white noise relating to delta function selfcorrelations of the environment.
3. We find that on subhorizon scales the density matrix violates positivity, a symptom of the fact that the master equation is not in Lindblad – the manifestly positivity preserving – form. This shows that on subhorizon scales the assumptions that go into deriving the master equation cannot be completely consistent. However, positivity is recovered on superhorizon scales when the state becomes very squeezed. We also see the same behaviour at the level of the Wigner function, which satisfies a FokkerPlanck equation exhibiting unphysical “negativediffusion” on subhorizon scales which becomes unimportant as modes exit the horizon. In this sense we provide a complete picture of how decoherence and positivity emerge as modes exit the horizon leading to the interpretation of the density matrix as a classical Probability Density Function (PDF) after it becomes both diagonal and positive.
4. It is a wellknown fact that in realistic systems, entanglement between spacelike separated regions is extremely sensitive to noise from the local environment BP (), it is therefore important to ask how vulnerable entanglement produced during inflation is to decoherence. The final part of this paper uses the quantum master equation and its solutions to study the effects of environmental noise on entanglement produced during inflation. The principle measure we use is quantum discord OZ (); HV (). Our analysis reveals that the environment does indeed weaken the strength of discord with the corrections being of order \lambda/H, a ratio characterising the relative strength of entangled particle pair production and the interaction rate with the environment.
The structure of the remainder of the paper is as follows. Beginning in sec. II, we discuss open quantum systems in general terms using the paradigm of quantum Brownian motion to describe various concepts, which is particularly fitting as the modebymode master equation of cosmological perturbations assumes the same schematic form. Then, in sec. III, we introduce the simplified model used to study an open system of primordial perturbations and review the derivation of the corresponding master equation outlined in Boyanovsky:2015jen ().
Once this groundwork is laid, in sec. IV we provide a detailed analysis of the master equation, drawing comparisons with quantum Brownian motion BP (); CL (); H (), where we identify the usual fluctuationdissipation terms, in addition to the unitary coherent evolution associated to the isolated squeezed state. We find that, although the master equation is not in Lindblad form, positivity of the density matrix is established on superhorizon scales, providing a natural emergence of classicality at late times. We also write the master equation in terms of the Wigner function which we then use to construct a Langevintype equation. This provides an equivalent description of the model – in the sense of expectation values – in terms of a stochastic theory of classical trajectories. In this sense we provide a firstprinciples derivation in the context of the model (3), of a classical stochastic equation which reproduces the dynamics of the quantum master equation on superhorizon scales. It should be noted that this is different from the stochastic inflation paradigm Star1 (); SY (), where the noise is described by subHubble modes in contrast to the present setup, where the noise arises from another field, conformally coupled so that it can mimic the subHubble environment in a computational tractable model Boyanovsky:2015jen ().
Finally, in sec. V we shall address the question of whether quantum correlations (generated by entangled pair creation) are robust to environmental decoherence. Recently, a novel approach to entanglement in inflationary cosmology has been explored by Lim Lim () and Vennin and Martin VM (). The idea is to harness the power of correlations in the system of perturbations using quantum discord OZ (); HV () which relies on the difference between classical and quantum correlations to quantify the extent to which a particular system exhibits quantum features.
It is natural to expect entanglement between different regions of a system to be subdued due to entanglement with a local environment. Indeed it has already been noted in other contexts that an external environment weakens discord BA (). The detailed analysis of decoherence provided by our quantum master equation and its solutions allows us to examine precisely this question. We find that environmental decoherence has a tendency to reduce the size of the discord and offer an interpretation of this in terms of the competing effects of particle pair creation in de Sitter (which tends to increase the strength correlations) and random production of particles due to the environmental scatterings, which weakens the strength of correlations. In section VI we offer our conclusions and some suggestions for future work.
II The quantum master equation
Since our central goal to is understand how a quantum state can end up classical, the nature of the quantumtoclassical transition deserves some comment. We take a pragmatic view here: classical mechanics is manifestly an emergent phenomenon of quantum systems in certain circumstances and the two key issues are: (i) what are those circumstances? And, (ii) what are the phenomenological rules that describe how the classical state emerges from an underlying quantum state?
We view the emergence of classical mechanics in the same way one views any effective theory of some underlying microscopic system. The resulting description will always be phenomenological and only appropriate to a coarsegrained perspective. By “coarsegrained” we mean that the effective theory is couched in terms of observables that have some realistic finite resolution scales, e.g. the dynamics of fluctuations on superhorizon scales in the cosmological setting. For observers probing these scales, an evolving quantum state can be equivalent to an ensemble of classical stochastic trajectories in the sense that quantum expectation values are captured by stochastic averaging. The final step is to postulate that an individual classical stochastic trajectory is real for the coarsegrained observer.
II.1 Analogy with quantum Brownian motion
A useful toy model which has many parallels with the present investigation is the muchdiscussed subject of quantum Brownian motion H (); BP (); CL (); JZ (). The system consists of a particle, which, in the simplest setup, moves in one dimension, interacting with an environment. The environment can be modelled in many different ways leading to the same universal behaviour for the particle. The first level of coarse graining is to take a perspective – defined by a set of observables – which act on the particle alone, e.g. its position and momentum. In these circumstances, one can trace out the environment and work with the reduced density matrix of the particle, \rho=\text{Tr}_{\mathcal{E}}\Psi\rangle\langle\Psi. Here, \ket{\Psi} is the state vector of the joint environmentparticle system whose evolution is unitary. In general, \rho does not satisfy an autonomous dynamical equation; however, there are circumstances for which one can derive an approximate autonomous equation for \rho known as the “master equation”. In going from the unitary evolution of the total state, to nonunitary evolution encoded in the quantum master equation, there are essentially three stages of approximation which must be made and justified:

The first is the Born approximation, in which one replaces the total (pure) state of the system \Psi\rangle by a product of density matrices for the particle and environment \rho_{S}\otimes\rho_{\mathcal{E}}. The idea is that the environment is so large compared with the system (the particle) that in every instant the system is interacting with a fresh bit of the environment and the overall environment is barely affected by the interaction with the much smaller system. Any entanglement that previously builds up between the particle and the environment is encoded in the von Neumann entropy of the reduced density matrix of the particle which increases over time.

The second approximation is the Markov approximation that allows one to replace terms that depend on the history of the system by the instantaneous state \rho. This approximation is justified when the dynamics of \mathcal{E} is suitably fast compared with the slow dynamics of the macroscopic particle.

Finally, the resulting equation for \rho has terms that oscillate rapidly on the fast environmental time scale. Temporal coarse graining amounts to removing these terms by hand in a form of rotating wave approximation BP ().
If all these steps are followed, a master equation for \rho can emerge with solutions that are physically consistent, i.e. which preserve the positivity of the spectrum of \rho. The latter is needed for a consistent density matrix which must have a spectrum \{p_{n}\} for which p_{n}\geq 0 and \sum_{n}p_{n}=1 (the normalization condition). The necessary and sufficient condition for this is that the master equation can be written in Lindblad form, that is as
\begin{split}\displaystyle\dot{\rho}=\frac{1}{i\hbar}[H,\rho]+\sum_{i}\alpha_{% i}\Big{(}2A_{i}\rho A_{i}^{\dagger}A^{\dagger}_{i}A_{i}\rho\rho A^{\dagger}_% {i}A_{i}\Big{)}\ ,\end{split}  (5) 
with \alpha_{i}\in\mathbb{R}>0. It is the second term here that encodes the interaction with the environment. For quantum Brownian motion, the form of the master equation in the position basis \rho(x,x^{\prime})=\langle x\rhox^{\prime}\rangle, is the CaldeiraLeggett master equation CL ().
\begin{split}\displaystyle\dot{\rho}&\displaystyle=\overbracket{\frac{1}{i% \hbar}[H_{\text{eff}},\rho]}^{\text{unitary evol.}}\overbracket{\gamma(xx^{% \prime})\Big{(}\frac{\partial}{\partial x}\frac{\partial}{\partial x^{\prime}% }\Big{)}\rho}^{\text{dissipation}}\\ &\displaystyle\qquad\qquad\underbracket{\frac{2m\gamma k_{B}T}{\hbar^{2}}(xx% ^{\prime})^{2}\rho}_{\text{decoherence}}\ .\end{split}  (6) 
Here, H_{\text{eff}} is the Hamiltonian of the particle with some additional additive terms that arise from the coupling of the particle and \mathcal{E}. In fact, strictly speaking, this master equation is not in Lindblad form. To render it so, one must take
\begin{split}\displaystyle A=x+i\beta p\ ,\qquad\alpha=\frac{2m\gamma k_{B}T}{% \hbar^{2}}\ ,\qquad\beta=\frac{\gamma}{\alpha}\ ,\end{split}  (7) 
in the Lindblad form (5), and this gives rise to extra terms quadratic in momentum. However, in realistic situations where T is not too small, the departure from positivity is small at the level of solutions, so the lack of terms needed for positivity in (6) can be neglected. The characteristic features of the master equation are the dissipation (or relaxation) and fluctuation (decoherence) terms.
Now that we have an explicit example of an open quantum system, we can ask the allimportant question: what makes a state of this system “classical”? To start with the word “classical” is attached to quantum states in different ways. One answer is that a classical state is one which is highly localised in position and momentum, like a coherent state (i.e. a state which provides the minimal level of uncertainty allowed by the Heisenberg relation). If the spread \delta x is much smaller than the scale over which the potential energy V(x) varies then Ehrenfest’s theorem guarantees that the quantum state follows the classical trajectory. Unfortunately, this reasoning, although often invoked in textbooks, is completely unrealistic. Even in free systems, wave functions spread over time. This spreading is slow for a macroscopic particle but in any chaotic system this spreading becomes remarkably fast, order t_{c}\log(I/\hbar), where t_{c} is the characteristic (classical) chaos time scale and I is the characteristic (macroscopic) action scale of the system. In summary, we cannot rely on Ehrenfest’s theorem to describe the quantumtoclassical transition (see the discussion in Zurek:1998mj ()).
The second use of the word “classical” is in “semiclassical”, or WKB states. These kinds of state are pure and are smeared along classical trajectories of the system in phase space and so are not the kind of state that can describe the quantumtoclassical transition Berry ().
Finally, we have the notion of a classical state that is meaningful for the quantumtoclassical transition. These states are very spread out in both position and momentum. This seems paradoxical: how can such states be classical? It will emerge that the answer is that these states give quantum expectation values that can be captured by a classical stochastic process.
A state which is very spread out in position and momentum has coherence lengths \ell_{x} and \ell_{p} (the scales over which the offdiagonal elements of the density matrix fall off in position and momentum space, respectively) which are very small in the sense that
\begin{split}\displaystyle\ell_{x}\ell_{p}\ll\hbar\ .\end{split}  (8) 
Such states must necessarily be mixed states rather than pure and can be given a useful phase space representation by using, for example, the Wigner function. We choose the latter for its ubiquity even though there is nothing that really singles it out from other alternatives (e.g. the Husimi function):
\begin{split}\displaystyle W(x,p)=\frac{1}{\pi}\int_{\infty}^{\infty}dx^{% \prime}\,\rho\Big{(}x\frac{x^{\prime}}{2},x+\frac{x^{\prime}}{2}\Big{)}e^{2% ipx^{\prime}/\hbar}\ .\end{split}  (9) 
In general, there is nothing to guarantee the positivity of W(x,p) – unlike the Husimi function which is manifestly positive – but for a “classical state” spread out in phase space over an area \gg\hbar, W(x,p) becomes positive.
The master equation (6) can be written in terms of the Wigner function as
\begin{split}\displaystyle\dot{W}&\displaystyle=\overbracket{\{H_{\text{eff}},% W\}}^{\text{Poisson bra.}}+\overbracket{2\gamma\partial_{p}(pW)}^{\text{% dissipation}}+\overbracket{2m\gamma k_{B}T\partial_{p}^{2}W}^{\text{% decoherence}}\\ &\displaystyle+\sum_{n\geq 1}\underbracket{\frac{(1)^{n}\hbar^{2n}}{2^{2n}(2n% +1)!}\partial_{x}^{2n+1}V\partial_{p}^{2n+1}W}_{\text{quantum corrections}}\ .% \end{split}  (10) 
Note that the quantum corrections above arise from nonlinear terms in the dynamics, that is, on the derivatives of the potential V(x) of order three and greater.
We can now consider the effect of the various terms in the master equation. Firstly, with only the Poisson bracket term on the righthand side, Liouville’s theorem would ensure that the area of the support of W is preserved. Consider an initial state with minimal support allowed by the uncertainty principle A\sim\hbar – a coherent state. The nonlinear interactions in the Poisson bracket term have the tendency to make W spread out in phase space but at the same time vary on small scales so that the area remains fixed A\sim\hbar. Generally the area is squeezed in some directions and expands in other directions in such a way that the perimeter grows rapidly. Now consider the effect of the quantum corrections in (10). These generally destroy the probability density interpretation of W and encode non trivial quantum interference. Specifically, these terms will lead to a violation of the positivity of the Wigner function as it “interferes with itself”.
However, the decoherence terms coming from the interactions of the particle with the environment can rescue a probabilistic interpretation. In particular the decoherence term, causes the Wigner function to spread out but also become smooth on small scales. The positivity violations are alleviated and the area of support is driven to A\gg\hbar. This spreading and smoothing goes hand in hand with coherence lengths that are driven to the regime (8). This lack of purity can also be quantified by the entanglement entropy. Note that the smoothing out caused by decoherence also implies that the different phase space representations of a state, Wigner, Husimi, etc., essentially become indistinguishable.
The fact that the Wigner function becomes smoothed out means that the derivatives \partial_{p}^{n}W are small and so the quantum corrections in (10) are suppressed Zurek:1998mj (). Effectively one can describe the dynamics of the Wigner function by the much simpler equation
\begin{split}\displaystyle\dot{W}\approx\{H_{\text{eff}},W\}+2\gamma\partial_{% p}(pW)+2m\gamma k_{B}T\partial_{p}^{2}W\ .\end{split}  (11) 
Now comes the key point for the quantumtoclassical transition. This simplified master equation has precisely the form of a FokkerPlanck equation where the decoherence term corresponds to diffusion and the dissipation term to friction. In other words, the Wigner function can now be given the interpretation of a classical probability density function for a set of classical stochastic trajectories described by the Langevin equations
\begin{split}\displaystyle\dot{x}=\{x,H_{\text{eff}}\}\ ,\qquad\dot{p}=\{p,H_{% \text{eff}}\}2\gamma p+\sigma\ ,\end{split}  (12) 
where \sigma(t) is a Gaussian white noise term with characteristics
\begin{split}\displaystyle{\mathbb{E}}(\sigma(t))=0\ ,\qquad{\mathbb{E}}(% \sigma(t)\sigma(t^{\prime}))=4m\gamma k_{B}T\delta(tt^{\prime})\ .\end{split}  (13) 
The final step in the quantumtoclassical transition is to interpret the Langevin equation as describing classical trajectories that have emerged out of the underlying quantum dynamics at the level of the coarsegrained perspective (after tracing over the environment). Note that the equivalence to a classical stochastic process means that quantum expectation values can be captured by averages over classical stochastic trajectories:
\begin{split}\displaystyle\text{Tr}\big{(}\rho(t){\cal O}(x,p)\big{)}\approx% \sum_{\{x(t^{\prime}),p(t^{\prime})\}}{\cal O}(x(t),p(t))\ ,\end{split}  (14) 
where the average on the righthand side is over the stochastic trajectories up to time t with random initial conditions that correspond to the initial Wigner function.
Let us summarize the phenomenological approach to the quantumtoclassical transition:
We shall now see how this same chain of steps emerges during inflation, with environmental interactions leading to diffusion/decoherence terms in the master equation that are then interpreted as stochastic noise acting on the emergent classical trajectories. The reduced dynamics is given by tracing out this environment with an additional level of coarse graining provided by restricting to superHubble modes – the appropriate scales for analysing the CMB.
It is worth pointing out at this stage, that what we have described above is a way to think about the quantumtoclassical transition and how classical trajectories can be thought to have emerged from the quantum system in an interacting system. However, the simplest cosmological model is free, i.e. V(x) is quadratic and there is no environment. In that case, the master equation for the Wigner function obviously takes the form of a FokkerPlanck equation immediately and therefore is clearly equivalent to a classical stochastic process where the trajectories in this case satisfy the classical equations of motion (there is no noise in this case because there is no diffusion term) but with initial conditions that are defined randomly with respect to a PDF that is the initial Wigner function (this is shown quite explicitly in KP (); Kiefer:2008ku (); Martin:2012pea ()). Note that this interpretation can be made even though the quantum state remains pure. It is really a moot point as to whether this as a valid description of the quantumtoclassical transition because in any realistic system, including the cosmological one, we expect that there will be interactions that cannot be ignored.
III Master equation for the reduced density matrix
In order to study an open system of the scalar cosmological perturbations, we introduce a simplified model considered by Boyanovsky Boyanovsky:2015jen () (whose notation we follow) by coupling the scalar perturbations \phi to an external scalar field \varphi via the action
\displaystyle S  \displaystyle=\int d^{3}x\,dt\,a^{3}\,\Big{\{}\frac{1}{2}\dot{\phi}^{2}\frac{% (\nabla\phi)^{2}}{2a^{2}}\frac{1}{2}\left(M_{\phi}^{2}+\xi_{\phi}\mathcal{R}% \right)\phi^{2}  
\displaystyle+\frac{1}{2}\dot{\varphi}^{2}\frac{(\nabla\varphi)^{2}}{2a^{2}}% \frac{1}{2}\left(M_{\varphi}^{2}+\xi_{\varphi}\mathcal{R}\right)\varphi^{2}% \lambda\,\phi\,:\varphi^{2}:\Big{\}},  (15) 
where
\mathcal{R}=6\left(\frac{\ddot{a}}{a}+\frac{\dot{a}^{2}}{a^{2}}\right),  (16) 
is the Ricci scalar and \xi=0,1/6 correspond to the cases of minimal and conformal coupling, respectively. The interaction is normalordered so that
:\varphi^{2}:=\varphi^{2}\braket{\varphi^{2}},  (17) 
where \langle\cdots\rangle denotes a QFT correlator. We shall work in conformally rescaled fields \phi=\chi/a and \varphi=\psi/a and introduce conformal time \eta=1/aH, in which the action becomes
\begin{split}\displaystyle S&\displaystyle=\int d^{3}x\,d\eta\,\Big{\{}\frac{1% }{2}\big{[}\chi^{\prime 2}(\nabla\chi)^{2}\mathcal{M}^{2}_{\chi}(\eta)\chi^{% 2}\big{]}\\ &\displaystyle+\frac{1}{2}\big{[}\psi^{\prime 2}(\nabla\psi)^{2}\mathcal{M}^% {2}_{\psi}(\eta)\psi^{2}\big{]}+\frac{\lambda}{H\eta}\,\chi:\psi^{2}:\Big{\}},% \end{split}  (18) 
where
\mathcal{M}_{\chi,\psi}(\eta)=\Big{[}\frac{\mathcal{M}_{\chi,\psi}^{2}}{H^{2}}% +12\Big{(}\xi_{\chi,\psi}\frac{1}{6}\Big{)}\Big{]}\frac{1}{\eta^{2}},  (19) 
provides an effective, timedependent mass.
In order to provide a simplified model of the cosmological perturbations, \chi plays the rôle of the MukhanovSasaki variable and so is a minimally coupled (\xi_{\chi}=0) massless field. In fact, the bare mass will be nonvanishing in order to absorb a UV divergence. On the other hand, for the environment field, choices can be made. We follow Boyanovsky:2015jen () and take \psi to be a conformally coupled field (\xi_{\psi}=1/6). The intuition for this choice is that the field \psi can either be viewed as a field in its own right or as a proxy for the high frequency components of \chi whose wavelength remain inside the horizon during inflation. This mimics a nonlinear selfinteraction that couples modes of different momentum. In this sense, the simplified model attempts to model an IRUV split of the degrees of freedom of a single (interacting) field.
It is convenient – at least to start with – to work in the interaction picture where the observables evolve in the free Heisenberg picture and the density matrix evolves according to the von Neumann equation with the interaction Hamiltonian
\frac{d\rho_{I}(\eta)}{d\eta}=i\left[H_{I}(\eta),\rho_{I}(\eta)\right].  (20) 
where
H_{I}(\eta)=\frac{\lambda}{H\eta}\int d^{3}x\,\chi({\bm{x}},\eta):\psi^{2}({% \bm{x}},\eta):,  (21) 
and where \chi({\bm{x}},\eta) and \psi({\bm{x}},\eta) are fields in the free Heisenberg representation. It is worth pointing out at this stage that since the interaction (3) is time dependent, naïve perturbation theory is only valid on timescales satisfying \lambda/H\lesssim\eta.
The von Neumann equation (20) has solutions
\rho_{I}(\eta)=\rho_{I}(\eta_{0})i\int^{\eta}_{\eta_{0}}d\eta^{\prime}\left[H% _{I}(\eta^{\prime}),\rho_{I}(\eta^{\prime})\right].  (22) 
Feeding this back into (20) gives the leading terms in the perturbative expansion describing the evolution of \rho_{I}:
\begin{split}\displaystyle\frac{d\rho_{I}(\eta)}{d\eta}&\displaystyle=i\left[% H_{I}(\eta),\rho_{I}(\eta_{0})\right]\\ &\displaystyle\int^{\eta}_{\eta_{0}}d\eta^{\prime}\left[H_{I}(\eta),\left[H_{% I}(\eta^{\prime}),\rho_{I}(\eta^{\prime})\right]\right].\end{split}  (23) 
In order to derive an autonomous equation for the evolution of the density matrix of the field \chi – the master equation – some familiar approximations are necessary. Firstly, the Born approximation assumes that the total density matrix can effectively be factorised as
\rho_{I}(\eta)=\rho_{I\chi}(\eta)\otimes\rho_{I\psi}(\eta_{0}).  (24) 
As described in sec. II, this is justified when the environment is large enough such that in each time interval, the system effectively interacts (weakly) with a fresh part of the environment. Essentially, the environment is unaffected by the interaction, while the system’s entanglement with the environment is encoded in the reduced density matrix \rho_{I\chi}. The reason why factorisation may be justified in the present context is due to the fact that a single mode of the field \chi couples to a pair of modes of the environmental scalar field and in a de Sitter background there are no particle thresholds. Brownian motion, described in sec. II, provides a familiar example which exploits these same approximations. In the present context, this approximation is described in greater detail in BP (); Boy1 (); Boyanovsky:2015jen ().
Combining the approximations (23) and (24) leads to the following equation Boyanovsky:2015jen () for the reduced interactionpicture density matrix \rho_{I\chi}(\eta) (in the following we will now denote the reduced density matrix of \chi: \rho_{I}\equiv\rho_{I\chi}):
\displaystyle\frac{d\rho_{I}(\eta)}{d\eta}  \displaystyle=\frac{\lambda^{2}}{H^{2}\eta}\int^{\eta}_{\eta_{0}}\frac{d\eta^% {\prime}}{\eta^{\prime}}\int d^{3}x\int d^{3}y\Big{\{}\chi({\bm{x}},\eta)\,% \chi({\bm{y}},\eta^{\prime})\,\rho_{I}(\eta^{\prime})\,G^{>}({\bm{x}}{\bm{y}}% ,\eta,\eta^{\prime})+\rho_{I}(\eta^{\prime})\,\chi({\bm{y}},\eta^{\prime})\,% \chi({\bm{x}},\eta)\,\,G^{<}({\bm{x}}{\bm{y}},\eta,\eta^{\prime})  
\displaystyle\,\chi({\bm{x}},\eta^{\prime})\,\rho_{I}(\eta^{\prime})\,\chi({% \bm{y}},\eta)\,\,G^{<}({\bm{x}}{\bm{y}},\eta,\eta^{\prime})\,\chi({\bm{y}},% \eta)\rho_{I}(\eta^{\prime})\chi({\bm{x}},\eta^{\prime})\,\,G^{>}({\bm{x}}{% \bm{y}},\eta,\eta^{\prime})\Big{\}},  (25) 
where
\begin{split}\displaystyle G^{>}({\bm{x}}{\bm{y}},\eta,\eta^{\prime})&% \displaystyle=\text{Tr}\left[:\psi^{2}({\bm{x}},\eta)::\psi^{2}({\bm{y}},\eta)% :\rho_{I\psi}(\eta_{0})\right],\\ \displaystyle G^{<}({\bm{x}}{\bm{y}},\eta,\eta^{\prime})&\displaystyle=\text{% Tr}\left[:\psi^{2}({\bm{y}},\eta)::\psi^{2}({\bm{x}},\eta):\rho_{I\psi}(\eta_{% 0})\right],\end{split}  (26) 
are environmental correlators.
Our next step is to make the Markov approximation by taking \rho_{I}(\eta^{\prime})\rightarrow\rho_{I}(\eta) in the integral expression in eq. (III). This is consistent with the weak coupling limit, as can be seen with an integration by parts of (III), noting that d\rho_{I}/d\eta\propto\lambda^{2}/H^{2}\ll 1 and neglecting the resulting term of \mathcal{O}(\lambda^{4}). Going to momentum space, we write
\displaystyle G^{\lessgtr}({\bm{x}}{\bm{y}},\eta,\eta^{\prime})  \displaystyle=\frac{1}{V}\sum_{{\bm{k}}}\mathcal{K}^{\lessgtr}(k,\eta,\eta^{% \prime})e^{i{\bm{k}}\cdot({\bm{x}}{\bm{y}})},  (27) 
where V is the quantisation volume, an IR regulator.
In order to proceed, we must compute the objects K^{\lessgtr}(p,\eta,\eta^{\prime}). As shown in Boyanovsky:2015jen (), after specialising to a BunchDavies (BD) vacuum BD () for the initial conditions of the field \psi, \rho_{I\psi}(\eta_{0})=\tensor{[}]{\ket{0}}{{}_{\psi}}\!\!\tensor{[}_{\psi}]{% \bra{0}}{}, we find
\displaystyle\mathcal{K}^{>}(q,\eta,\eta^{\prime})\equiv K(q,\eta,\eta^{\prime% })=2\int\frac{d^{3}k}{(2\pi)^{3}}v(k,\eta)v^{*}(k,\eta^{\prime})v(p,\eta)v^{*}% (p,\eta^{\prime})  (28) 
and {\cal K}^{<}(q,\eta,\eta^{\prime})=K^{*}(q,\eta,\eta^{\prime}), where q={\bm{k}}+{\bm{p}} and v(k,\eta), etc., are the mode functions associated to the environment field \psi. They satisfy
v^{\prime\prime}+\Big{[}k^{2}\frac{1}{\eta^{2}}\Big{(}\nu_{\psi}^{2}\frac{1}% {4}\Big{)}\Big{]}v=0,  (29) 
with \nu^{2}_{\psi}=9/4(M_{\psi}^{2}/H^{2}+12\xi_{\psi}). The steps above lead to a master equation
\begin{split}\displaystyle\frac{d\rho_{I}(\eta)}{d\eta}=\frac{\lambda}{H^{2}% \eta}\int^{\eta}_{\eta_{0}}\frac{d\eta^{\prime}}{\eta^{\prime}}\sum_{{\bm{k}}}% &\displaystyle\Big{\{}\chi_{{\bm{k}}}(\eta)\,\chi_{{\bm{k}}}(\eta^{\prime})% \rho_{I}(\eta)K(k,\eta,\eta^{\prime})+\rho_{I}(\eta)\,\chi_{{\bm{k}}}(\eta^{% \prime})\,\chi_{{\bm{k}}}(\eta)K^{*}(k,\eta,\eta^{\prime})\\ &\displaystyle\chi_{{\bm{k}}}(\eta)\rho_{I}(\eta)\,\chi_{{\bm{k}}}(\eta^{% \prime})\,K^{*}(k,\eta,\eta^{\prime})\chi_{{\bm{k}}}(\eta)\rho_{I}(\eta)\,% \chi_{{\bm{k}}}(\eta)\,K^{*}(k,\eta,\eta^{\prime})\Big{\}},\end{split}  (30) 
where \chi_{\bm{k}}(\eta) are the Fourier modes of the field satisfying \chi^{\dagger}_{{\bm{k}}}=\chi_{{\bm{k}}}.
Since \psi is massless and conformally coupled, with BD vacuum initial conditions, it has simple mode functions
v(k,\eta)=\frac{e^{ik\eta}}{\sqrt{2k}}.  (31) 
In this case,
K(p,\eta,\eta^{\prime})=\frac{i}{8\pi^{2}}e^{ip(\eta\eta^{\prime})}\mathcal% {P}\left[\frac{1}{\eta\eta^{\prime}}\right]+\frac{1}{8\pi}\delta(\eta\eta^{% \prime}),  (32) 
where \mathcal{P} stands for the principal part, which can be written as
\displaystyle\mathcal{P}\left[\frac{1}{\eta\eta^{\prime}}\right]=\frac{1}{2}% \frac{d}{d\eta}\log\left[\frac{(\eta\eta^{\prime})^{2}+\varepsilon^{2}}{(% \bar{\eta})^{2}}\right].  (33) 
Here \varepsilon\to 0 is a UV regulator and \bar{\eta} is an arbitrary scale to make the argument of the logarithm dimensionless, and, as we shall see shortly, it also acts as a renormalisation scale. Inserting eqs. (32) and (33) into eq. (30) leads to terms of the form
\displaystyle\int\frac{d\eta^{\prime}}{\eta^{\prime}}\chi_{{\bm{k}}}(\eta^{% \prime})K[q;\eta,\eta^{\prime}]  
\displaystyle\equiv\frac{1}{2}\int^{\eta}_{\eta_{0}}\!\!d\eta^{\prime}\,\,% \chi_{{\bm{k}}}\frac{e^{ik(\eta\eta^{\prime})}}{\eta^{\prime}}\frac{d}{d% \eta^{\prime}}\log\left[\frac{(\eta\eta^{\prime})^{2}+\varepsilon^{2}}{(\bar% {\eta})^{2}}\right].  (34) 
Integrating by parts yields
\begin{split}&\displaystyle\left.\frac{1}{2}\,\,\chi_{{\bm{k}}}(\eta^{\prime% })\frac{e^{ik(\eta\eta^{\prime})}}{\eta^{\prime}}\log\left[\frac{(\eta\eta^% {\prime})^{2}+\varepsilon^{2}}{(\bar{\eta})^{2}}\right]\right_{\eta^{\prime}% =\eta_{0}}^{\eta^{\prime}=\eta}\\ &\displaystyle+\int^{\eta}_{\eta_{0}}\!\!d\eta^{\prime}\,\,\frac{d}{d\eta^{% \prime}}\left[e^{i(\eta\eta^{\prime})}\frac{\chi_{{\bm{k}}}(\eta^{\prime})}% {\eta^{\prime}}\right]\log\left[\frac{\eta\eta^{\prime}}{\bar{\eta}}\right].% \end{split}  (35) 
As explained in Boyanovsky:2015jen () by identifying the RG scale with \eta_{0} and by taking k\eta_{0}\gg 1 we can neglect the lower end of the surface term. The upper end gives a UV divergence which is removed by an additive mass renormalisation of \chi.
These steps lead to a master equation consisting of local and nonlocal terms
\frac{d\rho_{I}(\eta)}{d\eta}=\mathcal{D}_{L}[\rho_{I}]+\mathcal{D}_{NL}[\rho_% {I}],  (36) 
where
\begin{split}&\displaystyle\mathcal{D}_{L}[\rho_{I}]=\frac{i\lambda^{2}}{8\pi% ^{2}H^{2}\eta}\log\Big{[}\frac{\varepsilon}{\eta_{0}}\Big{]}\sum_{{\bm{p}}}% \left[\chi_{\bm{q}}(\eta)\chi_{{\bm{p}}}(\eta),\rho_{I}(\eta)\right]\\ &\displaystyle\frac{\lambda^{2}}{16\pi H^{2}\eta^{2}}\sum_{{\bm{p}}}\Big{\{}% \chi_{{\bm{p}}}(\eta)\chi_{{\bm{q}}}(\eta)\rho_{I}(\eta)\\ &\displaystyle+\rho_{I}(\eta)\chi_{\bm{p}}(\eta)\chi_{{\bm{p}}}(\eta)2\chi_{% {\bm{p}}}(\eta)\rho_{I}(\eta)\chi_{{\bm{p}}}(\eta)\Big{\}},\end{split}  (37) 
and,
\begin{split}\displaystyle\mathcal{D}_{NL}[\rho_{I}]&\displaystyle=\frac{i% \lambda^{2}}{8\pi H^{2}\eta^{2}}\sum_{{\bm{p}}}\Big{\{}\chi_{{\bm{p}}}(\eta)X_% {{\bm{q}}}(\eta)\rho_{I}(\eta)\\ &\displaystyle\rho_{I}(\eta)\overline{X}_{{\bm{p}}}(\eta)\chi_{{\bm{p}}}(% \eta)+\chi_{{\bm{p}}}(\eta)\rho_{I}(\eta)\overline{X}_{{\bm{p}}}(\eta)\\ &\displaystyleX_{{\bm{p}}}(\eta)\rho_{I}(\eta)\chi_{{\bm{p}}}(\eta)\Big{\}},% \end{split}  (38) 
with
\begin{split}\displaystyle X_{{\bm{k}}}(\eta)=\int^{\eta}_{\eta_{0}}\!\!d\eta% ^{\prime}\,\,\frac{d}{d\eta^{\prime}}\Big{[}e^{ik(\eta\eta^{\prime})}\frac{% \chi_{{\bm{k}}}(\eta^{\prime})}{\eta^{\prime}}\Big{]}\log\Big{[}\frac{\eta% \eta^{\prime}}{\eta_{0}}\Big{]}\end{split}  (39) 
and
\begin{split}\displaystyle\overline{X}_{{\bm{k}}}(\eta)=\int^{\eta}_{\eta_{0}% }\!\!d\eta^{\prime}\,\,\frac{d}{d\eta^{\prime}}\Big{[}e^{ik(\eta\eta^{\prime}% )}\frac{\chi_{{\bm{k}}}(\eta^{\prime})}{\eta^{\prime}}\Big{]}\log\Big{[}\frac% {\eta\eta^{\prime}}{\eta_{0}}\Big{]}.\end{split}  (40) 
The \varepsilon divergence in the local term is removed by mass renormalisation of \chi. Notice that the nonlocal contribution arises as a result of integrations over longtime environment correlations, and is therefore associated with memory effects. In contrast, the local terms are Markovian, and in the stochastic interpretation that we will develop, corresponds to white noise. This distinction will prove important in our subsequent analysis. This completes our derivation of the master equation and review of the main steps in Boy1 (); Boyanovsky:2015jen ().
IV Interpreting the master equation
We begin by rewriting the master equation in the Schrödinger picture, which allows us to draw close parallels with quantum Brownian motion. To do this, we rewrite the nonlocal contribution in terms of Heisenberg operators at time \eta by writing the modes \chi_{{\bm{k}}}(\eta^{\prime}) in the integral kernel as
\chi(\eta^{\prime})=f(\eta^{\prime},\eta)\chi(\eta)+k^{1}g(\eta^{\prime},\eta% )\Pi(\eta),  (41) 
where f and g are functions to be determined. For now we have suppressed the momentum labels {\bm{k}}. The factor k^{1} is there on dimensional grounds and will drop out of our calculations at a later stage. Notice in particular that the Heisenberg equations of motion for \chi(\eta^{\prime}) imply that f and g must satisfy the mode equation with respect to \eta^{\prime}; namely
f^{\prime\prime}+\left(k^{2}\frac{a^{\prime\prime}}{a}\right)f=0,  (42) 
and similarly for g. Furthermore, the relation
\Pi=\chi^{\prime}\frac{a^{\prime}}{a},  (43) 
means that
\Pi(\eta^{\prime})=\Big{(}f^{\prime}\frac{a^{\prime}}{a}f\Big{)}\chi(\eta)+k^% {1}\Big{(}g^{\prime}\frac{a^{\prime}}{a}g\Big{)}\Pi(\eta).  (44) 
Noting that a^{\prime}/a=1/\eta, we see that f and g must satisfy the boundary conditions
\begin{split}&\displaystyle f(\eta^{\prime}=\eta)=1,\quad\partial_{\eta^{% \prime}}f(\eta^{\prime}=\eta)=1/\eta,\\ &\displaystyle g(\eta^{\prime}=\eta)=0,\qquad\partial_{\eta^{\prime}}g(\eta^{% \prime}=\eta)=k.\end{split}  (45) 
Solving explicitly for the mode functions, we find
\begin{split}\displaystyle g\left(\eta^{\prime},\eta\right)&\displaystyle=\Big% {(}\frac{1}{k\eta^{\prime}}\frac{1}{k\eta}\Big{)}\cos\big{[}k\left(\eta\eta^% {\prime}\right)\big{]}\\ &\displaystyle\Big{(}1+\frac{1}{k^{2}\eta\eta^{\prime}}\Big{)}\sin\big{[}k% \left(\eta\eta^{\prime}\right)\big{]},\end{split}  (46) 
and
\displaystyle f(\eta^{\prime},\eta)  \displaystyle=\frac{\sin\left[k\left(\eta\eta^{\prime}\right)\right]}{k\eta^{% \prime}}+\cos\left[k\left(\eta\eta^{\prime}\right)\right].  (47) 
With these definitions we find that the nonlocal objects X_{{\bm{k}}}(\eta) in (38)(40) can be written as a linear combination of Heisenberg operators \chi(\eta) and \Pi(\eta), multiplied by integrals over the functions f(\eta^{\prime},\eta) and g(\eta^{\prime},\eta), i.e.,
X_{{\bm{k}}}(\eta)=k\,\chi_{{\bm{k}}}(\eta)F(\eta)+\Pi_{{\bm{k}}}(\eta)G(\eta)  (48) 
where
\displaystyle F(\eta)  \displaystyle=\frac{1}{k}\int_{\eta_{0}}^{\eta}d\eta^{\prime}\,\,\frac{d}{d% \eta^{\prime}}\Big{[}f(\eta^{\prime},\eta)\frac{e^{ik(\eta\eta^{\prime})}}{% \eta^{\prime}}\Big{]}\log\Big{[}\frac{\eta\eta^{\prime}}{\eta_{0}}\Big{]},  (49)  
\displaystyle G(\eta)  \displaystyle=\frac{1}{k}\int_{\eta_{0}}^{\eta}d\eta^{\prime}\,\,\frac{d}{d% \eta^{\prime}}\Big{[}g(\eta^{\prime},\eta)\frac{e^{ik(\eta\eta^{\prime})}}{% \eta^{\prime}}\Big{]}\log\Big{[}\frac{\eta\eta^{\prime}}{\eta_{0}}\Big{]},  (50) 
and correspond to longtime correlations in the environment. We can now perform the switch to the Schrödinger picture. Using the result
U(\eta,\eta_{0})\chi_{\bm{k}}(\eta)U^{1}(\eta,\eta_{0})=\chi_{\bm{k}}  (51) 
and the relation between the interaction density matrix and the (timedependent) Schrödinger picture density operator
\rho_{I}(\eta)=U_{0}^{1}(\eta,\eta_{0})\rho(\eta)U_{0}(\eta,\eta^{\prime}),  (52) 
as well as
\frac{dU_{0}(\eta,\eta_{0})}{d\eta}=i\left[H_{0}(\eta),\rho\right],  (53) 
we find that after appropriate insertions of UU^{1} in the standard way, one obtains a Schrödinger representation of the master equation
\begin{split}\displaystyle\frac{d\rho(\eta)}{d\eta}&\displaystyle=i[H_{0},% \rho(\eta)]\frac{\lambda^{2}}{16\pi H^{2}\eta^{2}}\sum_{{\bm{k}}}\Big{\{}\chi% _{{\bm{k}}}\chi_{{\bm{k}}}\rho(\eta)+\rho(\eta)\chi_{{\bm{k}}}\chi_{{\bm{k}}% }2\chi_{{\bm{k}}}\rho(\eta)\chi_{{\bm{k}}}\Big{\}}\\ &\displaystyle+\frac{ik\lambda^{2}}{8\pi^{2}H^{2}\eta}\sum_{{\bm{k}}}\Big{\{}F% (\eta)\left[\chi_{{\bm{k}}}\chi_{{\bm{k}}}\rho(\eta)\chi_{{\bm{k}}}\rho(% \eta)\chi_{{\bm{k}}}\right]F^{*}(\eta)\left[\rho(\eta)\chi_{{\bm{k}}}\chi_{{% \bm{k}}}\chi_{{\bm{k}}}\,\rho(\eta)\chi_{{\bm{k}}}\right]\Big{\}}\\ &\displaystyle+\frac{i\lambda^{2}}{8\pi^{2}H^{2}\eta}\,\sum_{{\bm{k}}}\Big{\{}% G(\eta)\left[\chi_{\bm{k}}\Pi_{{\bm{k}}}\rho(\eta)\Pi_{{\bm{k}}}\rho(\eta)% \chi_{{\bm{k}}}\right]G^{*}(\eta)\left[\rho(\eta)\Pi_{{\bm{k}}}\chi_{{\bm{k}% }}\chi_{{\bm{k}}}\,\rho(\eta)\Pi_{{\bm{k}}}\right]\Big{\}},\end{split}  (54) 
The Hamiltonian in (54) is
\begin{split}\displaystyle H_{0}=\frac{1}{2}\sum_{\bm{k}}\Big{[}\Pi_{{\bm{k}}}% \Pi_{{\bm{k}}}+k^{2}\chi_{{\bm{k}}}\chi_{{\bm{k}}}+\frac{a^{\prime}}{a}\left% (\chi_{{\bm{k}}}\Pi_{{\bm{k}}}+\Pi_{{\bm{k}}}\chi_{{\bm{k}}}\right)\Big{]}.% \end{split}  (55) 
The full density matrix can then be factorised as
\rho(\eta)=\prod_{{\bm{k}}}\otimes\,\rho({\bm{k}},{\bm{k}})+\mathcal{O}(% \lambda^{3}),  (56) 
where the modes of different wavelengths only mix at \mathcal{O}(\lambda^{3}). Under this assumption (remembering to pick out the +{\bm{k}} and {\bm{k}} terms from the mode sum in eq. (54)) we find, dropping the {\bm{k}} subscripts, the following master equation for the twomode density matrix \rho({\bm{k}},{\bm{k}})
\begin{split}\displaystyle\frac{d\rho}{d\eta}&\displaystyle=\overbracket{i% \Big{[}H_{\text{eff}},\rho\Big{]}}^{\text{unitary evol.}}+\overbracket{\frac{i% \lambda^{2}}{8\pi^{2}H^{2}\eta}\,\text{Re}\,G\,\Big{\{}\left[\chi,\{\Pi^{% \dagger},\rho\}\right]+\left[\chi^{\dagger},\{\Pi,\rho\}\right]\Big{\}}}^{% \text{dissipation/relaxation}}\\ &\displaystyle\underbracket{\frac{\lambda^{2}}{16\pi^{2}H^{2}\eta^{2}}\left(% \pi+2\,k\,\eta\,\,\text{Im}\,F\right)\Big{\{}\left[\chi^{\dagger},\left[\chi,% \rho\right]\,\right]+\left[\chi,\left[\chi^{\dagger},\rho\right]\,\right]\Big{% \}}}_{\text{fluctuation/decoherence}}\frac{\lambda^{2}}{8\pi^{2}H^{2}\eta}\,% \text{Im}\,G\,\Big{\{}\left[\chi,\left[\Pi^{\dagger},\rho\right]\right]+\left[% \chi^{\dagger},\left[\Pi,\rho\right]\right]\Big{\}},\end{split}  (57) 
where the renormalised Hamiltonian H_{\text{eff}}=H_{0}+\delta H contains the contributions
\begin{split}\displaystyle\delta H&\displaystyle=\frac{k\lambda^{2}}{4\pi^{2}% H^{2}\eta}\,\,\text{Re}\,F\,\chi^{\dagger}\chi,\\ \displaystyle H_{0}&\displaystyle=k^{2}\chi^{\dagger}\chi+\Pi^{\dagger}\Pi+% \frac{a^{\prime}}{2a}\left(\chi\Pi^{\dagger}+\Pi\chi^{\dagger}+\chi^{\dagger}% \Pi+\Pi^{\dagger}\chi\right).\end{split}  (58) 
This reveals striking similarities between the dynamics of inflation and quantum Brownian motion as described in sec. II. This provides a deeper insight into the master equation first presented in Boy1 () and will prove important in subsequent analysis. It also paints a richer and more detailed picture than the simple master equation presented in BLM () (see for instance their eq. (10)). However, as in BLM () we too shall see that the master equation is dominated by a simple diffusion term on superhorizon scales. The first term on the righthand side describes the unitary dynamics (with renormalised Hamiltonian) and the second and third terms correspond respectively to dissipation (i.e. a friction term) and fluctuations/decoherence of the system due to the environmental noise. The final term, which does not fit the quantum Brownian motion scheme, will be discussed more shortly. These first three terms are easily identified with corresponding terms in the CaldeiraLeggett master equation of quantum Brownian motion (6). Whilst the similarities are striking, an important difference is that in the cosmological setting the couplings of the various terms are timedependent.
The interpretation of the dissipationfluctuation terms becomes more apparent if we write the master equation in terms of matrix elements in the field basis \rho(\chi_{{\bm{k}}},\tilde{\chi}_{{\bm{k}}})\equiv\bra{\chi_{{\bm{k}}}}{\rho}% ({\bm{k}},{\bm{k}})\ket{\tilde{\chi}_{{\bm{k}}}}. Given any operator M, the action of \Pi is given by
\bra{\chi}{M}{\Pi}^{\dagger}\ket{\tilde{\chi}}=i\frac{\partial}{\partial\tilde% {\chi}}\bra{\chi}{M}\ket{\tilde{\chi}},\qquad\bra{\chi}{M}{\Pi}\ket{\tilde{% \chi}}=i\frac{\partial}{\partial\tilde{\chi}^{*}}\bra{\chi}{M}\ket{\tilde{\chi% }}.  (59) 
Theses relations follow from the commutation relations [\chi,\Pi^{\dagger}]=1 and the fact that \exp(i\chi\Pi^{\dagger}) is the generator of translations. The action of other terms involving \Pi can be inferred from these expressions by appropriate complex conjugation. This gives rise to the following equation for the matrix elements:
\begin{split}\displaystyle\frac{d\rho}{d\eta}&\displaystyle=i\Big{(}k^{2}% \frac{k\lambda^{2}}{4\pi^{2}H^{2}\eta}\,\,\text{Re}\,F\Big{)}\left(\chi^{2}% \tilde{\chi}^{2}\right)+i\left(\frac{\partial^{2}\rho}{\partial\chi\partial% \chi^{*}}\frac{\partial^{2}\rho}{\partial\tilde{\chi}\partial\tilde{\chi}^{*}% }\right)\frac{a^{\prime}}{a}\left(\chi\frac{\partial\rho}{\partial\chi}+\chi^% {*}\frac{\partial\rho}{\partial\chi^{*}}+\tilde{\chi}\frac{\partial\rho}{% \partial\tilde{\chi}}+\tilde{\chi}^{*}\frac{\partial\rho}{\partial\tilde{\chi}% ^{*}}+2\rho\right)\\ &\displaystyle+\frac{\lambda^{2}}{4\pi^{2}H^{2}\eta}\,\text{Re}\,G\left[\left(% \chi\tilde{\chi}\right)\left(\frac{\partial\rho}{\partial\chi}\frac{\partial% \rho}{\partial\tilde{\chi}}\right)+\left(\chi^{*}\tilde{\chi}^{*}\right)\left% (\frac{\partial\rho}{\partial\chi^{*}}\frac{\partial\rho}{\partial\tilde{\chi% }^{*}}\right)\right]\\ &\displaystyle\frac{\lambda^{2}}{8\pi^{2}H^{2}\eta^{2}}\left(\pi+2\,k\eta\,\,% \text{Im}\,F\right)\left\chi\tilde{\chi}\right^{2}\,\rho\\ &\displaystyle+\frac{i\lambda^{2}}{4\pi^{2}H^{2}\eta}\,\text{Im}\,G\left[\left% (\chi\tilde{\chi}\right)\left(\frac{\partial\rho}{\partial\chi}+\frac{% \partial\rho}{\partial\tilde{\chi}}\right)+\left(\chi^{*}\tilde{\chi}^{*}% \right)\left(\frac{\partial\rho}{\partial\chi^{*}}+\frac{\partial\rho}{% \partial\tilde{\chi}^{*}}\right)\right],\end{split}  (60) 
which can be compared with the Brownian motion case in equation (6). One can also use this master equation to study the phasespace of the theory as captured by the Wigner function, cf. (9),
\begin{split}&\displaystyle W(\chi,\Pi)=\\ &\displaystyle\frac{1}{\pi^{2}}\int^{\infty}_{\infty}d^{2}\chi^{\prime}\rho% \left(\chi\frac{\chi^{\prime}}{2},\chi+\frac{\chi^{\prime}}{2}\right)e^{i2% \chi^{\prime}_{R}\Pi_{R}+i2\chi^{\prime}_{I}\Pi_{I}},\end{split}  (61) 
which provides a quantum analogue of a phasespace probability distribution. By Wigner transforming the master equation (60) we can derive the following equation for W, mirroring the case of quantum Brownian motion (10),
\begin{split}\displaystyle\frac{\partial W}{d\eta}&\displaystyle=\{{H_{\text{% eff}},W\}}\frac{\lambda^{2}}{4\pi^{2}H^{2}\eta}\,\text{Re}\,G\left(\partial_{% \Pi}\,(\,\Pi W)+\partial_{\overline{\Pi}}(\bar{\Pi}W)\right)\\ &\displaystyle+\frac{\lambda^{2}}{8\pi^{2}H^{2}\eta^{2}}\left(\pi+2k\eta\,% \text{Im}F\right)\partial_{\Pi}\partial_{\bar{\Pi}}W\\ &\displaystyle\frac{\lambda^{2}}{8\pi^{2}H^{2}\eta}\,\text{Im}\,G\left(% \partial_{\chi}\partial_{\bar{\Pi}}W+\partial_{\bar{\chi}}\partial_{\Pi}W% \right),\end{split}  (62) 
This has the usual terms: decoherence (the diffusion terms involving twoderivatives), the dissipation/friction term (single derivatives) and the renormalisation of the Hamiltonian. Note that in our model, the field \chi has no selfinteractions and so there is no analogue of the quantum terms in (10).
The diffusion term on the righthand side can be written in the form \mathcal{D}_{ij}\partial_{i}\overline{\partial}_{j}W, where the derivatives are defined by \partial_{1}=\partial_{\chi} and \partial_{2}=k\partial_{\Pi}, with the diffusion matrix defined by
\mathcal{D}=\frac{\lambda^{2}}{8\pi^{2}H^{2}k\eta}\left(\begin{array}[]{cc}0&% \text{Im}\,G\\ \text{Im}\,G&\quad(\pi+2k\eta\,\text{Im}F)/(k\eta)\end{array}\right).  (63) 
This has eigenvalues
\begin{split}\displaystyle\mathcal{\mu}_{\pm}&\displaystyle=\frac{\lambda^{2}}% {16\pi^{2}H^{2}k^{2}\eta^{2}}(\pi+2k\eta\,\text{Im}F)\\ &\displaystyle\pm\frac{\lambda^{2}}{16\pi^{2}H^{2}k^{2}\eta^{2}}\sqrt{\left(% \pi+2k\eta\,\text{Im}F\right)^{2}+4\left(k\eta\text{Im}\,G\right)^{2}}.\end{split}  (64) 
We see immediately that the eigenvalue \mu_{}<0. In fact, the plot in Fig. 1 shows that these eigenvalues remain comparable in magnitude, so that positive and negative diffusion, at first sight appear to be of equal importance. This can also be seen from the asymptotic form of F and G for \eta k\ll 1
\begin{split}\displaystyle F(\eta)&\displaystyle\simeq\frac{1}{k\eta}\Big{(}1% \log\frac{\eta}{\eta_{0}}\Big{)},\\ \displaystyle G(\eta)&\displaystyle\simeq\frac{i}{2k\eta}+\frac{1}{3}\Big{(}1% \log\frac{\eta}{\eta_{0}}\Big{)},\end{split}  (65) 
which reveal that \mu_{\pm} scale as \pm 1/\eta^{2} on superhorizon scales. This “antidiffusion” has no classically meaningful analogue, and can be traced to the \text{Im}\,G term in master equation (60).
The other violation of positivity comes from the behaviour of the coefficients in the master equation, which exhibit oscillatory behaviour, fluctuating sign at early times. Similar behaviour is seen for coefficients of the master equation in flat space Boy1 (), where the oscillations are due to the nonlocal environment correlations corresponding to memory effects. These effects are usually dealt with by temporal coarsegraining in the form of a kind of rotating wave approximation. This approximation is key to writing the master equation in Lindblad, and hence manifestly positivity preserving form.
At this stage, our situation may seem a little perilous: we have a master equation which is not in Lindblad form and exhibits several features, e.g. antidiffusion and oscillatory coefficients, which might jeopardise the positivity of the density matrix. In fact, although the unitary evolution of the total state ensures positivity of the full density matrix, the series of simplifying approximations outlined in sec. II and carried out explicitly in sec. III provide no guarantee that the spectrum of the reduced density matrix remains positive. Whilst Hermiticity and normalisation follow trivially, the most important property – that \rho(\eta) be positive definite – remains to be established. Indeed, even the CaldeiraLeggett model can lead to a small violation of positivity H (); BP (). In fact, we shall see the operator form of the master equation and the oscillatory behaviour of the coefficients mean that the master equation manifestly violates positivity on subhorizon scales.
However the master equation alone does not tell the full story. The real question is whether physically relevant solutions of the master equation maintain positivity as a function of time. This shall be the subject of the following sections.
IV.1 Solving the master equation
In accordance with our previous boundary conditions, we shall take Bunch Davis initial conditions which imply a Gaussian ansatz for the density matrix
\rho(\chi,\tilde{\chi})=\mathcal{N}\exp\left\{\Omega\chi^{2}\Omega^{*}\tilde% {\chi}^{2}\frac{1}{2}\xi\chi\tilde{\chi}^{2}\right\},  (66) 
and initial conditions
\Omega(\eta_{0})=1,\qquad\xi(\eta_{0})=0,\qquad\mathcal{N}(\eta_{0})=\frac{2% \Omega_{R}}{\pi},  (67) 
where \Omega=\Omega_{R}+i\Omega_{I} is complex and \xi is real, ensuring the Hermiticity of the density matrix. Notice that although we use the notation \Omega_{R} for the Gaussian width, it does not describe a simple isolated squeezed state, and, in particular, is not a simple expression which can be written in terms of squeezing parameters r_{k} and \varphi_{k} of e.g. Lim (); VM (). On the contrary, \Omega contains \mathcal{O}(\lambda^{2}) corrections and in general must be solved numerically, as we shall see below.
Substituting the ansatz (66) into eq. (60), and comparing the coefficients of \chi^{2}, \tilde{\chi}^{2}, \chi\tilde{\chi}^{*} and \tilde{\chi}\tilde{\chi} gives the evolution equations for \Omega and \xi and \mathcal{N}:
\displaystyle\frac{d\Omega}{d\eta}  \displaystyle=ik^{2}+\frac{2}{\eta}\Omegai\Omega_{R}\xii\Omega^{2}\frac{i% \lambda^{2}}{4\pi^{2}H^{2}\eta}\Big{(}k\,\text{Re}\,F\,\text{Im}\left(G\Omega% \right)\Big{)},  (68)  
\displaystyle\frac{d\xi}{d\eta}  \displaystyle=2\xi\Big{(}\Omega_{I}+\frac{1}{\eta}\Big{)}+\frac{\lambda^{2}}{4% \pi^{2}H^{2}\eta^{2}}\Big{(}\pi+2\eta k\,\text{Im}\,F+2\eta\,\text{Re}\,\left(% G\Omega\right)+2\eta\,\text{Re}\,G\xi\Big{)},  (69)  
\displaystyle\frac{d\mathcal{N}}{d\eta}  \displaystyle=2\,\mathcal{N}\Big{(}\Omega_{I}+\frac{1}{\eta}\Big{)}.  (70) 
Note that F and G are in fact dimensionless, and functions only of the dimensional combination k\eta, as can be seen from eqs. (49) and (50).
We immediately recognise that \xi governs decoherence since it controls the suppression of the offdiagonal elements of \rho, corresponding to the decay of quantum interference. Notice also that by taking the real part of eq. (68) and using the normalisation condition (70) we get
\frac{1}{\Omega_{R}}\frac{d\Omega_{R}}{d\eta}=\frac{2}{\eta}+2\Omega_{I}=\frac% {1}{\mathcal{N}}\frac{d\mathcal{N}}{d\eta},  (71) 
whose solutions give the expected normalisation condition \mathcal{N}=2\Omega_{R}/\pi for a (complex) Gaussian, providing a consistency check for the evolution equations.
Using these equations, and the asymptotic behaviour for F and G given in eq. (65) we find the following superhorizon behaviour for \Omega_{R}, \Omega_{I} and \xi in the k\eta\ll 1 limit:
\displaystyle\Omega_{R}  \displaystyle=k^{3}\eta^{2}\exp\left[\frac{\lambda^{2}}{12\pi^{2}H^{2}}\Big{(}% \log^{2}\frac{\eta}{\eta_{0}}\frac{4}{3}\log\frac{\eta}{\eta_{0}}\Big{)}\right]  
\displaystyle+O\left(\eta^{3}\right),  (72)  
\displaystyle\Omega_{I}  \displaystyle=\frac{\lambda^{2}}{12\pi^{2}H^{2}}\frac{1}{\eta}\left(\log\frac{% \eta}{\eta_{0}}\frac{2}{3}\right)k^{2}\eta+O\left(\eta^{3}\right),  (73)  
\displaystyle\xi  \displaystyle=\frac{\lambda^{2}}{12\pi^{2}H^{2}}\frac{1}{\eta}+O\left(\eta% \right).  (74) 
These expressions will be used throughout the remainder of this paper to study the superhorizon behaviour of fluctuations.
It is noteworthy that the expression for \Omega_{R} in (72) has the form of a resummed quantity in perturbation theory because it is \Omega_{R} that determines the allimportant power spectrum that is a key observable for the CMB. This resummed expression was first derived in Boyanovsky:2015jen () and in order to make contact with that analysis, let us focus in more detail on the equation for \Omega_{I}. Including the bare mass term for \chi and the UV divergence, this equation takes the form
\begin{split}\displaystyle\frac{d\Omega_{I}}{d\eta}&\displaystyle=k^{2}+\frac{% {\cal M}_{\chi}^{2}}{H^{2}\eta^{2}}+\frac{\lambda^{2}}{4\pi^{2}H^{2}\eta^{2}}% \log\frac{\varepsilon}{\eta_{0}}\\ &\displaystyle+\frac{2}{\eta}\Omega_{I}+\frac{\lambda^{2}}{4\pi^{2}H^{2}\eta^{% 2}}\Big{(}1\log\frac{\eta}{\eta_{0}}\Big{)}+\cdots.\end{split}  (75) 
Note that this expression makes clear that the subtraction point we chose to be \eta_{0} is in fact arbitrary. Let us set the renormalized mass to zero but include a specific finite counterterm:
\begin{split}\displaystyle{\cal M}_{\chi}^{2}=\frac{\lambda^{2}}{4\pi^{2}}% \log\frac{\varepsilon}{\eta_{0}}\frac{\lambda^{2}}{6\pi^{2}}.\end{split}  (76) 
Solving (75) and then integrating (71) to find \Omega_{R} (and fixing the constant of integration appropriately), we find
\begin{split}\displaystyle\Omega_{R}&\displaystyle=k^{3}\eta^{2}\exp\Big{[}% \frac{\lambda^{2}}{12\pi^{2}H^{2}}\Big{(}\log^{2}(k\eta)\\ &\displaystyle\qquad\qquad2\log(k\eta)\log(k\eta_{0})\Big{)}\Big{]}.\end{split}  (77) 
This quantity gives directly the power spectrum of the scalar perturbations:
\begin{split}\displaystyle{\cal P}(k,\eta)=\frac{k^{3}H^{2}\eta^{2}}{2\pi^{2}}% \cdot\langle\chi_{\bm{k}}\chi_{\bm{k}}^{\dagger}\rangle=\frac{k^{3}H^{2}\eta^{% 2}}{(2\pi)^{2}}\cdot\frac{1}{\Omega_{R}}=\frac{H^{2}}{(2\pi)^{2}}\exp\Big{[}% \frac{\lambda^{2}}{12\pi^{2}H^{2}}\Big{(}2\log(k\eta_{0})\log(k\eta)\log^{2% }(k\eta)\Big{)}\Big{]}.\end{split}  (78) 
This is the result quoted in Boyanovsky:2015jen () for the correction to the power spectrum. Note that the coefficient of the \log(k\eta) in the exponent can be altered by changing the finite part of the counterterm while the \log^{2}(k\eta) term is universal.
IV.2 Positivity on super horizon scales and the nature of decoherence
Positivity is guaranteed if, and only if, \rho has nonnegative eigenvalues. To find these, we note that the density matrix (elements) factorise into a product of two identical Gaussians of the real and imaginary parts of \chi, i.e.,
\rho(\chi,\tilde{\chi})=\rho_{I}(\chi_{I},\tilde{\chi}_{I})\rho_{R}(\chi_{R},% \tilde{\chi}_{R}),  (79) 
where
\begin{split}&\displaystyle\rho_{R}(\chi_{R},\tilde{\chi}_{R})\\ &\displaystyle=\sqrt{\frac{2\Omega_{R}}{\pi}}\exp\Big{[}\Omega\chi_{R}^{2}% \Omega^{*}\tilde{\chi}_{R}^{2}\frac{\xi}{2}(\chi_{R}\tilde{\chi}_{R})^{2}% \Big{]},\end{split}  (80) 
with an identical expression for \rho(\chi_{I},\tilde{\chi}_{I}) given by replacing \chi_{R}\rightarrow\chi_{I}. The eigenstates of \rho are then simply tensor products \rho_{I} and \rho_{R} eigenstates, so that positivity of \rho is then equivalent to positivity of \rho_{I} and \rho_{R}. Following JZ (), we can write \rho_{R} in diagonal form as
\rho_{R}=\sum_{n}p_{n}\ket{\varphi_{n}}\bra{\varphi_{n}}.  (81) 
The eigenstates \ket{\varphi_{n}} and their eigenvalues p_{n} satisfy (dropping the R subscript on \chi)
p_{n}\,\varphi_{n}(\chi)=\int d\chi^{\prime}\rho_{R}(\chi,\chi^{\prime})% \varphi_{n}(\chi^{\prime}),  (82) 
where \varphi_{n}(\chi)\equiv\braket{\chi}{\varphi_{n}} gives the position representation of the \left\{\ket{\varphi_{n}}\right\}. The Gaussian form of \rho_{R} means that \varphi_{n}(x) take the same form as harmonic oscillator wave functions, i.e.,
\varphi_{n}(\chi)=NH_{n}(\alpha\chi)\exp(\beta\chi^{2}),  (83) 
where H_{n} are the Hermite polynomials, N is a normalisation factor, and \alpha and \beta are functions to be determined. Substituting this into (82) and using the integral identity in GR () we find
\begin{split}\displaystyle\alpha^{2}&\displaystyle=2\left[\Omega_{R}(\Omega_{R% }+\xi)\right]^{1/2},\\ \displaystyle\beta&\displaystyle=i\Omega_{I}+\sqrt{\Omega_{R}(\Omega_{R}+\xi)}% .\end{split}  (84) 
We can also read off
p_{n}=p_{0}(1p_{0})^{n},  (85) 
where
\displaystyle p_{0}  \displaystyle=\frac{2}{1+\sqrt{1+\xi/\Omega_{R}}}.  (86) 
Positivity requires 0\leq p_{0}\leq 1, which, assuming \Omega_{R}>0, holds if and only if \xi>0. From the numerical solution for \xi in Fig. 2, one can see that at early times positivity is violated, but that \xi becomes positive as the modes exit the horizon.
This can be understood by examining the equation (69) for \xi, which can be rewritten in the form
\begin{split}\displaystyle\frac{d\xi}{d\eta}&\displaystyle=\Big{[}\frac{% \lambda^{2}}{2\pi^{2}H^{2}\eta}\,\text{Re}\,G+\frac{1}{\Omega_{R}}\frac{d% \Omega_{R}}{d\eta}\Big{]}\xi\\ &\displaystyle+\frac{\lambda^{2}}{4\pi^{2}H^{2}\eta^{2}}\left[\pi+2k\eta\,% \text{Im}\,F+2\eta\,\text{Re}\,(G\Omega)\right].\end{split}  (87) 
One notes that on superhorizon scales (k\eta\ll 1) the nonlocal (memory) effects in the first square braket [...]\xi corresponding to \text{Re}(G) become subdominant, giving way to the leading behaviour from \Omega_{R}^{1}{d\Omega_{R}}/{d\eta}\simeq 2/\eta. Similarly the “forcing term” in the second squarebraket is dominated by the first term (which is ultimately traced to the local contribution in the master equation) with the nonlocal effects from \text{Im}F and \text{Re}(G\Omega) becoming negligible. Thus on superhorizon scales we find
\frac{d\xi}{d\eta}\simeq\underbracket{\,\,\frac{2}{\eta}\,\,\xi}_{\text{% squeezing}}+\underbracket{\,\,\frac{\lambda^{2}}{4\pi H^{2}\eta^{2}}\,\,}_{% \text{interaction}},  (88) 
The second term in (88) drives the positive growth of \xi. It originates from the local term in the first line of (54) and is associated to the decoherence term in (60). Equivalently, the relevant superhorizon dynamics in the evolution equation for \xi can be reproduced by replacing the environment correlators according to
K(k,\eta,\eta^{\prime})\longrightarrow\frac{1}{8\pi^{2}}\delta(\eta^{\prime}% \eta),  (89) 
back in eq. (30); in other words, long time environmental correlations become unimportant for decoherence at late times, which is driven only by, what will become in the stochastic interpretation, white noise. This is a common feature of many systems, where they exhibit essentially memoryless (Markov) behaviour at late times, when long time environmental correlations become negligible. This is precisely the case for the flat space master equation considered in Boy1 (), where initially oscillatory coefficients in the master equation, induced by shorttime environment correlations lead to positivity on longtimes. Hence, we see that the oscillatory behaviour in the interaction terms associated to the kernels F and G, which can be traced to nonlocal/memory effects, leads to an initial violation of positivity. However, at late times, this gives way rather elegantly to the dominance of effects driven by local terms and their white noise leading to
\xi\simeq\frac{\lambda^{2}}{12\pi H^{2}}\frac{1}{\eta}\quad\text{for}\quadk% \eta\ll 1.  (90) 
It should be noted in passing however, that by contrast, the object \Omega_{R} is dominated by nonlocal contributions to the master equation on superhorizon scales. In this sense, the power spectrum receives corrections due to memory effects – see eq. (78). We mention this purely to point out that although \xi becomes dominated by local effects, there remain other physically interesting quantities which are dominated by memory effects in the master equation. Indeed, the resummation of nonlocal contributions and their secular log terms is one of the key results in Boyanovsky’s paper Boyanovsky:2015jen ().
The growth in \xi can be traced back to the fact the interaction strength is effectively timedependent and grows with time, as can be seen from eq. (3), where one could in principle absorb the timedependence into the coupling constant by taking \lambda^{2}\rightarrow\lambda^{2}(\eta)\equiv\lambda^{2}/H^{2}\eta^{2}. It is therefore interesting to ask how generic this growth in \xi is, and what happens when one looks at other interactions, e.g. marginal couplings such as \mathcal{H}_{int}=\sqrt{g}\phi^{2}\varphi^{2} as well as other relevant couplings to see how modeldependent the scaledependence of decoherence is. In addition this might save us from some of the problems associated with having a coupling constant that grows in time, thus invalidating perturbativity at late times. We leave such investigations for future work.
At this point we also highlight an important subtlety in the evolution equation for \xi. Let us consider the first term in eq. (88) of the form 2/\eta which corresponding to squeezing. This can be traced to the \Omega_{R}^{1}\Omega_{R}^{\prime} term in (87) or equivalently to the a^{\prime}/a term in the Hamiltonian (58). Notice that this term gives a negative contribution to d\xi/d\eta. At first sight, this makes things look as though the squeezing is slowing the increase in \xi and therefore resisting the suppression of offdiagonal elements. This would certainly run counter to the party line that squeezing makes a state more classical. However, such reasoning is misleading as we now explain.
This can be seen most clearly in the free case with \lambda=0. Explicitly, the evolution equations become
\begin{split}\displaystyle\frac{d\xi}{d\eta}&\displaystyle=\frac{1}{\Omega_{R}% }\frac{d\Omega_{R}}{d\eta}\xi,\\ \displaystyle\frac{d\Omega}{d\eta}&\displaystyle=ik^{2}+\frac{2\Omega}{\eta}i% \Omega_{R}\xii\Omega^{2}.\end{split}  (91) 
One finds that the solutions satisfy \xi\propto\Omega_{R}. Since \Omega_{R}(\eta_{0})=1, solutions are given by \xi=\xi_{0}\Omega_{R}. Again, this makes it very tempting to suggest \Omega_{R} drives \xi to zero. However, \rho takes the form
\rho_{\text{free}}(\chi,\tilde{\chi})=\frac{2\Omega_{R}}{\pi}\exp\left(\Omega% _{R}\left[\chi+\tilde{\chi}+\frac{\xi_{0}}{2}\chi\tilde{\chi}^{2}\right% ]\right).  (92) 
From this we see that \Omega_{R} does decrease the size of \xi, and therefore increases the offdiagonal extent of the density matrix, but only in the sense that as time passes, it rescales the entire ellipse
\chi+\tilde{\chi}+\frac{\xi_{0}}{2}\chi\tilde{\chi}^{2}=\Omega_{R}^{1},  (93) 
corresponding to surfaces of \rho(\chi,\tilde{\chi})=\text{constant}. This is illustrated in the figure 3. This clearly shows how although \Omega_{R}\rightarrow 0 does cause \xi to shrink (as offdiagonal elements grow in size), the relative size of diagonal and offdiagonal elements stays the same, as can be seen from the fact the contours retain their shape. So in actual fact, squeezing only decreases \xi in a trivial sense in that \Omega_{R} scales the size of all density matrix elements and so contrary to initial impressions, squeezing does not, in any real sense, resist the onset of decoherence. Explicitly, in the free case, the Wigner area and entanglement entropy depend only on the ratio \xi/\Omega_{R} which is constant as explained above. Therefore these quantities remain constant – consistent with the idea that squeezing cannot make the state any “lessdecohered” since it is associated to unitary evolution alone. Similarly the purity \text{Tr}[\rho^{2}] remains constant due to the unitary evolution of the von Neumann equation. When interactions are switched on, we again have the \Omega_{R}^{1}\Omega_{R}^{\prime} term, but the situation is exactly the same – this acts only as an overall rescaling and comes from the unitary dynamics.
IV.3 Langevin equation
We now show how it is possible to derive a classical Langevin equation which reproduces the behaviour of the density matrix \rho(\chi,\tilde{\chi}) on superHubble scales. By this, we mean that, mirroring the case of quantum Brownian motion, the master equation for the Wigner function W(\chi,\Pi) can be interpreted as a FokkerPlanck equation associated to an auxiliary classical stochastic system whose dynamics is described by a Langevin equation. The interpretation is then to go beyond the mere idea that the classical system is auxiliary to the idea that it describes the classical (stochastic) trajectories that have emerged from the quantum system.
To study the superhorizon behaviour of solutions, we return to equations (68) and (69). From the asymptotic behaviour of G and F, (65) and eqs. (72)(74), we see that on superhorizon scales the contribution to solutions from \,\text{Im}\,G is subleading because it is multiplied by \Omega_{R} which goes to 0 like \eta^{2}. Thus, at the level of the solutions, squeezing ensures that the antidiffusion term turns out to be negligible. The latetime behaviour is therefore captured by the FokkerPlancklike equation
\displaystyle\frac{\partial W}{d\eta}  \displaystyle=\{{H_{\text{eff}},W\}}\frac{\lambda^{2}}{4\pi^{2}H^{2}\eta}\,% \text{Re}\,G\left(\partial_{\Pi}\,(\,\Pi W)+\partial_{\overline{\Pi}}(\bar{\Pi% }W)\right)  
\displaystyle+\frac{\lambda^{2}}{8\pi^{2}H^{2}\eta^{2}}\left(\pi+2k\eta\,\text% {Im}\,F\right)\partial_{\Pi}\partial_{\bar{\Pi}}W,  (94) 
valid in the region k\eta\to 0, which, by virtue of the latetime behaviour of F, has a positive diffusion matrix.
We now go the extra step and argue that from the perspective of the superHubble physics, a classical stochastic dynamics emerges from the underlying quantum system whose FokkerPlanck equation is (94). The individual trajectories of the classical system are solutions of the Langevin equations for the classical variable \chi:
\begin{split}\displaystyle\chi^{\prime}&\displaystyle=\Pi\frac{\chi}{\eta},\\ \displaystyle\Pi^{\prime}&\displaystyle=\Big{(}\frac{1}{\eta}+\frac{k\lambda^{% 2}}{4\pi^{2}H^{2}\eta}\,\text{Re}\,G\Big{)}\Pi\\ &\displaystyle\Big{(}k^{2}+\frac{{\cal M}_{R}^{2}}{H^{2}\eta^{2}}\frac{k% \lambda^{2}}{4\pi^{2}H^{2}\eta}\,\text{Re}\,F\Big{)}+\sigma,\end{split}  (95) 
where we have allowed for a renormalised mass. The above can be rendered as a single equation for \chi:
\begin{split}&\displaystyle\chi^{\prime\prime}\frac{k\lambda^{2}}{4\pi^{2}H^{% 2}\eta}\,\text{Re}\,G\,\chi^{\prime}+\Big{(}k^{2}+\frac{{\cal M}_{R}^{2}}{H^{2% }\eta^{2}}\\ &\displaystyle\frac{2}{\eta^{2}}\frac{k\lambda^{2}}{4\pi^{2}H^{2}\eta}\,% \text{Re}\,F\frac{k\lambda^{2}}{4\pi^{2}H^{2}\eta^{2}}\,\text{Re}\,G\Big{)}% \chi=\sigma,\end{split}  (96) 
where the white noise \sigma satisfies
\mathbb{E}(\sigma(\eta))=0,\qquad\mathbb{E}(\sigma(\eta)\sigma(\eta^{\prime}))% =\frac{\lambda^{2}}{8\pi H^{2}\eta^{2}}\delta(\eta\eta^{\prime}).  (97) 
This provides a set of stochastic equations capable of reproducing the same superHubble dynamics of the quantum master equation in the sense that we have explained in sec. II. It should be noted that this is quite different from the kind of stochastic equations encountered in the usual stochastic inflation paradigm Star1 (); SY () where FokkerPlanck equations emerge even in the absence of interactions, where subHubble modes are interpreted as noise. By contrast, our noise arises from a genuine coupling of different Fourier modes between \chi and \psi so that there is no noise in the \lambda\rightarrow 0 limit.
It is interesting that the alternative formalism for dealing with open quantum systems, the influence action approach leads directly to a Langevin equation. For the scalar perturbations, the influenceaction Langevin equation was derived in Boyanovsky:2015jen (). This equation appears to be very different from the one we have written down here in that it features a memory integral and also Gaussian coloured noise. In fact, one can easily show that if one substitutes the zeroth order expression for \chi(\eta^{\prime}) from (41) in the memory integral on the righthand side of eq. (3.41) of Boyanovsky:2015jen () one gets precisely the homogeneous terms in our Langevin equation (96). Therefore, the difference between the Langevin equations is simply the fact that our equation (96) has only the local (white) noise component of the coloured noise term of the Langevin equation in Boyanovsky:2015jen (). But we have argued that this is a good approximation for superhorizon modes.
We now solve the equations (68)(70) numerically and examine the evolution of various quantities which characterise the emergence of classical stochastic behaviour as modes exit the horizon.
Since from now on we shall work modebymode, we shall work in units where k=1 for the remainder of the paper, which greatly simplifies notation by effectively removing any k from most expressions. Dimensionful expressions can always be restored by appropriate insertions of k. Throughout what follows we shall take the indicative values
\frac{\lambda}{H}=0.1,\qquad k\eta_{0}(=\eta_{0})=20.  (98) 
Obviously the initial value for k\eta_{0} is unrealistically small, but the results are rather insensitive to its actual value.
IV.4 Entropy
One of the most natural ways to capture the effects of decoherence is to compute the entropy of an open system. In the present context this is the entanglement entropy between the field \chi (or rather its modes labelled by k) and the environment field \psi, S^{k}_{\text{ent}}[\chi,\psi]. This entanglement entropy is, of course, identified with the von Neumann entropy of the reduced density matrix \rho of the \chi field S_{\text{vN}}=\text{Tr}(\rho\log\rho). The growth of entropy (Fig. 4) provides an important test for the emergence of classicality and the process of decoherence as the system evolves from an initially pure to a mixed state.
Given that the state factorizes \rho=\rho_{R}\rho_{I}, as in (79), the total entropy is the sum of identical contributions from \rho_{R} and \rho_{I}, i.e.
S^{k}_{\text{ent}}[\chi,\psi]\equiv S_{\text{vN}}(\rho)=S(\rho_{I})+S(\rho_{R}% )=2S(\rho_{R}).  (99) 
where
\displaystyle S(\rho_{R})  \displaystyle=\sum_{n=0}^{\infty}p_{n}\log p_{n},  (100) 
with p_{n} given by eqs. (85) and (86). We find
\begin{split}\displaystyle S(\rho_{R})&\displaystyle=\log\frac{2}{\sqrt{1+\xi% /\Omega_{R}}+1}\\ &\displaystyle\frac{\sqrt{1+\xi/\Omega_{R}}1}{2}\log\frac{\sqrt{1+\xi/\Omega% _{R}}1}{\sqrt{1+\xi/\Omega_{R}}+1},\end{split}  (101) 
with the corresponding plot shown in Fig. 4.
Notice that S(\rho) depends on the ratio \xi/\Omega_{R}, so that the emergence of classical stochastic behaviour (as the state becomes more mixed), as characterised by an increase in entropy, depends on the relative strength of squeezing and decoherence. In the approximation of KP (); KPS2000 (); Kiefer:2006je (), \xi was assumed to be constant, and so this ratio grows as \xi/\Omega_{R}\sim\eta^{2}, whereas in reality, our analysis reveals it scales as \eta^{3} due to the timedependence of \xi\sim\eta^{1}. Hence our analysis reveals that the growth in entropy is even faster. The same is true for the Wigner area discussed below. This demonstrates the importance of a proper analysis of the evolution of the density matrix in realistic systems.
IV.5 Wigner function and phase space
The tradeoff between squeezing and decoherence can be better understood by exploring the phase space portrait of the state as provided by the Wigner function. Explicitly, by performing the integral in eq. (61) with the Gaussian ansatz (66), we get
\begin{split}\displaystyle W(\chi,\Pi)=\frac{4}{\pi^{2}}\frac{\Omega_{R}}{\xi+% \Omega_{R}}\exp\Big{[}\frac{2\left\Pi+\Omega_{I}\chi\right^{2}}{\Omega_{R}+% \xi}2\Omega_{R}\left\chi\right^{2}\Big{]}.\end{split}  (102) 
Notice that this factorises into two identical Wigner functions for W(\chi_{R},\Pi_{R}) and W(\chi_{I},\Pi_{I}). The evolution of W(\chi_{R},\Pi_{R}) is plotted in figures 5 and 6.
Using this, one can compute the spread in the “position” and “momentum” of the field as captured by the correlators
\braket{\chi_{R}^{2}}=\frac{1}{4\Omega_{R}},\qquad\braket{\Pi_{R}^{2}}=\frac{1% }{4\Omega_{R}}(\Omega^{2}+\xi\Omega_{R}).  (103) 
Notice that since we set k=1, we are working in a set of coordinates in which \chi and \Pi effectively have the same dimension. To restore dimensions one simply reinserts factors of k. The size of the area in phase space explored by the system is characterised by the Wigner ellipse defined by the contour (dropping the R subscript)
1=\frac{2\left(\Pi+\Omega_{I}\chi\right)^{2}}{\Omega_{R}+\xi}+2\Omega_{R}\chi^% {2},  (104) 
This describes an ellipse with a “tilt” of angle \alpha relative to the \chi axis shown in figure 7. The ellipse is aligned with a rotated set of axes \chi^{\prime} and \Pi^{\prime}, in which its equation is
\frac{\chi^{\prime 2}}{a^{2}}+\frac{\Pi^{\prime 2}}{b^{2}}=1,  (105) 
with semimajor and minor axes a and b, respectively. The extreme squeezing of the ellipse amounts gives a onetoone relation between the values of \Pi and \chi of the form \Pi\simeq\tan\alpha\chi. The two sets of coordinates are related by a rotation through an angle \alpha:
\left(\begin{array}[]{c}\chi^{\prime}\\ \Pi^{\prime}\end{array}\right)=\left(\begin{array}[]{cc}\cos\alpha&\sin\alpha% \\ \sin\alpha&\cos\alpha\end{array}\right)\left(\begin{array}[]{c}\chi\\ \Pi\end{array}\right).  (106) 
Inserting this expression into (105) and comparing with (104) gives, after a little algebra,
\begin{split}&\displaystyle a^{2},b^{2}=\\ &\displaystyle\frac{\Omega_{R}+\xi}{\Omega^{2}+\Omega_{R}\xi+1\mp\sqrt{(% \Omega^{2}+\Omega_{R}1)^{2}+4\Omega_{I}^{2}}},\end{split}  (107) 
and a rotation angle
\tan 2\alpha=\frac{2\Omega_{R}}{1\Omega^{2}\Omega_{R}\xi}.  (108) 
The area of this ellipse is given by
\mathcal{A}=\pi ab=\frac{\pi}{2}\left(1+\frac{\xi}{\Omega_{R}}\right)^{1/2}  (109) 
From fig. 7, we see that the minor axis is larger in the interacting case due to the system diffusing in phase space. This is reflected by a growth in the Wigner area (fig. 8), whose size is determined by the ratio \xi/\Omega_{R}. The growth in the Wigner area is characteristic of open systems and indicates that the state is becoming more “classical”. Indeed, the inverse of the area is roughly the product of the position and momentum correlation lengths and increasing area means that the product of the correlation lengths is becoming smaller than \hbar(=1). In the associated classical stochastic system, the widening of the ellipse is due to diffusion created by the noise in the Langevin equation (96). The broadening of the Wigner ellipse is also reflected in the momentum correlator \braket{\Pi_{R}^{2}} plotted in fig. 9. Indeed from eq. (103) we see the main contribution to \braket{\Pi_{R}^{2}} comes from \xi on superhorizon scales. Thus the broadening of the ellipse corresponds to local effects. By contrast, the length of the ellipse, which characterizes the size of \braket{\chi^{2}}, appears to be roughly the same in the free and interacting cases. The reason for this is that the corrections to \braket{\chi^{2}} are derived from nonlocal (memory) contributions, and are therefore much weaker in comparison the localcontributions driving the momentum diffusion and broadening of the ellipse.
V Sensitivity of discord to decoherence
Entanglement is a defining feature of quantum systems and is a subtle form of correlation that cannot be mimicked by a classical system. However, it is not an indication of nonlocal interaction: correlation is not interaction! So entanglement, and its ensuing correlations, are the calling card of a quantum system. Given the notorious sensitivity of entanglement to environmental noise, it is of particular interest to investigate the question of how robust the entanglement created during inflation is to the effects of decoherence. Examining this question will be the objective of the remainder of this paper.
In the present context, by writing the Hamiltonian (58) in terms of creation and annihilation operators as
H_{0}=\frac{k}{2}\left[\,a_{{\bm{k}}}a_{{\bm{k}}}^{\dagger}+a_{{\bm{k}}}a_{{% \bm{k}}}^{\dagger}\right]+i\frac{a^{\prime}}{a}\left[a_{{\bm{k}}}^{\dagger}a_{% {\bm{k}}}^{\dagger}+a_{{\bm{k}}}a_{{\bm{k}}}\right],  (110) 
it becomes apparent that the squeezing term involving a^{\prime}/a leads to the creation of entangled particle pairs with opposite momenta. These field quanta are then separated outside the horizon leading to spacelike field correlations on superHubble scales, whose strength is related to the expansion scale H. The purpose of the last part of this paper is to examine how quantum correlations are affected by decoherence, which intuitively one expects to render an initially quantum system more classical.
V.1 Overview of quantum discord
The entanglement entropy is the most familiar way to capture the strength of quantum correlations between two regions A and B of a bipartite system, and provides a simple way of quantifying the degree of entanglement within a system. However, the entanglement entropy only provides such a measure when the state of the total system is pure. The question is how can one measure quantum entanglement between A and B when the total state is not pure. In recent years, a new measure – quantum discord OZ (); HV () – has been proposed as a means of capturing quantum correlations in just such a situation. Quantum discord exploits the difference in classical and quantum posterior probabilities and the discrepancy, at the quantum level, between two classically equivalent definitions of mutual information, defined for classical systems as
I(A:B)=H(A)+H(B)H(A,B),  (111) 
where H is the Shannon entropy. Classically, this is equivalent to
J(A:B)=H(A)H(AB),  (112) 
where H(AB) is the Shannon entropy conditional on B defined by
H(AB)=\sum_{b}P(B=b)H(AB=b),  (113) 
with
H(AB=b)=\sum_{a}p(A=aB=b)\log p\left(A=aB=b\right).  (114) 
The equivalence of these two quantities is ensured by Bayes’ theorem
\begin{split}\displaystyle p(AB=b)=\frac{p(A,B=b)}{p(B=b)}.\end{split}  (115) 
However, at the quantum level, defining a posterior probability is more subtle. Indeed, the act of knowing and measuring the state of B may affect the state of the whole system \rho_{AB}! It is this subtlety which lies at the heart of quantum discord. Instead, in the quantum picture, one must first define a notion of conditional entropy, and construct a conditional density matrix. To do this, we must choose a positive operator valued measure (POVM) \left\{\Pi^{B}_{n}\right\} – this provides a set of measurement outcomes for the system B on which to condition the state of A. After a measurement outcome of \Pi_{n} on B the state is projected according to
\rho_{AB}\rightarrow\frac{\Pi^{B}_{n}\rho_{AB}\Pi^{B}_{n}}{P_{n}},  (116) 
where
P_{n}=\text{Tr}_{AB}\rho_{AB}\Pi^{B}_{n}  (117) 
is the probability that the measurement outcome \Pi^{B}_{n} occurs and ensures the correct normalisation of \rho_{AB=\Pi^{B}_{i}}. Hence the reduced density matrix of A, conditional on the measurement outcome \Pi_{n}^{B}, is
\rho_{AB=\Pi^{B}_{i}}=\text{Tr}_{B}\frac{\rho_{AB}\Pi^{B}_{n}}{P^{B}_{n}},  (118) 
where we assumed that the POVM is projective in the sense that \Pi_{n}^{2}=\Pi_{n}. We can then define a conditional entropy over all possible measurement outcomes in B:
S(AB=\{\Pi^{B}_{n}\})=\sum_{n}P_{n}S(\rho_{AB=\Pi^{B}_{n}}).  (119) 
With this translation to the quantum picture, the natural definition of J(A:B)_{\Pi^{B}_{n}} becomes
\mathcal{J}(A:B)=S(A)S(AB=\{\Pi^{B}_{n}\}).  (120) 
In fact, this is the Holevo information from A to B ZZ (). The generalisation of the mutual information is more straightforward and consists simply of replacing Shannon by von Neumann entropy:
\mathcal{I}(A:B)=S(A)+S(B)S(A,B).  (121) 
The quantum discord is then defined by the difference
\begin{split}\displaystyle\delta(A:B)_{\Pi_{n}^{B}}&\displaystyle=\mathcal{I}(% A:B)\mathcal{J}(A:B)_{\Pi^{B}_{n}}\\ &\displaystyle=S(B)S(A,B)+\sum_{i}P_{i}S(\rho_{AB=\Pi_{n}}).\end{split}  (122) 
Notice the appearance of the POVM \{\Pi_{n}\} in the definition, which reminds us how conditional probability is intimately related to measurement at the quantum level.
Note that the discord depends implicitly on the choice of POVM and so one can define the discord by minimizing over the choice of POVM. In general, this minimizing is difficult to perform and below we shall content ourselves with calculating a discord for a particular physically motivated choice of POVM. Our POVM dependent discord is therefore only an upper bound on the discord.
V.2 Decoherence and Discord in inflation
Recently, discord has been applied in the context of cosmological inflation in Lim () to probe correlations between primordial fluctuations and environment degrees of freedom and in VM () to probe correlations amongst an isolated system of primordial fluctuations. Our aim is to build on the approach of VM () and ask what happens to these correlations once the system is exposed to an external decohering environment. This is illustrated by the cartoon in Fig. 10. We shall follow the approach of VM () and choose the most natural partitioning of the (open) system by splitting the Hilbert space into +{\bm{k}} and {\bm{k}} modes. Physically, this corresponds to studying entanglement between field quanta created backtoback by the de Sitter expansion as they journey outside the horizon and experience decoherence. Since we have an additional environment, we are considering a tripartition of the overall system: the +{\bm{k}} mode, {\bm{k}} mode and environmental field \psi.
We begin by writing the density matrix in a form which makes this partitioning more transparent, expressing it in terms of the occupation number basis as
{\rho}(+{\bm{k}},{\bm{k}})=\sum_{n,m}C_{n,m,n^{\prime},m^{\prime}}\ket{n_{% \bm{k}},m_{{\bm{k}}}}\bra{n^{\prime}_{\bm{k}},m^{\prime}_{{\bm{k}}}},  (123) 
where C_{n,mn^{\prime},m^{\prime}}=\bra{n_{\bm{k}},m_{{\bm{k}}}}\rho\ket{n^{\prime% }_{\bm{k}},m^{\prime}_{{\bm{k}}}} are timedependent coefficients. We must now decide on a choice POVM in order to compute the conditional density matrix (118). We follow VM () and choose \Pi^{{\bm{k}}}_{n}=\ket{n_{{\bm{k}}}}\bra{n_{{\bm{k}}}}, which leads to
\begin{split}\displaystyle\rho_{AB=n_{{\bm{k}}}}=\frac{1}{p_{n}}\text{Tr}_{% {\bm{k}}}\left[\rho\Pi^{{\bm{k}}}_{n}\right]=\frac{1}{p_{n}}\sum_{m}p_{nm}% \ket{m_{{\bm{k}}}}\bra{m_{{\bm{k}}}},\end{split}  (124) 
where
p_{nm}=\bra{m_{\bm{k}}n_{{\bm{k}}}}\rho\ket{m_{\bm{k}}n_{{\bm{k}}}},  (125) 
and, in analogy with eq. (117),
\begin{split}\displaystyle p_{n}&\displaystyle=\text{Tr}_{(+{\bm{k}},{\bm{k}}% )}\left[\rho\,\Pi^{{\bm{k}}}_{n}\right]\\ &\displaystyle=\sum_{m}\bra{m_{\bm{k}}n_{{\bm{k}}}}\rho\ket{m_{\bm{k}}n_{{% \bm{k}}}}\equiv\sum_{m}p_{nm}.\end{split}  (126) 
After a series of somewhat lengthy calculations described in appendix A, we find
\displaystyle p_{n}=4\Omega_{R}\frac{\!\!\!\left[12\Omega_{R}+\Omega^{2}+% \xi\Omega_{R}\right]^{m}}{\,\,\left[1+2\Omega_{R}+\Omega^{2}+\xi\Omega_{R}% \right]^{m+1}},  (127) 
\begin{split}&\displaystyle p_{nm}=\frac{2(m+n)!}{m!\,n!}\frac{2\Omega_{R}\,% \xi^{m+n}}{\left[1+2\Omega_{R}+\Omega^{2}+\xi\Omega_{R}+\xi\right]^{m+n+1}}% \\ &\displaystyle\times\,\,\tensor{[}_{2}]{F}{{}_{1}}\Big{[}\!m,\!n,\!m\!\!n,% \,\,\frac{1}{\xi^{2}}\left(12\Omega_{R}+\Omega^{2}+\xi\Omega_{R}\xi\right% )\left(1+2\Omega_{R}+\Omega^{2}+\xi\Omega_{R}+\xi\right)\Big{]},\end{split}  (128) 
where we remember that earlier in sec. IV.3 we switched to units where k=1. To restore dimensionality to (127) and (128) one simply reinserts appropriate factors of k. Notice there is a nonzero probability for states to have different occupation numbers in the \pm{\bm{k}} modes, which cannot happen in the noninteracting case as we begin in the vacuum and particles are created only in pairs with opposite momenta. However, the \chi\psi^{2} interaction violates \chi number in such a way as to allow differing numbers in the \pm{\bm{k}} modes. The state of the reduced density \rho({\bm{k}}) is
\rho({\bm{k}})=\sum_{n}p_{n}\ket{n_{{\bm{k}}}}\bra{n_{{\bm{k}}}}.  (129) 
This is a thermal state with temperature
\beta_{k}=\log\left[\frac{12\Omega_{R}+\Omega^{2}+\xi\Omega_{R}}{1+2\Omega% _{R}+\Omega^{2}+\xi\Omega_{R}}\right].  (130) 
The corresponding entropy S({\bm{k}}) is given by
S({\bm{k}})=(1+\braket{n_{k}})\log(1+\braket{n_{k}})\braket{n_{k}}\log\braket% {n_{k}}  (131) 
where
\braket{n_{\bm{k}}}=\frac{1}{e^{\beta_{\bm{k}}}1},  (132) 
is nothing more than the BoseEinstein distribution.
Returning to the conditional entropy (124) we have
\displaystyle S\left(AB=\{\Pi^{n}_{{\bm{k}}}\}\right)  \displaystyle=\sum^{\infty}_{m=0}\frac{p_{mn}}{p_{n}}\log\Big{(}\frac{p_{mn}}% {p_{n}}\Big{)}.  (133) 
After a short calculation we find
\mathcal{J}=2\sum^{\infty}_{n=0}p_{n}\log p_{n}+\sum^{\infty}_{m,n=0}p_{mn}% \log p_{mn}.  (134) 
On the other hand S({\bm{k}})=S({\bm{k}})=\sum_{n}p_{n}\log p_{n} and S({\bm{k}},{\bm{k}})\equiv S^{k}_{\text{ent}}[\chi,\psi] is identified with the entanglement entropy of the modes of the field \chi and \psi calculated in sec. IV.4; hence
\begin{split}\displaystyle{\cal I}=2\sum_{n=1}^{\infty}p_{n}\log p_{n}S^{k}_% {\text{ent}}[\chi,\psi]\end{split}  (135) 
and then (121) and (122) give the POVMdependent discord
\begin{split}\displaystyle\delta&\displaystyle\equiv\delta({\bm{k}}:{\bm{k}})% _{\Pi_{n}^{{\bm{k}}}}\\ &\displaystyle=S^{k}_{\text{ent}}[\chi,\psi]\sum^{\infty}_{m,n=0}p_{mn}\log p% _{mn}.\end{split}  (136) 
For the noninteracting case we find that p_{nm}=\delta_{nm}p_{n}, which leads to
\delta_{\text{free}}=\sum_{n}p^{{\text{free}}}_{n}\log p_{n}^{{\text{free}}}% \equiv S_{\text{free}}({\bm{k}}).  (137) 
This is nothing more than the entanglement entropy of \rho_{\text{free}}({\bm{k}}). Indeed, it is a wellknown fact (see e.g. AD ()) that the discord of a bipartite division of a system in a pure state is equal to the entanglement entropy. On the other hand, when there is interaction with a third system – in this case the \psi field – it is not possible to define the entanglement entropy between the +{\bm{k}} and {\bm{k}} subsystems and then one needs quantities like the discord to provide a measure of the quantum correlations between the +{\bm{k}} and {\bm{k}} sectors. We should emphasize again that what we have calculated is a discord that depends upon the choice of POVM, albeit that it is a particularly natural choice, and as such it provides an upper bound on the true discord. Nevertheless it does provide a measure of sorts of the quantum correlation between our sub systems.
The discord is plotted in Fig. 11. Note that the “free discord” (which is really just the entanglement entropy) dominates the discord in the interacting case, meaning that the \mathcal{O}\left(\lambda^{2}/H^{2}\right) corrections to the discord cannot be seen in fig. 11. The difference between the two is seen more clearly in fig. 12. Thus, as expected, we see decoherence does have a tendency to weaken discord, though the effect is perturbatively small in our case. One can see this more explicitly by counting powers of \lambda. Firstly, we note that since the entropy of a pure state is zero, it follows that
\begin{split}\displaystyle S^{k}_{\text{ent}}[\chi,\psi]=\mathcal{O}\left(% \lambda^{2}/H^{2}\right).\end{split}  (138) 
Furthermore, we have
p_{mn}=\begin{cases}p^{\text{free}}_{n}+\mathcal{O}\left(\lambda^{2}/H^{2}% \right),&\text{if}\ m=n\\ \mathcal{O}\left(\lambda^{2}/H^{2}\right),&\text{if}\ m\neq n\end{cases}.  (139) 
Therefore, from eqs. (122) and (137) we see
\displaystyle\delta  \displaystyle=\delta_{\text{free}}+\mathcal{O}\left(\lambda^{2}/H^{2}\right).  (140) 
It has already been noted in other contexts that discord is weakened by the presence of an environment. For example, in BA (), the effect of various quantum channels acting on spin chains was considered, and it was found that discord in the system was not only weakened by decoherence but actually decayed. This confirms the intuition, that, in a broad sense, decoherence should reduce discord.
For our setup, the discord does not decrease in time, but its growth rate is slowed relative to the free case as can be seen from fig. 12.
An interpretation of this might go as follows. In the free case, the discord is just the entanglement entropy between the +{\bm{k}} and {\bm{k}} modes. It grows with time as more and more entangled pairs are created, leading the +{\bm{k}} and {\bm{k}} sectors to become strongly correlated as time goes on. If one measures n particles in the {\bm{k}} state, this pairrule ensures that there have to be n particles in state +{\bm{k}}.
However, in the interacting case, the \chi\psi^{2} coupling allows the numbers in the +{\bm{k}} and {\bm{k}} modes to fluctuate independently beyond pair production due to \psi\psi\leftrightarrow\chi scattering processes. Therefore, if a measurement of the +{\bm{k}} mode reveals an occupation number n, the state of the {\bm{k}} mode need not contain the same number of particles. This is reflected by the fact that when interactions are switched on, p_{mn} can be nonzero for n\neq m. Hence the occupation numbers in the \pm{\bm{k}} modes are more weakly correlated in the interacting case due to random environmentinduced fluctuations.
Thus two processes – pairproduction and environmental scattering – are competing to increase and decrease the discord respectively. In the present setup, pairproduction wins, leading to an increasing discord albeit at a slower rate in the interacting theory.
One can gain a little insight into this last point with some power counting. The pair production rate of a massless field on dimensional grounds must be of order H. Similarly, the scattering rate between a massless field \chi and a conformally coupled field \psi is order \lambda. Therefore we expect the correlations between +{\bm{k}} and {\bm{k}} sectors to increase if the pair production rate is greater than the scattering rate, i.e., if 1>\lambda/H, but this is precisely the perturbativity condition. Hence for 1>\lambda/H we expect the discord to increase. In the \lambda>H case we might expect the opposite, but that would lie beyond the perturbation theory used in our present analysis. One might question how these conclusions would be affected in the case of a different interaction Hamiltonian. We leave such questions for future work.
It is important to realise that, in general, many quantum features – including entanglement – are not robust to the effects of decoherence. In typical “laboratory” setups, even the weakest of interactions lead to rapid decoherence and the decay of discord, as shown in BA (). Therefore, it is far from obvious that quantum correlations should survive the effects of environmental noise. This reveals the importance of the full analysis of decoherence, discord and entanglement presented here. The fact that quantum correlations do survive in the present case is due to the continual production of entanglement by cosmological expansion which counteracts the effects of decoherence. It is worth noting that such a mechanism is quite novel to the cosmological scenario.
VI Discussion
In this work we studied the process of decoherence in inflation. By considering an explicit interaction between scalar perturbations and a conformally coupled environment field, we derived the master equation for the reduced density matrix of perturbations first presented in Boy1 (). We have provided a deeper analysis and interpretation of the master equation by explicitly solving for the density matrix with a Gaussian ansatz and Bunch Davies initial conditions. This revealed a scaledependent nature of decoherence, with the rate at which offdiagonal elements of the density matrix are driven to zero growing in strength on superHubble scales. We also found that nonlocal contributions (i.e. memory effects) from environmental correlators have a negligible effect on the superhorizon strength of decoherence which is dominated by local (white noise) from the environment. By contrast, the Gaussian width receives contributions from nonlocal memory effects due to longtime environment correlations, which feed into the power spectrum. This confirms from a different perspective the power spectrum calculation of Boyanovksy Boyanovsky:2015jen (). We also showed that potential positivityviolating effects are suppressed on superHubble scales leading to a well behaved density matrix for modes that cross the horizon. This allowed us to construct a set of classical stochastic equations. The resulting Langevin and FokkerPlanck equations provided an emergent classical description of the same dynamics as the quantum master equation on superHubble scales. It should be noted that these are due to genuine interaction with an environment due to direct couplings between different Fourier modes of the system and environment field and are therefore quite different to the effective stochastic equations encountered in the usual stochastic inflation paradigm Star1 (); SY ().
In the remainder of the paper we used these tools to examine the question of the quantumtoclassical transition and the emergence of classicality from an initially pure quantum state. The most natural quantity to consider was the entropy, which grew in time as decoherence took effect. In particular we saw that entropy was determined solely by the ratio \xi/\Omega_{R}, so that the extent to which the state became mixed was controlled by the relative strength of squeezing and decoherence. We also examined the Wigner function and found that the area of the Wigner ellipse \mathcal{A}=\pi/2(1+\xi/\Omega_{R})^{1/2} can be captured by the same ratio \xi/\Omega_{R}, revealing that decoherence has a tendency to increase the size of phasespace explored by the system due to diffusion. Remarkably, however, the power spectrum itself was only affected in a very weak way via nonlocal memory effects.
Finally, we turned our attention to the question of quantum correlations in the system and how robust they are to environmental noise. We used quantum discord as a measure of the strength of quantum correlations by partitioning the twomode density matrix into \pm{\bm{k}} modes. Physically, this corresponds to the correlations between particle pairs created backtoback by the de Sitter expansion. It is these correlations which lead to the spacelike field correlations \braket{\chi_{\bm{k}}\chi_{{\bm{k}}}} on large scales. By looking at the discord between these modes, we were able to see the effect of decoherence on the strength of quantum correlations between \pm{\bm{k}} modes.
For isolated states, pair creation induces strong correlations between the \pm{\bm{k}} sectors, but environmental interactions lead to local fluctuations in the numbers in the \ket{n_{\bm{k}}} and \ket{n_{{\bm{k}}}} states, reducing the strength of correlations between the two. As a result, decoherence reduces the strength of quantum discord, and in this particular sense, renders the system more classical. However, these corrections were of \mathcal{O}(\lambda^{2}/H^{2}), which for our weakly interacting system was perturbatively small and discord remained dominated by the entanglement due to pairproduction. This suggests that for weakly interacting systems, at least in the pure de Sitter phase of expansion, pair production dominates over decoherence effects provided the reaction rate with the environment is less than the pairproduction rate characterised by H. For this reason it is important to reproduce our analysis of discord and decoherence with other interactions.
There are various avenues for developing this work. Firstly, given the timedependent nature of the interaction Hamiltonian, perturbativity breaks down at late times, meaning that one can only push the analysis used here and in Boy1 () so far outside the Horizon. One resolution to this problem might be to consider other kinds of interactions, e.g. a marginal coupling of the form \lambda\int d^{4}x\sqrt{g}\phi\varphi^{3}=\lambda\int d^{4}x\chi\psi^{3}, which does not suffer from the same scaling. This would also provide an interesting way to probe the modeldependence of the growth in decoherence strength observed for the particular interaction considered here.
In this paper we have only considered pure de Sitter expansion, rather than a complete slowroll analysis which allows inflation to end. Therefore, a more developed analysis of decoherence should include a generic inflationary potential which would allow one to see how decoherence changes as one exists the inflationary era. In particular, one might expect that once the rate of particle pair creation has dropped after inflation, decoherence may persist and wipe out any quantum discord.
Finally, there are many parallels with pair creation in black hole physics. Indeed, both de Sitter and black hole spacetimes possess horizons and exhibit Hawking radiation, which is described by the reduced density matrix after tracing out the remaining member of the pair – see our eq. (130). Hence the analysis presented here might provide an interesting springboard for studying decoherence and open quantumsystems in the context of black holes.
Acknowledgments
JIM would like to thank Gianmassimo Tasinato for first sparking his interest in the question of open systems in inflation, Markus Müller for many useful discussions on quantum information, quantum optics and decoherence and Davide De Boni for helpful conversations about Wigner functions, squeezing and decoherence. This work was supported by the STFC grants ST/K502376/1 and ST/L000369/1.
Appendix A Discord calculations
We begin by describing how to compute the matrix elements of the reduced density matrix \rho(+{\bm{k}}) in the number basis, defined by
\displaystyle\bra{n}\rho(+{\bm{k}})\ket{m}=\sum_{\ell}\bra{n,\ell}\rho\ket{m,\ell}  (141) 
where we have suppressed the {\bm{k}} subscripts. Inserting a complete set of position eigenstates on the right hand side gives
\displaystyle\bra{n}\rho(+{\bm{k}})\ket{m}=\int d^{2}x\int d^{2}y\,\ \rho(y,x)% \sum_{\ell}\psi_{m\ell}(x)\psi_{n\ell}^{*}(y).  (142) 
where \psi_{ml}(x)=\braket{x}{m\ell} is the position space wavefunction of the number state. This can be written in the usual way, in terms of ladder operators acting on the ground state wavefunction, i.e.
\psi_{m,\ell}(x)=\frac{\left(a^{\dagger}\right)^{m}}{\sqrt{m!}}\frac{\left(% \bar{a}^{\dagger}\right)^{\ell}}{\sqrt{\ell!}}\sqrt{\frac{2}{\pi}}e^{x^{2}}  (143) 
where the action of the ladder operators in position space is given by
a^{\dagger}=\frac{x\bar{\partial}}{\sqrt{2}},\qquad a=\frac{x+\bar{\partial}}% {\sqrt{2}}  (144) 
and the bar on the ladder operator denotes complex conjugation, and comes from the ladder operators associated to the {\bm{k}} modes. Notice that since k will drop out of any quantities at the end of the calculation, we have effectively set k=1 from the outset. From this it follows that
\psi_{m\ell}(x)\psi^{*}_{n\ell}(y)=\frac{2}{\pi}\frac{1}{\ell!}\frac{1}{\sqrt{% n!}\sqrt{m!}}(a^{\dagger}_{x})^{m}(\bar{a}^{\dagger}_{x})^{\ell}(\bar{a}^{% \dagger}_{y})^{n}(a_{y}^{\dagger})^{\ell}e^{x^{2}y^{2}}.  (145) 
Noting that
(\bar{a}^{\dagger}_{x}a_{y}^{\dagger})^{\ell}e^{x^{2}y^{2}}=(\sqrt{2}% \bar{x}\sqrt{2}y)^{\ell}e^{x^{2}y^{2}},  (146) 
and inserting this expression into the sum over \ell gives
\sum_{\ell}\psi_{m\ell}(x)\psi_{n\ell}^{*}(y)=\frac{2}{\pi}\frac{1}{\sqrt{m!}% \sqrt{n!}}(a^{\dagger}_{x})^{m}(\bar{a}^{\dagger}_{y})^{n}e^{x^{2}y^{2}+% 2\bar{x}y}\ .  (147) 
Next one has to consider the action of the remaining derivatives. Using \bar{a}_{y}=(\bar{y}\partial_{y})/\sqrt{2} gives
\begin{split}\displaystyle\sum_{\ell}\psi_{m\ell}(x)\psi_{n\ell}^{*}(y)&% \displaystyle=\frac{2}{\pi}\frac{1}{\sqrt{m!}\sqrt{n!}}(a^{\dagger}_{x})^{m}% \left[\sqrt{2}\left(\bar{y}\bar{x}\right)\right]^{n}e^{x^{2}y^{2}+2\bar% {x}y}\\ &\displaystyle=\frac{1}{\sqrt{m!}\sqrt{n!}}(a^{\dagger}_{x})^{m}\left.\frac{d^% {n}}{d\lambda^{n}}e^{x^{2}y^{2}+2\bar{x}y+\sqrt{2}\lambda(\bar{y}\bar{x% })}\right_{\lambda=0}.\end{split}  (148) 
We then act once again with a_{x}^{\dagger}=(x\bar{\partial}_{x})/\sqrt{2}, giving rise to
\sum_{\ell}\psi_{m\ell}(x)\psi_{n\ell}^{*}(y)=\frac{2}{\pi}\frac{1}{\sqrt{m!}% \sqrt{n!}}\frac{d^{n}}{d\lambda^{n}}\left\{\left[\sqrt{2}(xy)+\lambda\right]^% {m}e^{x^{2}y^{2}+2\bar{x}y+\sqrt{2}\lambda(\bar{y}\bar{x})}\right\}_{% \lambda=0}.  (149) 
Now we invoke the Leibniz rule for generalised products and take the \lambda\rightarrow 0 limit giving
\sum_{\ell}\psi_{m\ell}(x)\psi_{n\ell}^{*}(y)=\frac{2}{\pi}\frac{1}{\sqrt{m!}% \sqrt{n!}}\sum^{n}_{k=0}\left(\begin{array}[]{c}n\\ k\end{array}\right)\frac{m!}{(mk)!}\left(\sqrt{2}(xy)\right)^{mk}\left(% \sqrt{2}\lambda(\bar{y}\bar{x}\right)^{nk}e^{x^{2}y^{2}+2\bar{x}y}\ .  (150) 
Inserting the matrix elements \rho(y,x) gives
\begin{split}\displaystyle\bra{n}\rho(+{\bm{k}})\ket{m}&\displaystyle=\frac{2}% {\pi}\frac{1}{\sqrt{m!}\sqrt{n!}}\sum^{n}_{k=0}\frac{n!m!}{(nk)!(mk)!}\frac{% \sqrt{2}^{m+n2k}}{k!}\frac{2\Omega_{R}}{\pi}\\ &\displaystyle\times\int d^{2}x\int d^{2}y\,\,\left(xy\right)^{mk}\left(\bar% {y}\bar{x}\right)^{nk}\exp\Big{(}(1+\Omega)y^{2}(1+\Omega^{*})x^{2}% \frac{\xi}{2}xy^{2}+2\bar{x}y\Big{)}.\end{split}  (151) 
The (complex) Gaussian integral is only nonvanishing for m=n which gives
\begin{split}&\displaystyle\bra{n}\rho(+{\bm{k}})\ket{m}=\delta_{m,n}\frac{4% \Omega_{R}}{\pi^{2}}\sum^{n}_{k=0}\frac{n!}{[(nk)!]^{2}}\frac{2^{nk}}{k!}\\ &\displaystyle\times\int d^{2}x\int d^{2}y\,\,\left(x\bar{y}+y\bar{x}x^{2}% y^{2}\right)^{nk}\exp\Big{[}(1+\Omega)y^{2}(1+\Omega^{*})x^{2}\frac{% \xi}{2}xy^{2}+2\bar{x}y\Big{]}.\end{split}  (152) 
We can write the Gaussian moments in terms of derivatives, so that
\begin{split}&\displaystyle\bra{n}\rho(+{\bm{k}})\ket{m}=\delta_{mn}\,\frac{4% \Omega_{R}}{\pi^{2}}\sum^{n}_{k=0}\frac{n!}{[(nk)!]^{2}}\frac{2^{nk}}{k!}\\ &\displaystyle\times\left.\frac{d^{nk}}{d\,z^{nk}}\int d^{2}x\int d^{2}y\,% \exp\Big{[}(1+\Omega+z)y^{2}(1+\Omega^{*}+z)x^{2}\frac{\xi}{2}xy^{2}% +(2+z)\bar{x}y+z\bar{y}x\Big{]}\right_{z=0}.\end{split}  (153) 
Performing the multidimensional complex Gaussian integral gives a contribution \pi^{2}/\det A, where A is the covariance matrix appearing in the exponent. Its determinant is given by
\begin{split}\displaystyle\det A=1+2\Omega_{R}+\Omega^{2}+\xi\Omega_{R}+z2% \Omega_{R}\equiv\alpha+\beta z.\end{split}  (154) 
From this it follows that
\begin{split}\displaystyle\bra{n}\rho(+{\bm{k}})\ket{m}&\displaystyle=\delta_{% mn}4\Omega_{R}\sum^{n}_{k=0}\frac{n!}{[(nk)!]^{2}}\frac{2^{nk}}{k!}\left.% \frac{d^{nk}}{d\,z^{nk}}\left[\alpha+\beta z\right]^{1}\right_{z=0}\\ &\displaystyle=\delta_{m,n}4\Omega_{R}\sum^{n}_{k=0}\frac{n!}{(nk)!}\frac{2^{% nk}}{k!}(1)^{nk}(\beta)^{nk}\alpha^{1(nk)}\ .\end{split}  (155) 
Performing the binomial sum, we get
\bra{n}\rho(+{\bm{k}})\ket{m}=\delta_{mn}4\Omega_{R}\frac{\left[12\Omega_{R}+% \Omega^{2}+\xi\Omega_{R}\right]^{n}}{\left[1+2\Omega_{R}+\Omega^{2}+\xi% \Omega_{R}\right]^{n+1}}\equiv\delta_{mn}p_{n}\ .  (156) 
Notice this satisfies the correct normalisation \sum_{n}p_{n}=1. The procedure for computing p_{mn}=\bra{mn}\rho\ket{nm} is similar, leading, after a lengthy calculation, to
\begin{split}&\displaystyle p_{nm}=\frac{2(m+n)!}{m!\,n!}\frac{2\Omega_{R}\,% \xi^{m+n}}{\left[1+2\Omega_{R}+\Omega^{2}+\xi\Omega_{R}+\xi\right]^{m+n+1}}% \\ &\displaystyle\times\,\,\tensor{[}_{2}]{F}{{}_{1}}\Big{[}\!m,\!n,\!m\!\!n,% \,\,\frac{1}{\xi^{2}}\left(12\Omega_{R}+\Omega^{2}+\xi\Omega_{R}\xi\right% )\left(1+2\Omega_{R}+\Omega^{2}+\xi\Omega_{R}+\xi\right)\Big{]}.\end{split}  (157) 
References
 (1) A. H. Guth, Phys. Rev. D23, 347 (1981).
 (2) A. H. Guth and SY. Pi, Phys. Rev. Lett. 49, 1110 (1982).
 (3) P. Ade et al. (Planck Collaboration), Astron. Astrophys. 571, A1 (2014), [astroph/1303.5062].
 (4) R. Adam et al. (Planck) (2015), [astroph/1502.01582].
 (5) P. A. R. Ade et al. (Planck) (2015), [astroph/1502.01589].
 (6) P. A. R. Ade et al. (Planck) (2015), [astroph/1502.02114].
 (7) C. Kiefer, D. Polarski, Annalen Phys. 7 (1998) 137158 [grqc/9805014].
 (8) C. Kiefer and D. Polarski, Adv. Sci. Lett. 2, 164 (2009) [astroph/0810.0087].
 (9) J. Martin, V. Vennin and P. Peter, Phys. Rev. D 86, 103524 (2012) [hepth/1207.2086].
 (10) H. P. Breuer and F. Petruccione “The theory of open quantum systems”, OUP (2002).
 (11) C. Kiefer , D. Polarski, A. A. Starobinsky, Phys. Rev. D 62 (2000) 043518 [grqc/9910065].
 (12) C. Kiefer, I. Lohmar, D. Polarski and A. A. Starobinsky, Class. Quant. Grav. 24, 1699 (2007) [astroph/0610700].
 (13) R. Brandenberger, R. Laflamme, M. Mijić, Physica Scripta. T36, 265268, 1991.
 (14) C. Kiefer, D. Polarski, A. A. Starobinsky, Int. J. Mod. Phys. D7 (1998) 455462, [grqc/9802003].
 (15) D. Boyanovsky, New J. Phys. 17 (2015) no. 6, 063017, [hepph/1503.00156].
 (16) D. Boyanovsky, Phys. Rev. D 93, 043501 (2016) [astroph.CO/1511.06649 ].
 (17) M. A. Sakagami, Prog. Theor. Phys. 79, 442 (1988).
 (18) F. C. Lombardo and D. Lopez Nacir, Phys. Rev. D 72, 063506 (2005) [grqc/0506051].
 (19) D. Campo and R. Parentani, Phys. Rev. D 78, 065045 (2008) [hepth/0805.0424].
 (20) P. Martineau, Class. Quant. Grav. 24, 5817 (2007) [astroph/0601134].
 (21) C. P. Burgess, R. Holman and D. Hoover, Phys. Rev. D 77, 063534 (2008) doi:10.1103/PhysRevD.77.063534 [astroph/0601646].
 (22) D. Polarski and A. A. Starobinsky, Class. Quant. Grav. 13, 377 (1996) [grqc/9504030].
 (23) J. Lesgourgues, D. Polarski and A. A. Starobinsky, Nucl. Phys. B 497, 479 (1997) [grqc/9611019].
 (24) C. Kiefer, D. Polarski and A. A. Starobinsky, Int. J. Mod. Phys. D 7, 455 (1998) [grqc/9802003].
 (25) C. Kiefer, J. Lesgourgues, D. Polarski and A. A. Starobinsky, Class. Quant. Grav. 15, L67 (1998) [grqc/9806066].
 (26) A. Albrecht, P. Ferreira, M. Joyce and T. Prokopec, Phys. Rev. D 50, 4807 (1994) [astroph/9303001].
 (27) C. P. Burgess, R. Holman, G. Tasinato and M. Williams, JHEP 1503, 090 (2015) [hepth/1408.5002].
 (28) C. P. Burgess, R. Holman and G. Tasinato, JHEP 1601, 153 (2016) [grqc/1512.00169].
 (29) A. Matacz, Phys. Rev. D 55, 1860 (1997) [grqc/9604022].
 (30) E. Nelson, JCAP 1603 (2016) 022, [grqc/160103734].
 (31) H. Ollivier and W. H. Zurek, Phys. Rev. Lett. 88, 017901 (2001), [quantph/0105072].
 (32) L. Henderson and V. Vedral, Journal of Physics A: Mathematical and General 34, 6899 (2001), [quantph/0105028].
 (33) A. O. Caldeira, A. J. Leggett, Physica 121A, 587616 (1983).
 (34) J. J. Halliwell, J. Phys. A40, 12, [quantph/0607132].
 (35) A. A. Starobinsky, Lect.Notes Phys. 246 (1986) 107â126.
 (36) A. A. Starobinsky, J. Yokoyama, Phys. Rev. D50 (1994) 63576368 , [astroph/9407016].
 (37) E. A. Lim, Phys. Rev. D91, 083522 (2015) [hepth/14105508].
 (38) J. Martin, V. Vennin Phys. Rev. D93, 023505 (2016) [astroph/1510.04038].
 (39) I. Bose, A. K. Pal, Int. J. Mod. Phys. B27, 1345042 (2013), [quantph/12051300].
 (40) E. Joos, H. D. Zeh, Z. Phys. B59 223243 (1985)
 (41) W. H. Zurek, Phys. Scripta T 76, 186 (1998) [Acta Phys. Polon. B 29, 3689 (1998)] [quantph/9802054].
 (42) M.V. Berry, Phil. Trans. Soc. 287 237271 (1997)
 (43) T. Bunch, and P. Davies, Proc. R. Soc. London A360, 117 (1978). 34, 6899 (2001).
 (44) I. S. Gradshteyn, I.M. Ryzhik, “Table of Integrals, Series, and Products”, p837, Academic Press Inc. (1980).
 (45) M. Zwolak and W. H. Zurek, Scientific Report 3 1729 (2013), [quantph/13034659].
 (46) G. Adesso, A. Datta, Phys. Rev. Lett. 105, 030501 (2010), [quantph/10034979].