Zero-dimensional models for gravitational and scalar QED decoherence

We investigate the dynamics of two quantum mechanical oscillator system-bath toy models obtained by truncating to zero spatial dimensions linearized gravity coupled to a massive scalar field and scalar QED. The scalar-gravity toy model maps onto the phase damped oscillator, while the scalar QED toy model approximately maps onto an oscillator system subject to two-photon damping. The toy models provide potentially useful insights into solving for open system quantum dynamics relevant to the full scalar QED and weak gravitational field systems, in particular operational probes of the decoherence for initial scalar field system superposition states.


I. INTRODUCTION
The non-existence of macroscopic mass system quantum superposition states under everyday conditions is commonly understood to be due to interactions with the system environment; air molecules, photons, and internal system defects cause the rapid decoherence of position and energy superposition states into apparent mixtures of either/or alternatives that are indistinguishable from a classical statistical distribution [1][2][3]. By placing the system in ultrahigh vacuum, shielding it from external electromagnetic radiation, and cooling the system down to its ground state, quantum mechanics would in principle allow for macroscopic system superposition states to be prepared and measured. However, there is one environment that cannot be screened out-gravity, as expressed dynamically at the classical level through Einstein's field equations [4][5][6][7][8][9][10].
From a fundamental perspective, it is interesting to try to quantify the effect of the gravitational environment on macroscopic mass/energy superposition states; even if the predicted gravitationally induced decoherence times are much longer than for everyday, electromagnetic environments, having a good quantitative understanding of the former would allow us to place in principle, unavoidable bounds on the coherence times of macroscopic superposition states, and furthermore help point the way towards possible future experiments to probe the role of gravity in enforcing macroscopic classicality.
Under terrestrial or space-based laboratory conditions corresponding to weak spacetime curvature [11], it should be sufficient to work with linearized gravity [12], where the matter system-gravitational environment action is quadratic in metric deviations h µν from Minkowski spacetime η µν [= diag(−+++)]: g µν = η µν +κh µν , where κ = √ 32πG (with natural units = c = 1). Furthermore, modeling the matter system through quantum excitations of a massive scalar field φ, a "first-principles" starting point for investigating gravitational decoherence * qidong.xu.gr@dartmouth.edu † miles.p.blencowe@dartmouth.edu is the following action: where the system, environment, and interaction actions are respectively: and with T µν (φ) the scalar field energy-momentum tensor, U µνρσ (φ) a quadratic in φ tensor [13], and h = h µ µ . Quantization might then proceed through the derivation of a master equation for the density matrix of the scalar matter system, with the (assumed for simplicity) thermal gravitational environmental degrees of freedom traced out [4][5][6]. Alternatively, a quantum Langevin equation might be derived for the scalar matter field operator, again with the gravitational environment integrated out. One route to obtaining such effective equations is the closed time path integral approach, which is particularly suited to field systems [14].
However, as a coupled system-environment field theory with a non-quadratic interaction and a gauge symmetry (i.e., general coordinate invariance), the derivation of the quantum gravitational decoherence dynamics presents additional challenges beyond the usual systemenvironment models considered in non-relativistic quantum mechanics [15][16][17]. One challenge involves the necessity for making various approximations in order to solve for the reduced system dynamics. For example, in the usual open quantum systems analyses, it is assumed that the system+environment is initially in a product state, e.g., the system is in a superposition of two distinct wavepacket or energy states and the environment is arXiv:2005.02554v2 [quant-ph] 6 Dec 2022 in a thermal state. Such a product state can result in an initial "burst" of decoherence that depends on the upper cut-off physics of the environment, which in the case of gravity is unknown. Furthermore, a Born and possibly Markovian approximation is made [15][16][17], where the influence of the environment on the system is treated perturbatively to lowest non-trivial order, while the environment is assumed to respond rapidly relative to the system dynamics timescale.
Another challenge concerns requiring gauge invariance of the calculated decoherence rates for them to be meaningful, in particular when assuming a finite temperature environment that comprises gauge degrees of freedom (e.g., photons or gravitons) [16]. A common, direct approach [15][16][17] to obtaining decoherence rates for open quantum systems, either with or without gauge degrees of freedom, is to examine the time evolution of the offdiagonal matrix elements of the system density operator in the state basis of interest (e.g., energy eigenstates, position eigenstates etc.). However, the density operator is not a gauge invariant quantity.
A more consistent approach is to extract the decoherence rates through an operational procedure, i.e., involving an in principle measurement that can be ascribed to a particular expectation value of an observable. One such example is the particle detection number density in an atom interference set-up. A minimal way to get scalar matter field quanta initially in superposition states to interfere is by spatially trapping the field quanta in a three-dimensional harmonic confining potential [18]; the system action (2) is then supplemented by the term where Ω is the characteristic oscillation frequency and the potential center coincides with the spatial origin r = 0 in the rest frame of the confining potential. Referring to Fig. 1, we might then consider a thought experiment where an initial (t = 0) N -nucleon state corresponding to being in a superposition of two collective coherent states with coordinate parameters r = (x, 0, 0), x = x 01 > 0, x 02 < 0 (and hence superposition separation x 01 − x 02 ) undergoes gravitational decoherence. Once during every oscillation period, the coherent states in the superposition pass through each other in the region centered at x = 0, resulting in an interference pattern for the x-dependence of the local particle detection probability as indicated in the figure. A measure of the degree of coherence is the so-called "interferometric visibility", defined below in Eq. (30); a reduction in visibility over time is interpreted as a signature of gravitationally induced decoherence or dephasing.
In particular, suppose that we have a particle detector with center of mass worldline (t, r) in the vicinity of r = 0 and described by the local field observable V −1 V drφ(r, 0) 2 , where φ(r, 0) is the scalar field operator (in the Schrödinger picture), and V is a coordinate averaging volume (assumed very small) that reflects the FIG. 1: Scheme for operationally defining gravitational decoherence. An initial spatial superposition of N nucleon-oscillator coherent states gives rises to a spatial interference pattern in the particle detection probability whenever the particle wavefunctions pass through each other at x = 0. The particle detector is indicated centered at some given x-location. A suppression of the x-dependent interference pattern in the particle detection probability is interpreted as gravitational decoherence.
fact that a real detector is not pointlike, but rather occupies some nonzero volume in space. The visibility can then be obtained in terms of the following expectation value: where the density matrix ρ(t) describes the N nucleons initially in a coherent superposition state and interacting with a thermal graviton bath environment. The expectation value (6) gives a measure of the spatial particle number density (smeared over the small volume V ) and is to be viewed as the field-theoretic counterpart to the configuration space probability density V −1 V dr r|ρ ho (t)|r for a single, non-relativistic quantum harmonic oscillator described by the evolving density matrix ρ ho (t).
The just described set-up shares features of atom and molecular wave interferometry experiments [19][20][21], but utilizing optical traps [22][23][24][25]. The latter enables the two wavefunction components making up the superposition to interfere multiple times as they oscillate through each other, rather than just once as in most matter wave interferometry set-ups. Furthermore, no additional coupled degrees of freedom such as spins manipulated by external magnetic fields in a Stern-Gerlach-type apparatus [24] are required in the system-environment action order to implement the interferometer; a full, relativistic quantum field theoretic analysis can be applied with just the addition of a harmonic confining potential for the scalar field. We must emphasize however that our set-up should not be viewed necessarily as a possible way to feasibly measure gravitational decoherence, but rather as an in-principle operational procedure to quantify the decoherence via the above-defined visibility measure.
With the above-describe challenges in mind, in the present paper we shall consider as a first step, two toy system-environment models that are in turn closely related through dimensional reduction to the above scalar field-gravity system and to scalar field quantum electrodynamics (QED) [26]. The Lagrangian for scalar QED is where φ is a complex-valued scalar field, D µ = ∂ µ − ieA µ is the covariant derivative, and F µν = ∂ µ A ν − ∂ ν A µ is the electromagnetic field strength tensor. We have also included a three-dimensional harmonic confining potential [c.f. Eq. (5)] in order to facilitate operational probes of (de)coherence for initial scalar field spatial quantum superposition states as discussed above.
The models presented in Sec. II below are "toys" in the sense that there is no spatial coordinate-just a time coordinate-and hence are formally zero-dimensional (0d) field models. Our motivation here is to utilize the toy models in order to validate the above-described operational interferometric approach to decoherence as well as certain standard approximation methods, thus giving confidence in eventually applying a similar approach to quantifying actual gravitationally induced decoherence; the 0d model was in fact utilized in Ref. 4 to lend support for an initial gravitational decoherence derivation.
As zero-dimensional field systems, the toy models lack any gauge symmetry, however. It is for this reason that full scalar field QED is also useful for investigating decoherence and verifying that the considered decoherence measures are gauge invariant. In particular, what constitutes a gauge invariant observable is conceptually clearer in scalar QED than in gravity and thus the former also serves as a useful pedagogical stepping stone towards addressing gravitational decoherence.
In Sec. II, we introduce the 0d toy model Lagrangians. Section III analyzes the quantum dynamics of the scalarweak gravity toy model, by utilizing an exact solution to the full system-environment Schrödinger equation assuming an initial system-environment product state, with the oscillator system state expressed in a number state basis and environment in a thermal state. These solutions are then utilized to determine the decoherence dynamics of initial superpositions of system oscillator coherent states through an operational interference fringe visibility analysis that is the single particle counterpart to that described above. Section IV analyzes both the classical and quantum dynamics of the scalar QED toy model. In particular, both classical and quantum Langevin equations as well as a quantum master equation are derived for the system oscillator interacting with its oscillator bath. By making various approximations, the 0d model is shown to map onto that of a simpler oscillator sys-tem with 'two-photon' damping. The master equation is numerically solved to determine the decoherence dynamics of initial superpositions of system oscillator coherent states, again utilizing the operational interference fringe visibility approach. Section V gives some concluding remarks.

II. 0D TOY MODELS
We consider in turn two distinct oscillator systemenvironment models described by the following Lagrangians: and Both model Lagrangians describe an oscillator system with mass M and bare frequency Ω that is coupled to a bath of independent oscillators with assumed identical masses m and frequencies ω i . The two models differ in their system-bath couplings; in particular, the system oscillator couples via its energy to the bath oscillator coordinates in Lagrangian L grav , a 0d analogue of the T µν h µν coupling term in Eq. (4). On the other hand, the interaction term in Lagrangian L qed is obtained via a 0d analogue of the gauge principle of minimal coupling: ∂ µ → ∂ µ −ieA µ . Expanding out the kinetic energy part of Lagrangian (9) gives both cubic and quartic interaction terms, which are respectively linear and quadratic in the bath coordinates [c.f. Eq. (31)]; the full, scalar QED Lagrangian (7) possesses analogous nonlinear terms. Note that the coupling strength parameters λ i in Eqs. (8) and (9) have different dimensions. While Lagrangian (8) is obtained by truncating the spatial dimension from three down to zero for the quantum field theory description of the graviton-matter interaction, and we hence adopt Lagrangian (8) as our 0d gravitional decoherence model in this work, we note in passing that the Lagrangian also yields the standard Hamiltonian of an optomechanical system under the conditions of weak system-bath coupling, where a single optical mode furnishes the system oscillator degree of freedom, while the bath comprises a very large number of mechanical degrees of freedom. This is in contrast to usually-considered optomechanical systems [27], where only one or a few mechanical degrees of freedom are considered. With our focus being on the in-principle operational approach to gravitational and scalar QED decoherence in this work, we shall neglect other possible decoherence mechanisms for both 0D models. The optomechanical interpretation and experimental realizations of Eq. (8) are explored in a separate companion work [28]. As we shall see later in Sec. III, when placed in an initial superposition of coherent states, such a system mode undergoes dephasing without damping-the latter behavior a consequence of the fact that the interaction Hamiltonian commutes with the system Hamiltonian. The resulting, effective system dynamics coincides with that of the so-called 'phase damped' oscillator [17].
In Sec. IV, we show that Lagrangian (9) describes approximately an oscillator system subject to two-photon damping [29]. This is in contrast to the usual quantum Brownian oscillator model with single photon damping and results in qualitatively different decoherence dynamics from the latter for initial superpositions of coherent states.

III. SCALAR GRAVITY MODEL
A. Solving the model Starting with the Lagrangian L grav (8), the system and bath momentum coordinates are where we omit the 'grav' subscript from now on. The model Hamiltonian is For weak system-environment (bath) coupling as in the case of weak field gravity, we can Taylor expand the kinetic energy coupled bath term to obtain approximately Quantizing and expressing the Hamiltonian (13) in terms of the oscillator system and bath creation and annihilation operators, which are defined through the respec- , the scalar gravity model Hamiltonian is where the system-bath coupling is redefined as λ i = 2mωi λ i . We recognize in Eq. (14) the familiar form of the standard optomechanical Hamiltonian, but with a bath of mechanical oscillator modes (labelled by index i) in contrast to the usually considered situation of just one mechanical mode [27].
Solving for the quantum evolution, we will make the common assumption that the system and bath are in an initial product state ρ s ⊗ ρ bath . While the latter assumption facilitates solving for the quantum dynamics, it is not always justified experimentally, since it necessarily requires that the system quantum state can be sufficiently isolated and prepared quickly enough compared to the interaction time scale with the bath degrees of freedom. While such an approximation may be justified for an electromagnetic environment under certain conditions, a mass-energy system can never be isolated from its gravitational environment. Nevertheless, as we shall see, the ability to solve exactly for the scalar gravity model quantum dynamics will give insights into the consequences of assuming a product state.
It is convenient to work in a basis of system number states and bath coherent states |n, {α i } ; the time evolution for such a state can be written as: Following the analysis of Ref. [30], Eq. (15) can be evaluated as: Supposing the bath to initially be in a thermal state, we can express its initial density matrix in the coherent basis as follows [31]: where β −1 = k B T , with k B Boltzmann's constant and T the bath temperature. Decomposing the initial systemenvironment state in the number state basis: we have for the time evolution of the number state outer products after tracing out the bath degrees of freedom: As we shall show later below in Sec. III B, the time evolution of an arbitrary initial reduced oscillator system state can be obtained by decomposing in terms of the number state outer product solutions (19).
In order to carry out the sum over bath degrees of freedom in Eq. (19), we will assume an "Ohmic" bath spectral density with exponential cut-off frequency ω c : We note that the bath spectral density is in general determined by the microscopic coupling strength λ i and the spectral structure of the bath modes; here we consider the linear in ω Ohmic dependence since a thermal graviton bath is Ohmic [4]; a cut-off is necessary in order to tame infinities that would otherwise arise for ω → ∞ in the intermediate stages of the analysis; the exponential cut-off form is primarily motivated by calculational convenience, enabling exact solutions to be obtained for the oscillator system reduced dynamics in the number state basis. In the gravitational decoherence settings, the cut-off energy scale is related to the nature of the 'Planckian' physics, which is not known yet. From the perspective of the effective field theory approach to quantum gravity, this unknown high energy physics is assumed to decouple from the low energy dynamics, appearing only in renormalized parameters that are determined phenomenologically [12]; as we shall see below, the cut-off dependence can be absorbed in part through a frequency and non-linear Kerrtype self-interaction renormalization. However, the cutoff does affect the initial decoherence dynamics, which is tied to the artificial initial product state assumption.
Using Eq. (20) to replace the sum in Eq. (19) with an integral over the continuous variable ω, we obtain Note that for an optomechanical system realization, the environment would have a finite extent resulting in a nonzero, lower frequency cut-off ω 1 : 0 < ω 1 ω c ; here we set the lower frequency cut-off ω 1 = 0, reflecting the actual gravitational environment with effectively infinite spatial extent. In the case of a non-zero lower frequency cut-off, the system will not fully decohere in the long time limit [28], while for the ω 1 = 0 gravitational setting, the dephasing dynamics displays exponential decay with time as we shall see in the following. The integral in the above expression can be evaluated analytically to give Taking the limit β ω c → ∞ (i.e, upper cut-off frequency large compared to the bath temperature), we have With approximation (23), Eq. (21) becomes We now discuss the various terms appearing in Eq. (24).
First, note that the outer product is timeindependent for n = n , a consequence of the fact that the system oscillator Hamiltonian commutes with the system-bath interaction Hamiltonian. The first, pure imaginary term −iΩ(n − n )t in the argument of the exponential is just the free oscillator system evolution. The second pure imaginary term is upper cut-off dependent and comprises both a linear term in system number, which renormalizes the system oscillator frequency Ω, and a quadratic term in system number that is in fact of the same form as the free evolution of a Kerr nonlinear oscillator expressed in the number state basis: Thus, we should properly include a Kerr-type nonlinearity in our starting Hamiltonian (14), with the environmentally induced term i CΩ 2 π (n − n )(n + n )ω c t renormalizing the nonlinear interaction strength Λ kerr . The latter term may be thought of as somewhat analogous to the Newtonian gravitational self-interaction arising from the interaction of a matter system with its gravitating environment. Since we are primarily concerned with decoherence in the present work, we will neglect the quadratic in number term, supposing that it renormalizes an existing Kerr nonlinearity with resulting negligible renormalized coupling strength Λ kerr . For t ω −1 c , the tan −1 (ω c t) term in (25) tends to π/2; this term can be absorbed through a shift in the time coordinate: t →t = t − π/(2ω c ).
Taking into account the system frequency and Kerr nonlinearity renormalizations as just described, Eq. (24) simplifies to where we have redefined the coupling constant as C → C = CΩ 2 and dropped the tilde. The real term on the second line of the argument of the exponential results in decoherence, i.e., exponential decay of the outer product for n = n . In the high temperature (equivalently long time) limit corresponding to t β ω −1 c , Eq. (27) can be approximated asymptotically as From Eq. (28), it is clear that the outer product terms for n = n decay exponentially with rate given by (n − n ) 2 Ck B T / . Note however, that for early, 'Planckian' (by analogy with gravity) times t < ∼ ω −1 c , the rate of decoherence is governed by the upper cut-off ω c , resulting in the logarithm term appearing in Eq. (28); depending on the magnitude of the ratio ω c /k B T 1, there may already be a significant 'burst' of decoherence during the 'Planckian' regime before the later, high temperature exponential decoherence regime. The fact that the decoherence rate depends on the upper cut-off frequency ω c is a consequence of assuming an initial system-environment product state [32]. The latter assumption is tantamount to supposing that the system initial state can be prepared on time scales shorter than ω −1 c (or equivalently, the system-environment interaction is switched on over a time scale shorter than ω −1 c ). While this may be possible for low energy, solid state system environments (i.e., phonons), for an actual gravitational environment with corresponding characteristic Planck time scale, the system state cannot be similarly isolated from the gravitational environment; an analysis which accounts for the system remaining correlated with the environment while its state is being prepared on timescales that are long compared with ω −1 c , is expected to result in a subsequent decoherence rate that does not depend on the upper cutoff frequency of the environment.

B. Decoherence
In the following, we will use Eq. (27) to determine the decoherence dynamics of the oscillator system for initial coherent state superpositions of the form where |α 1 and |α 2 denote coherent states and N is the normalization constant. We consider coherent states since they describe most closely a cooled down, macroscopic oscillator center of mass system. We emphasise that we do not rely on any of the approximations that are often invoked in the study of open quantum system dynamics (beyond assuming an initial product state). In particular, the following analysis is valid for both short/long time scales and high/low temperatures. Note from the form of (the exact) Eq. (27), that the system will evolve into a classical mixture of number states with probability coefficients that are identical to the coefficients of the initial system state. In contrast to other types of system-bath interaction where the final steady state of the system is usually temperature dependent, for the present model the temperature only determines how fast the system decoheres-not its long time limit steady state. From Eq. (27), we also see that the decoherence rate is proportional to (n − n ) 2 for a superposition of two number states |n and |n . Thus, for a superposition of coherent states, we expect that the larger the average energy difference, the more rapid the decoherence. This trend is apparent in the oscillator system Wigner function [33] snapshots shown in Fig. 2. For the initial, example superposition state with α 1 = 3, α 2 = −7, the negative Wigner function regions disappear in the long time limit (signifying loss of quantum coherence). On the other hand, for the initial example superposition states with nearby coherent state parameter magnitudes: α 1 = −α 2 = 3, α 1 = 3 and α 2 = −5, negative Wigner function regions remain in the long time limit (signifying remaining quantum coherence), as is seen more clearly for the zoomed-in Fig. 3. Such trends are consistent with decoherence only resulting for initial spatial superpositions where the states making up the superposition have sufficiently distinct average energies; initial spatial superpositions with the same (or nearby) average energies for the states making up the superposition do not completely decohere. Note also from the Wigner function snapshots in Fig. 2 and Fig. 3 that the initial coherent superpositions phase-diffuse first into crescent-like regions and then eventually into rings. This is consistent with the fact that, as mentioned above, the final state is always a mixture of number states. The above findings are in accord with first investiga-tions on the gravitational decoherence of massive scalar quantum field initial superposition states [4,5], where it is found that superpositions comprising distinct energy states decohere.
Following from the discussion in Sec. I, an operational way (i.e., in principle measurement procedure) to quantify the coherence is through the system oscillator position detection probability density P (x, t) = x|ρ(t)|x [c.f. the full field-theoretic counterpart Eq. (6)] when the two (initially coherent) wavefunctions making up the superposition pass through each other at x = 0; these time instants are τ n = Ωt n = π(n + 1/2), n = 0, 1, 2, . . . for the initial coherent state superposition examples considered above, as can be seen for the early time snapshots in Fig 2. The presence of coherence is manifested in P (x, t) having an oscillatory dependence about x = 0. The latter operational approach corresponds to a two-slit inteference measurement, where the harmonic potential plays the role of the slits by (periodically) bringing the wavefunction components in the initial superposition together. Figure 4 shows the position probability distribution function in the long time limit, steady state for the various example initial coherent state superpositions; we can see that the probability density indicates interference fringes in the vicinity of x = 0 consistent with the presence of negative-valued Wigner function regions shown in Fig. 2; the snapshots can be interpreted as the marginal probability distributions obtained by integrating over the momentum coordinate Wigner function distributions. In particular, the interference remains for α 1 −α 2 , where the average energies of each coherent state making up the initial superposition are not too dissimilar. Note that the other, larger scale scale probability variations in Fig. 4 are due to the final, steady state being a mixture of different number states, as mentioned earlier above.
We adopt the commonly used 'visibility' as a measure of the size of the interference fringes, defined as where P max is the central maximum of the probability density P (x) at x = 0 , and P min is the first local minimum of the probability to the right (or left) of the central maximum. The decrease in visibility over time starting from the initial superposition state (29), provides an operational, quantitative measure of decoherence; Fig. 5 gives the visibility as a function of time for various example, initial coherent state, bath temperature, and systembath coupling parameters. As to be expected, the visibility decreases more rapidly the higher the temperature and the stronger the coupling. Also, the more dissimilar in magnitude α 2 < 0 is from α 1 > 0 (and hence the larger the average energy difference) in the initial coherent state superposition, the more rapid is the decrease in visibility.

IV. SCALAR QED MODEL
In contrast to the scalar-gravity 0d toy model, the scalar QED toy model does not admit an exact, analytical solution for its quantum dynamics. We will therefore utilize various approximation methods towards solving for its quantum dynamics. In particular, we consider both quantum Langevin and master equation approaches, and approximations within these approaches that take advantage of the assumed weak system-bath interaction to show that the model maps onto the oscillator system with two-photon damping.

A. Classical Langevin equation
We start with the Lagrangian L qed (9) and will first derive the Langevin equation that describes the classical oscillator system dynamics interacting with the oscillator bath following the approach of Ref. [34]. Expanding out the kinetic energy term of Eq. (9), we have where we omit the 'qed' subscript from now on. The system and bath momentum coordinates are and the model Hamiltonian is Hamiltonian (34) is to be compared with the gravity toy model Hamiltonian (13), which differs solely in the form of the system coordinate part of the interaction term; both models have in common a quadratic system coordinate coupling, to be contrasted with the usually studied oscillator system-oscillator bath model with interaction term that is linear in the coupled system and bath coordinates.
Hamilton's equations for the system and bath coordi- Formally integrating the equations of motion (35), (36) for the bath coordinates and expressing in terms of the system coordinates: where we have performed an integration by parts that allows to identify system renormalization and damping terms as we shall see below. Substituting the solution (39) for q i (t) into the equations of motion (37), (38) for the system coordinates leads to the following non-linear Langevin equations: where the renormalized system Hamiltonian is the bath memory kernel is and the bath random force function is In particular, the first term on the right hand side of the equals sign in the Langevin equations (40),(41) describes the Hamiltonian evolution, the second term describes nonlinear damping, and the third term the random force. Note that the bath induces a quartic anharmonic potential in the system Hamiltonian (42). Such a term is analogous to a Coulomb self-interaction potential in the scalar QED field system. After making the rotating wave approximation (RWA), the interaction term reduces to a Kerr-type nonlinearity [c.f. Eq. (26)]. Together with 'two-photon' damping [see Eq. (71) below], the resulting open system quantum dynamics can generate quantum states with negative-valued Wigner function regions in the long-time limit, starting from initial Gaussian states [35]. In the following, with decoherence dynamics our main subject of interest, we will neglect this induced potential energy term, supposing that it renormalizes an existing anharmonic potential with resulting negligible renormalized coupling strength.
Assuming a thermal equilibrium canonical ensemble distribution for the initial bath coordinates q i (0), p i (0) in Eq. (44), it can be shown that the fluctuation dissipation relation (FDR) between the memory kernel and the random force follows: where k B is Boltzmann's constant and T is the bath temperature. We shall assume that the bath responds rapidly on the time-scale of the system oscillator dynamics, so that memory kernel is approximated as K(t−τ ) = λ 2 k 0 δ(t − τ ), where k 0 is a constant. The Langevin equations (40), (41) then becomė with the FDR (45) taking the form The above delta function-approximated memory kernel can be obtained from a bath spectral density n(ω) with upper cut-off frequency ω c in the limit ω c → ∞. In particular, for a dense bath spectrum, we can approximate the sum over bath degrees of freedom with a bath spectral frequency integral: Assuming a Lorentzian spectral density the memory kernel (43) then becomes Taking the infinite limit ω c → +∞, we obtain the above delta function-approximated memory kernel: Note that we could equally well have assumed a spectral density with exponential cut-off function instead, as for the gravity toy model [c.f. Eq. (20)]; while the calculations are somewhat more straightforward for the Lorentzian spectral density, we do not expect any qualitative differences in the resulting system quantum dynamics. The motivation to use the Lorentzian spectral density here is purely calculational convenience. The classical, non-linear Langevin equations (46), (47) can be numerically solved as stochastic differential equations as we show in the following sections when comparing with the corresponding quantum dynamics.

B. Quantum Langevin equation
The quantum description is obtained through the correspondence principle where x, p and p i , q i become operators satisfying the canonical commutation relations: with all other commutators vanishing. From Eq. (34), the quantum Hamiltonian operator is where the interaction term on the second line is symmetrized in x and p in order that H is Hermitian. Formally integrating Heisenberg's equations of motion for the bath operators, we obtain the following quantum Langevin equations for the system position and momentum operators: where the force noise operator is given by Eq. (44) with the system/bath coordinates and momenta replaced by their corresponding operators, and K(t − τ ) is the memory kernel given by Eq. (43). It is convenient to express the quantum Langevin equations in terms of the system creation and annihilation operators which are defined through the usual relations x = 2M Ω (a + a † ), Under conditions of weak system-environment coupling, Eq. (57) can be simplified by applying the RWA as we now show. Making the substitution a(t) = A(t)e −iΩt in Eq. (57), we obtaiṅ Dropping fast rotating terms, neglecting time derivatives of A(τ ) (since A evolves at much slower rates than Ω), and setting A(τ ) = A(t) (Markov approximation), Eq. (58) becomes approximatelẏ where we have dropped fast oscillating terms and where γ = Ωλ 2 k 0 /2. For ω c Ω and neglecting the anharmonic interaction term, Eq. (60) simplifies tȯ Defining the noise operator as and utilizing Eqs. (44), (50) and the RWA, the usual noise operator (anti)commutation rules follow: where the Bose-Einstein thermal average occupation number is evaluated at twice the system oscillator frequency: n(2Ω) = e 2 Ω/k B T − 1 −1 . Finally, transforming back to the non-rotating frame, A(t) = a(t)e iΩt , we obtain our desired, RWA quantum Langevin equation: From Eq. (64), we see that the parameter γ has the dimensions of inverse time and characterizes the strength of a nonlinear damping term, while the third term is the nonlinear force noise operator. Equation (64) can be solved numerically as a quantum stochastic differential equation or approximately by first deriving the equations for the various moments in a, b, and their Hermitian conjugates and truncating at some order.

C. Quantum master equation
An alternative way to express the quantum dynamics is via the quantum master equation, where the time evolution is given by the oscillator system reduced density matrix. To second order in the interaction potential and assuming that the bath responds much more rapidly than the system oscillation timescale (Born-Markov approximation), the master equation for system density matrix ρ in the interaction picture is approximately [15][16][17], where ρ B is the initial thermal state of the bath, Tr B denotes the trace over the bath state and V (t) is the system-bath interaction Hamiltonian expressed in the interaction picture: where we have introduced creation/annihilation operators for the bath: In order to simplify the next steps, we introduce the following shorthand notation: Expanding out Eq. (65) and substituting in Eqs. (67) and (68), we obtain: where where n = n(2Ω) = e 2 Ω/k B T − 1 −1 . In Eq. (71), we recognize an oscillator subject to 'two-photon' damping. As a consistency check, we can obtain an equation for the expectation value of a starting either from the quantum Langevin equation (64) with a = Tr (a(t)ρ(0)) or from the master equation (71) with a = Tr (a(0)ρ(t)); both approaches coincide to give ȧ = −iΩ a − γ a † a 2 + 2γn(2Ω) a . (72)

D. Validity of the RWA and quantum vs classical dynamics
Starting with the 0d analogue scalar QED model Lagrangian (9), in the previous sections we derived a Markov approximated classical Langevin equation (46), (47), a Markov-RWA quantum Langevin equation (64), and a corresponding Markov-RWA quantum master equation (71). In the following, we will test the validity of the RWA at the classical level, as well as compare the classical versus RWA quantum dynamics for the averaged quantities a and a † a .
It is convenient to express the classical Langevin equations (46), (47) in terms of the complex coordinates (a, a * ) corresponding to the quantum annihilation/creation operators: In order to solve the non-RWA (73) and RWA (74) Langevin equations, we treat them as classical stochastic differential equations: where W is the standard Wiener process, i.e., a continuous-time random walk [36]. Figures 6 and 7 give numerical solutions to these classical stochastic equations as well as to the quantum master equation (71) (the latter solved using QuTiP [37]) for a range of damping parameters γ and bath temperatures T . These parameters are respectively expressed in terms of the dimensionless Q = Ω/γ factor and thermal average bath occupation number n. The quantum system is initially in a coherent state |α for which a|α = α|α , while the corresponding classical system is given an initial amplitude a(0) = α, in order to allow a direct comparison between the quantum and classical dynamics. From Fig. 6, it can be seen that increasing n and Q −1 both lead to faster decay of the amplitude, signalling the non-linear nature of the damping and noise terms in the system Langevin and quantum master equations. It can also be seen that the difference between non-RWA, RWA and classical vs quantum is barely visible with the chosen parameters. However, such differences clearly show up in Fig. 7 where we consider the time evolution of the average system number a † a . In particular, throwing away fast rotating terms due to the RWA results in smoothing of the oscillating behaviour of the non-RWA time evolution of a † a . Furthermore, the quantum simulation of a † a decays faster than the corresponding classical approximation.

E. Decoherence
In the following, we consider the evolution of system oscillator initial coherent state superpositions of the form where N is a normalization constant. Figure 8 displays the evolving state through its Wigner function representation [33] for a selection of α, n, and Q parameter values-obtained by numerically solving the master equation (71). Quantum coherence manifested in the presence of negative-valued Wigner function regions can survive longer than the amplitude damping time. This is to be contrasted with the commonly-investigated quantum Brownian oscillator model with single photon damping, described by the following master equation: For the latter master equation, decoherence proceeds more rapidly than amplitude damping. Note that the initial, even superposition state (77) is an eigenstate of the operator a 2 since a 2 |±α = α 2 |±α , so that the twophoton loss term in the master equation (71) preserves coherence [29]. In contrast, the even superposition state (77) flips to the odd superposition state N (|α − | − α ) under the action of a single annihilation operator a, hence the single photon loss term in the master equation (78) does not preserve coherence. Figure 9 gives snapshots of the system oscillator position probability density P (x, t) = x|ρ(t)|x when the two initial coherent state wavefunctions making up the superposition pass through each other at x = 0 (at time instants τ k = Ωt k = π(k + 1/2), k = 0, 1, 2, . . . ). These snapshots can be interpreted as the marginal probability distributions obtained by integrating over the momentum coordinate of Wigner function distributions that are similar to those shown in Fig. 8 (but for different parameter values). The presence of quantum coherence is manifested in P (x, t) having an oscillatory dependence about x = 0. In contrast to the gravity toy model (c.f., Fig. 4), the interference fringes survive longer than the initial coherent state peaks; even after 190 cycles a small amount of interference is still present, while the initial coherent states have decayed away.
Proceeding as in Sec. III B for the scalar gravity model, We can operationally quantify the decoherence of an initial superposition of coherent states by using the fringe visibility measure ν (30) for the position detection probability density. Figure 10 shows the time dependence of the visibility ν for a range of parameter choices. The rate at which the visibility is reduced increases with larger damping parameter and bath temperature as for the single photon damping case with master equation (78), but contrary to single photon damping the visibility reduction rate decreases with larger initial amplitude.

V. CONCLUDING REMARKS
In the present work, we have explored two 0d systembath models that share common features with a scalar field system weakly coupled to gravity, and also with scalar QED. The considered model systems comprise a single harmonic oscillator, with the gravitational and electromagnetic fields replaced by a bath of harmonic oscillators, in each case coupled to the oscillator system via non-quadratic interaction terms that resemble the respective scalar-weak field gravity and scalar QED interactions. We utilized these models as a test bed for an operational interference fringe visibility measure of decoherence, as well as for various standard open quantum systems approximation methods.
In particular, we have gained several insights working with the two models that may be of use for analyzing gravitational decoherence: (1) A relatively straightforward, interferometric operational approach for verifying decoherence dynamics can be analyzed that does not involve just extracting the off-diagonal terms of the system reduced density matrix (which is not a gauge invariant quantity in the full theory). (2) While the full scalar matter-gravitational system likely cannot be solved exactly as is the case for the corresponding model system, verified standard open quantum systems approximation methods (e.g., deriving a RWA quantum Langevin equation) may be applicable to the full system; the more involved closed time path integral approaches that are commonly applied to such dynamical quantum field system problems [14] can be guided by the simpler approximation method approaches that are common to nonrelativistic open quantum systems analyses. (3) 'Planckian', cut-off dependent terms can affect the initial decoherence dynamics. However, by being careful with the choice of initial system-environment state taking into account finite state preparation times, such cut-off dependence may be avoided.
The logical next step will be to apply the considered approximation methods to the scalar QED system, and verify that the interferometric observable quantities for probing decoherence are gauge invariant and accessible to analysis at low energy (i.e., 'table top' experiment) scales. We can then apply the lessons learned from the 0d models as well as the full scalar QED model to the scalar matter-weak gravity system. While the latter quantum field system is of course more challenging to analyze, the insights gained from the present work might nevertheless serve as a useful guide in developing an operational understanding of gravitational decoherence, just as the 0d model (8) proved valuable for the initial field theoretic investigation in Ref. [4].