Squeezing the Axion

We apply the squeezed state formalism to scalar field dark matter (e.g. axion) perturbations generated during inflation. As for the inflationary perturbations, the scalar field state becomes highly squeezed as modes exit the horizon. For as long as $H>m_\phi$ (with $H$ the Hubble rate and $m_\phi$ the scalar mass) the scalar field field does not interact during reheating, and we follow its evolution exactly as modes re-enter the horizon. We find that the quantum state remains squeezed after horizon re-entry during the hot big bang. This demonstrates a fact well-known in the theory of inflation: cosmological observables for scalar dark matter are accurately modelled by a classical stochastic field with a fixed phase. Our calculation covers all modes smaller than the present-day cosmic de Broglie wavelength. Larger scale modes mix gravitationally with the environment when $H


INTRODUCTION
The history of the Universe in the standard cosmological model is believed to proceed from an initial period of inflation (Guth 1981;Linde 1982;Albrecht and Steinhardt 1982), during which the Universe is in a quasi-de Sitter state with weakly broken scale invariance. During this epoch, the Friedmann-Robertson-Walker (FRW) scale factor, a, grows exponentially fast with respect to cosmic time, solving the horizon and flatness problems, as required by observations of the cosmic microwave background (CMB) anisotropies (Aghanim et al. 2020). Inflation also provides initial conditions for the small perturbations in the CMB, which arise from quantum vacuum fluctuations of the inflaton field, ϕ. However, the CMB sky is classical, and it is extremely challenging to look for cosmic signatures of our quantum past.
The explanation for this classicality is found within the inflation model itself. Quantum fluctuations are "squeezed" as they exit the horizon (Albrecht et al. 1994;Polarski and Starobinsky 1996). In the squeezed state, the mode functions can be treated as entirely real, with vanishing commutator. In fact, at the level of Gaussianity, all statistics can be matched exactly by classical stochastic fields, regardless of the level of squeezing, as long as the occupation number is large, i.e. nλ 3 1 where n is the number density and λ is the wavelength (see e.g. Sakurai 1967;Martin and Vennin 2016). Inflation creates the large particle number from the Bunch-Davies initial state ("pulling particles out of the vacuum"), and the act of squeezing is the dominance of the growing mode of perturbations over the decaying mode, which has the effect of fixing the phase of the classical stochastic model.
We note that the squeezing picture is somewhat simplistic: in the following we follow no interactions, and the Gaussian states are exactly equivalent to a classical j.kuss@stud.uni-goettingen.de david.j.marsh@kcl.ac.uk model. However the intuition imparted from squeezing, that off diagonal terms in the density matrix become vanishingly small and can be neglected, is essentially born out in the open effective field theory/Linblad equation description of Burgess et al. (2015).
Once all modes of interest have left the horizon during inflation, a further classicalization occurs. The inflaton field decays and "reheats" the other degrees of freedom, producing a thermal bath (Kofman et al. 1997). All particles, photons, baryons, cold dark matter (CDM), neutrinos, are then created out of this thermal bath. For the most part, the memory of the initial quantum state is carried in the classical stochastic correlations of density perturbations in these components (although see e.g. Maldacena 2016;Lim 2015;Martin and Vennin 2016, for discussions about accessing the quantum nature of inflation using cosmological observables).
An interesting possible exception to this story of quantum past and classical present is in the case of scalar field dark matter (DM) such as the axion (Peccei and Quinn 1977;Weinberg 1978;Wilczek 1978) and axion-like particles (see Marsh 2016, for a review). In the standard treatment of this model, the majority of the DM is created by non-thermal processes (Abbott and Sikivie 1983;Dine and Fischler 1983;Preskill et al. 1983). Crucially, the DM production requires no direct interactions with the thermal bath, and proceeds only from motion of the scalar field from a displaced initial state to the vacuum. It is described, at least initially, by a pure state.
The evolution of the classical mean field is generally thought to capture all of the relevant physics. The classical field exhibits a de Broglie wavelength, which manifests as a Jeans scale suppressing cosmic structure formation on small scales (Khlopov et al. 1985). The classical field undergoes Bose-Einstein condensation, and forms axion stars, or solitons (Schive et al. 2014;Levkov et al. 2018). In regions of coherent flows, such as in cosmic filaments, the classical field also displays interference patterns (Schive et al. 2014;Mocz et al. 2019). However, it has been an ongoing topic of research in the field of axion DM to ask whether this classical mean field is in fact the complete and correct description, or whether there are truly quantum and beyond mean field effects, i.e. quantum correlations, that must be accounted for (e.g Sikivie and Yang 2009;Guth et al. 2015;Berges and Jaeckel 2015;Sikivie and Todarello 2017;Hertzberg 2016;Dvali and Zell 2018;Lentz et al. 2019Lentz et al. , 2020.
In the following, we attempt to at least partially answer this question in a simple setting. We consider the case where the Peccei-Quinnn symmetry is broken during inflation. In this case, the axion field, φ, is linear for most of cosmic history, vastly simplifying the computation. In this case, the axion field is smoothed during inflation, giving rise to a background homogeneous condensate, φ 0 (t). This homogeneous field provides the axion DM relic density from vacuum realigment. Furthermore, such a homogeneous field is described exactly by mean field theory, with no exchange correlation (Lentz et al. 2019(Lentz et al. , 2020. The correlations of perturbations in the field are split into two modes: the adiabatic mode, and the isocurvature mode. The adiabatic axion fluctuations have asymptotically zero fluctuation in modes with k > 0 (Hlozek et al. 2015), while the k = 0 mode is described exactly by the classical mean field. All DM structure in this mode is inherited from the the source terms of the thermal bath. Therefore, because these modes are created thermally, like those of the standard cosmological model they will have lost memory of their quantum past.
We thus look to the isocurvature mode for the axion quantum state. Axion isocurvature fluctuations arise during inflation from vacuum fluctuations of the axion field (Axenides et al. 1983). Since the axion is so much lighter than the inflaton (typically m ϕ 10 14 GeV while m φ 1 eV), the axion behaves at all times as a spectator field (Gordon et al. 2000) (although see e.g. Takahashi et al. 2018;Graham and Scherlis 2018;Marsh and Yin 2021). The isocurvature modes leave the horizon just like the inflaton modes and, as we show below, become squeezed. The standard treatment of axion isocurvature perturbations (see e.g. Komatsu et al. 2009;Hertzberg et al. 2008;Hlozek et al. 2018) treats them like the inflaton and the adiabatic mode, with a quantum past, but a classical stochastic model of the structure formation phase during and after the hot big bang.
The fundamental difference between the axion and the inflaton (or indeed a more exotic model such as the curvaton  is that, by necessity if it is to be DM, the axion field is stable, and the vacuum realigment relic density is produced without need for any thermal interaction. During the hot big bang phase, modes that were super horizon after inflation re-enter the horizon. The evolution of the homogeneous background field splits the modes into two categories. The background field begins to oscillate when H ≈ m φ , and it is at this time that the scalar field DM is produced. Modes that enter the horizon prior to axion DM production, never interact or mix gravitationally with any external fields (in linear theory). Thus they cannot decohere, and "unsqueezing" of the initial pure state could occur. Large scale modes re-enter the horizon when H < m φ . The oscillation of the background field serves to "swicth on" the interaction with the gravitational potential, and consequent mixing with the thermal bath. This interaction is expected to decohere the large scale modes (Joos and Zeh 1985) (see also Allali and Hertzberg 2020, for decoherence inside galactic halos). The smaller scale modes (which are smaller than the cosmic de Broglie wavelength/Jeans scale), however, might retain memory of their quantum initial state. It is thus in axion isocurvature fluctuations that enter the horizon prior to axion particle production where the initial quantum state could be encoded in the evolution in the form of the phase, which is normally neglected.
In the present work we model this possibility and, in linear perturbation theory, evolve the axion field from the initial state during inflation, through reheating, and horizon re-entry during the hot big bang phase right up until the moment of axion DM production. We are thus able to compute the squeeze parameter of the axion field, including any possible unsqueezing of sub-de Broglie modes during horizon re-entry.
This paper is organised as follows. In Section 2 we introduce the basics of our cosmological model, in Section 3 we present the squeezing model of inflationary fluctuations, and in Section 4 we compute the squeeze parameter of the axion isocurvature mode and its evolution through horizon reentry. We conclude and discuss our results in Section 5.

MODEL
We employ a simple toy model for the background evolution consisting of single field slow-roll inflation with a quadratic potential V (ϕ) = 1 2 m 2 ϕ ϕ 2 . This choice is made purely for simplicity of the analytic background solution. Our conclusions concerning the axion quantum state are not strongly affected by this choice as long as 1) Peccei-Quinn symmetry is broken during inflation and not restored afterwards and 2) the axion remains a spectator field. Inflation is followed by instantaneous reheating leading to a radiation dominated epoch (i.e. the hot big bang). The FRW scale factor, a(t), is then given by 2mϕ . Reheating occurs for ϕ reheat = 1 √ 4πG , the scale factor is normalised to unity here, a reheat = 1. The mass of the inflaton is taken to be m ϕ = 10 13 GeV ≈ 10 −6 M P (Peter and Uzan 2009) (where M P = 1/ √ 8πG), the initial field value is chosen as ϕ i = 4 √ G yielding approximately one hundred e-folds of inflation. For the numerical treatment, it is convenient to use the natural logarithm ln(a) of the scale factor as the time variable. We can then express the Hubble parameter as The evolution of H and the comoving horizon is shown in Fig. 1. In this background model, the axion does not play any role: it is assumed to contribute only negligibly to the total energy density, which hence is dominated by the inflaton or by radiation respectively. The background evolution of the homogeneous axion field is given by: where we have assumed that the initial displacement of the axion field is less than the Peccei-Quinn symmetry breaking scale, φ < f a , such that the axion potential is approximately quadratic. Again, this choice is motivated purely by simplicity and our conclusions concerning the quantum state would not be strongly affected by using the full axion potential. We approximate the axion mass as constant, as for a generic axion-like particle. A temperature dependent mass, as for the QCD axion, would also not change the character of our conclusions, which apply in the limit H m φ . Within these approximations, the axion behaves like a damped oscillator and is essentially frozen at its initial value until the Hubble parameter becomes comparable to the axion mass, which we assume happens deep in the radiation dominated era.
For studying the perturbations we follow Gordon et al. (2000) in separating curvature and isocurvature perturbations. They consider the space spanned by two scalar fields χ 1 and χ 2 . By a rotation in this field space, a perturbation is split into a component along the the trajectory of the background field's evolution and one component perpendicular to it. It turns out that the former, δσ = cos(θ)δχ 1 + sin(θ)δχ 2 (θ being the angle between χ 1 and the background trajectory), is associated with curvature perturbations, R and the latter, δs = − sin(θ)δχ 1 + cos(θ)δχ 2 , with fluctuations in entropy, S = H δṗ p − δρ ρ . Identifying χ 1 with the inflaton and χ 2 with the axion, one finds that as long as the axion is frozen in by the Hubble damping, the angle θ is constantly zero. The inflaton is already the curvature field and the axion the isocurvature one. Further more, with θ = 0, the equation of motion for axion perturbations (in Fourier space) simplifies tö

SQUEEZED STATES
We follow Albrecht et al. (1994) and especially Polarski and Starobinsky (1996) in the application of the notion of squeezed states on primordial cosmic perturbations.
Here we briefly review the formalism based on Polarski and Starobinsky (1996) using the example of a generic massless field χ on an FRW background. It is convenient to use conformal time τ and define the variable y(τ ) ≡ a(τ )χ(τ ). The Lagrangian density for y can be written after integration by parts as where primes denote derivatives with respect to conformal time. With the momentum p(τ ) ≡ ∂˜ L ∂y = y − a a y one can construct the Hamiltonian (in Fourier space) as The time dependence of coordinate and momentum then are Now the field is quantised and one defines the common ladder operators: From the canonical commutation relation forŷ and p, the ladder operators inherit their usual commutator [â k (τ ),â † k ] = δ (3) (k − k ). One should note that the Heisenberg picture is employed here, where the timedependence is associated with the operators and the states remain time-independent.
To find the time evolution of the ladder operators, Eq. (8) is inserted into (7). This yields: An ansatz to solve these coupled differential equations iŝ which constitutes a Bogolubov transformation of the initial time ladder operators. As the commutation relations for the ladder operators shall hold for all times, one finds |α k (τ )| 2 − |β k (τ )| 2 = 1. This allows to parametrise the coefficient functions α k (τ ) and β k (τ ) as α k (τ ) = e −iθ k (τ ) cosh(r k (τ )), where r k , ϕ k and θ k are the squeezing parameters of interest.
To find the time evolution of the squeezing parameters, one first inserts Eq. (10) into Eq. (9). This yields The derivatives of the coefficient functions are Together these identities can be used in appropriate combinations [namely 1 2 (e iθ k α k + e −iθ k α * k ), 1 2 (e iθ k α k − e −iθ k α * k ) and 1 2 (e −i(θ k +2ϕ k ) β k + e i(θ k +2ϕ k ) β * k )] to derive the time dependence of the squeezing parameters. The result is With the mode function f k (τ ) = 1 √ 2k (α k (τ ) + β * k (τ )), one can write the field operator aŝ

SQUEEZING THE AXION
We now apply the formalism of squeezed states to the evolution of axion perturbations during inflation and early radiation domination. To do so, we derive the equations of motions for the squeeze parameter r k , ϕ k and θ k for the axion and solve them numerically.
The axion's equation of motion, Eq. (4), corresponds to the action Again, it is useful to define u = aφ. After integration by parts, the action becomes (17) To bring it to a form similar to Eq. (5) one subtracts the total derivative a a u 2 , yielding The corresponding Hamiltonian then iŝ withp k =û k − a aû k . One findsû =p + a aû and −p = m 2 effû + a ap . Introducing ladder operators as in Eq. (8), their equations of motion are: Now, it will be useful to define and we rewrite the expression a a + i 2 k − can apply a parametrisation as in Eq. (10) and (11) to obtain as the equation of motion for the squeezing parameters of the axion perturbations.
The initial conditions for the squeezing parameter r k , ϕ k and θ k can be inferred by comparing the mode function of the perturbations, to the Bunch-Davies vacuum (Baumann 2011 -Evolution of the squeezing parameters of axion perturbations for a given arbitrarily chosen k mode during inflation and radiation domination The mode in question leaves the horizon at ln(a) ≈ −19.7 and enters again at ln(a) ≈ 18.5. The state squeezes during inflation and continues becoming more and more squeezed after reheating. After reentering the horizon, the squeeze factor r k remains a large constant. The total phase also remains a constant after horizon re-entry. (There are numerical artifacts in the late time evolution of r k , since the rapidly growing ϕ k enters r k via cos. Early time numerical artifacts in the evolution of ϕ k occur due to limited temporal resolution.) In slight contrast to Albrecht et al. (1994), we find that the mode function coincides to the Bunch-Davies vacuum for: Now we numerically integrate equations (22). We assume the axion mass to be m φ = 1 meV, roughly at the upper end of the possible mass range for a QCD axion. The result is presented in Fig. 2 for a particular mode with k = 10 5 GeV. 1 One sees that r k grows (in modulus) when the corresponding mode leaves the horizon, i.e. the state becomes squeezed. Meanwhile φ k and θ k approach constant values. This pattern does not change with reheating. When the mode reenters the Horizon, φ k and θ k grow in modulus, but their sum remains constant. Accordingly, r k stays large, and the state does not become unsqueezed. The growth of the squeeze parameter to large values is consistent with the growth in comoving particle number density while the field remains frozen by Hubble friction.
That the axion field does not becomes unsqueezed upon horizon reentry in the radiation epoch is a result of physical interest. Squeezing initially becomes large when a mode exits the horizon, which happens due to the shrinking of aH. However, we observe that squeezing of superhorizon modes continues when the evolution of aH turns around at reheating. Subsequently, squeezing ceases, but does not reverse, when the mode reenters the horizon. The simple statements "modes squeeze on horizon exit", or "modes squeeze as the comoving horizon shrinks" cannot be naively time reversed. The radiation epoch and horizon reentry are not time reverses of inflation for squeezed states, despite superficial time reversed similarities.
The squeezed axion field can be modelled as a classical stochastic field with fixed phase at late times, and the commutator can be set to zero. This is because for large squeezing (large r k ) one can combine the two angles from Eq. (14) to get lim |r k |→∞ (θ k +ϕ k ) = 0, as both coth(2r k ) and tanh(r k ) approach 1 in this limit. Hence, one has θ k + ϕ k → δ k for a constant phase δ k (Fig. 2 shows this  clearly). For the function f k , it then follows (Polarski and Starobinsky 1996) The constant factor e −iδ k can be removed by a timeindependent phase rotation, hence f k can be made real for all times τ in the limit of large squeezing (r k → ∞).
Lastly we check that our formalism reproduces the expected results of approximate scale invariance for the isocurvature spectrum, and the subsequent classical evolution of the axion field. From the squeezing parameters, we calculated the power spectrum of the axion perturbations by P = k 3 Fig. 3 shows the time evolution of the mode k = 10 5 GeV. One sees that the power spectrum remains constant outside the horizon and decays upon reentry.
Evaluating the spectrum for different k at reheating (a = 1) yields approximate scale-invariance, as Fig. 4 demonstrates. The spectrum can be fit by the expected slow-roll result: where we choose the pivot scale k * = 1 GeV for which H * ≈ 4.9 · 10 13 GeV is the Hubble scale when this mode leaves the horizon (at a = e −31.5 , see Fig. 1), and the isocurvature spectral index is n I = 2 and = −Ḣ H 2 ≈ 0.016 is the slow roll parameter (see e.g. Komatsu et al. 2009). Hence, we recover the expected behaviour for isocurvature initial power spectrum in slow roll inflation. We also show the spectrum evaluated at a later time, a = e 20 , when small scale modes have reentered the horizon, showing the expected suppression and oscillatory features below the horizon scale found in the standard classical evolution (Khlopov et al. 1985;Hlozek et al. 2015;Hlozek et al. 2018).

CONCLUSION
We have computed the squeeze parameter of an axionlike scalar DM field during inflation and early radiation domination. The field picks up isocurvature perturbations from quantum vacuum fluctuations during inflation. These perturbations become squeezed when the wavelength becomes larger than the cosmological horizon, and as such the field and momentum operators can be rotated to a basis where the canonical commutators vanish. The scalar field is assumed to be stable and form the cosmic DM. Therefore, perturbations reenter  the horizon after reheating. In this regime mode functions begin to oscillate, however the overall phase remains constant, as does the squeeze parameter. Thus, the state of scalar field DM inflationary vacuum fluctuations is well described by a squeezed pure state on all scales that re-enter the horizon while the field is effectively massless, m φ H. We reiterate that the modes we have studied are expected to be the "most quantum" in axion cosmology. These modes are generated from vacuum fluctuations, and have not interacted with the thermal bath in any way. These modes re-enter the horizon before the scalar field mass becomes relevant, and correspond to scales below the cosmic de Broglie wavelength (Jeans scale).
Our calculation is valid as long as the background field is frozen by Hubble friction, that is, as long as the scalar field/axion mass is small compared to the Hubble parameter, m φ H. In our model, this happens around ln(a) = 28; the largest modes that reenter the horizon before this are at the order of 10 GeV, which redshifted to today correspond to ≈ 1 · 10 10 Mpc −1 When φ unfreezes and oscillates, interactions between the adiabatic and the isocurvature mode give rise to mixing between the axion perturbations and the thermal bath on all scales, mediated by gravity. In particular, this mixes isocurvature and curvature (Gordon et al. 2000) converting one into the other on large scales and at late times. Even in linear cosmological perturbation theory, this interaction mixes the axion field with the radiation bath (thus generating a non-zero curvature from the initial isocurvature state, see e.g. power series of isocurvature initial states in Bucher et al. 2000;Hlozek et al. 2018), and might be expected to give rise to decoherence of the isocurvature mode even before the non-linear regime (the same mechanism occurs in the curvaton model, Lyth and Wands 2002, except that unlike the curvaton the axion does not decay).
Scalar field DM fluctuations in the orthogonal adiabatic mode are generated solely from gravitational interactions with the radiation bath, and so must be described by a classical stochastic model. Furthermore, super-Jeans modes eventually grow and undergo nonlinear gravitational collapse into galactic DM halos. Allali and Hertzberg (2020) has shown that gravitational scattering of DM with hydrogen inside galaxies indeed leads to decoherence. If the field has non-gravitational interactions these also becomes relevant when m φ H (when the shift symmetry is broken), and will lead to additional decoherence. One possible application of our result is to provide the initial correlation state of the axion field in the preinflation symmetry breaking scenario. Another application of our result could be to provide the correct initial quantum state for the Peccei-Quinn field in the scenario when the U (1) PQ symmetry is broken after inflation, simply by replacing what we have termed the axion with two massless real fields. In this case, the subsequent evolution is highly non-linear, although governed by local interactions in the field potential rather than by gravity (see e.g. Berges and Jaeckel 2015;Gorghetto et al. 2018;Vaquero et al. 2019).
The present work can be seen as a quantum mechanical interpretation of well-known results in linear perturbation theory of scalar field/axion DM. It would be fruitful to further consider the quantum mechanical interpretation of non-linear phenomena, such as the axion bispectrum (e.g. Langlois and van Tent 2012)). Signals of quan-tum mechanics in cosmology are generically thought to arise as poles in higher-point correlation functions (Green and Porto 2020). However, these "cosmological collider" (Arkani-Hamed and Maldacena 2015) phenomena are generally present for particles with m H I , while the case we have studied applies in the opposite limit. While the present work was in preparation, Lu (2021) appeared considering such an "axion isocurvature collider". It would be interesting to develop the open-EFT/Linblad formalism of Ref. Burgess et al. (2015) applied to DM, which extends rigorously the intuition about classicalization derived from the squeezed state formalism.
Note Added : In the final stages of preparation, Kopp et al. (2021) appeared also discussing the squeezing of axion perturbations.