Cosmological abundance of the QCD axion coupled to hidden photons

We study the cosmological evolution of the QCD axion coupled to hidden photons. For a moderately strong coupling, the motion of the axion field leads to an explosive production of hidden photons by tachyonic instability. We use lattice simulations to evaluate the cosmological abundance of the QCD axion. In doing so, we incorporate the backreaction of the produced hidden photons on the axion dynamics, which becomes significant in the non-linear regime. We find that the axion abundance is suppressed by at most ${\cal O}(10^{2})$ for the decay constant $f_a = 10^{16}$ GeV, compared to the case without the coupling. For a sufficiently large coupling, the motion of the QCD axion becomes strongly damped, and as a result, the axion abundance is enhanced. Our results show that the cosmological upper bound on the axion decay constant can be relaxed by a few hundred for a certain range of the coupling to hidden photons.

Introduction -The axion is a pseudo Nambu-Goldstone boson in the Peccei-Quinn (PQ) mechanism [1][2][3][4], which provides a dynamical solution to the strong CP problem in the quantum chromodynamics (QCD). The axion is copiously produced in the early Universe, and it is one of the plausible candidates for dark matter (DM). There are many on-going and planned experiments aiming to detect the axion. See Refs. [5][6][7][8] for recent reviews on the QCD axion and axion-like particles.
The axion mass and coupling strengths are characterized by the decay constant, f a . The so-called axion window between the astrophysical and cosmological bounds is given by 4 × 10 8 GeV f a 10 12 GeV, where the lower bound is due to the neutrino burst duration of SN1987A [9][10][11] while the upper bound is due to the cosmological abundance of the axion produced by the misalignment mechanism [12][13][14]. The latter assumes the initial displacement of order f a from the CP-conserving minimum. The upper bound can be relaxed if one consideres entropy production in the early Universe [14][15][16][17][18], non-standard couplings to hidden monopoles [19][20][21] and a resonant mixing with axion-like particles [22].
where φ is the axion field, F Hµν = ∇ µ A Hν −∇ ν A Hµ is the field strength of the hidden U(1) H gauge field A Hµ and F µν Here and in what follows, we assume that there are no light hidden matter fields charged under U(1) H . The coupling constant β H ≡ g φγ f a is model-dependent, and as pointed out in Ref. [24] it can be much larger than unity if one considers the clockwork/aligned QCD axion model [27,28] (See also [29][30][31][32][33][34][35]).
Hidden photons are known to be explosively produced by tachyonic instability if the axion field evolves in time with a moderately large coupling β H [36]. As a result, the kinetic energy of the axion efficiently turns into hidden photons. In a context of axion inflation, this explosive production of the gauge fields affects the inflaton dynamics as well as the reheating processes [37]. In particular, the backreaction on the axion field is significant for a large coupling, and so, it is necessary to use lattice simulation for correctly describing the evolution of the inflaton and the reheating processes [38].
The coupling of the QCD axion to hidden photons (1) should also modify the evolution of the QCD axion, if the coupling β H is large enough. Recently, it was pointed out in Ref. [39] that such coupling can significantly reduce the energy density of coherent oscillations of the axion. (See a note added at the end of this Letter.) The main purpose of this Letter is to investigate how the coupling to hidden photons affects the cosmological abundance of the QCD axion including nonzero modes. Due to the explosive nature of the particle production, the system soon enters highly non-linear regime, which makes it dif-ficult to describe the entire evolution analytically. We, therefore, use lattice simulations to make a quantitative estimate of the axion abundance by incorporating the backreaction effects. To the best of our knowledge, this is the first lattice simulation for the QCD axion coupled to hidden photons.
We find that the QCD axion abundance is suppressed by at most O(10 2 ) for the decay constant f a = 10 16 GeV for a moderately large coupling β H 200. Then, the cosmological upper bound on f a can be relaxed by at most a factor of a few hundred. As we shall see, the reduction of the axion energy stops in a few Hubble times after the system enters non-linear regime, and both axion and hidden photons gradually asymptote to an equilibrium state. For a sufficiently large coupling, on the other hand, the motion of the QCD axion becomes strongly damped, which enhances the final axion abundance compared to the case without coupling to hidden photons. Therefore, the QCD axion abundance is not a monotonic function of the coupling to hidden photons: it is suppressed for β H 300 − 400, while it is enhanced for β H 300 − 400.
Axion electrodynamics -We start from the following Lagrangian of a system with the axion φ and a hidden photon A Hµ : The magnitude of the axion potential is given by the QCD topological susceptibility χ QCD , which is dependent on the temperature T . We adopt the result from the recent lattice QCD calculation in Ref. [40] and approximate χ QCD (T ) as with χ 0 = (7.6×10 −2 GeV) 4 , T c = 0.16 GeV and b = 8.2.
In the rest of this Letter, we normalize the scale factor a to unity when T = T c . In the temporal gauge, A Hµ = (0, A H ), the Lagrangian (2) leads to the following equations of motion: where the overdots represent derivatives with respect to the conformal time τ and H =ȧ/a is the conformal Hubble parameter. Note that the last equation is the constraint equation since the longitudinal component of the gauge field is not dynamical as in the usual electrodynamics. When the spatially homogeneous φ evolves in time, one can recast Eq. (5) intö in the Fourier space [36], where k ≡ |k| is the wave number, and the subscript ± indicates the helicity. One can see that one of the helicity components of the hidden photon A Hk,+ or A Hk,− satisfying k < β H |φ|/f a becomes tachyonic. Thus the kinetic energy of the axion zero mode is efficiently transferred to hidden photons by the tachyonic instability.
Once a large number of hidden photons are produced by tachyonic instabilities, they source the axion fluctuations through Eq. (4). Due to their large occupation number, such reflux can be so efficient that the axion zero mode is soon swamped by the nonzero modes. We have solved the equations of motion for the axion and the gauge field and confirmed that this is indeed the case. To make a quantitative estimate of the axion relic density, therefore, one must resort to lattice calculations.
Simulation -We have performed a simulation of the axion electrodynamics with a finite-difference method. Specifically, we have extended the Yee's algorithm [41] to accommodate the axion in the staggered grids, and adopted a periodic comoving box with a grid number N grid = 256 3 . The time integration is implemented by the leap-frog method, so that our algorithm is accurate up to the second order in both the time and spatial domains. The staggered grids are advantageous because the discretized constraint Eq. (6) can be satisfied automatically apart from round-off errors, which we have explicitly verified by monitoring the gradient ofȦ H − βH fa φ∇ × A H . Details of our simulation will be presented in a separate paper [42].
Axion abundance -We show in Fig. 1 the timeevolution of the energy densities (left) and number densities (right) of the axion and gauge field obtained in the lattice simulation with f a = 10 16 GeV, β H = 0 (red), 20 (green), 300 (blue) and the initial amplitude of the axion field φ ini = (π/3)f a . One can immediately see that for β H = 20, the axion number density gets suppressed as soon as the energy density of the gauge field becomes comparable to that of the axion. In this phase, the axion coherent oscillation is rapidly dissipated into the gauge field. At almost the same time, the axion nonzero modes are produced by the backreaction from the gauge field (see Fig. 2). After this event, the initial oscillation energy of axion is exhausted dominantly by the gauge field, and also fractionally by the axion nonzero modes.
For β H = 300 on the other hand, the gauge field increases more rapidly soon after the onset of the axion oscillation. In that case, the backreaction becomes significant before the axion completes a single oscillation  and acts as a friction on the axion dynamics. As a result, the axion rolls down the potential more slowly, which effectively delays the onset of the axion oscillation. This implies that when β H is further increased, the final abundance of the axion can be enhanced rather than suppressed, which will be shown explicitly later. The friction dominated regime is often applied to the axion warm inflation [37,38] (See also Refs. [43,44] for a similar effect in the relaxion model). Fig. 2 shows the energy spectra of the axion and gauge field for β H = 20 plotted at several different scale factors. We can see that, at a 5, the gauge field at low wave number is enhanced via the tachyonic instability. Subsequently, the axion nonzero modes are produced from the mode-mode coupling of the gauge field. Notably, the axion energy spectrum is flat in contrast to that of the gauge field (See the energy spectra at a = 5). This means that even the axion nonzero modes with momenta higher than the maximum tachyonic wave number, i.e. k β Hφ /f a , are produced. This is possible because these particles are produced from the mode-mode coupling of the gauge field between enhanced tachyonic mode (low momenta) and the vacuum fluctuations (high momenta). The modemode coupling of the axion nonzero modes and the gauge field sources the gauge field with different momentum as seen in the RHS of Eq. (5). When a 6, both the gauge field and the axion are populated over a wide range of momenta with a similar spectrum. This suggests that the axion nonzero modes and the gauge field are coupled with each other effectively in the nonlinear regime.
We finally examine how the suppression of the axion abundance depends on the coupling strength β H . Fig. 3 shows the ratio of the asymptotic abundance of the axion to that of the axion zero mode in the absence of the coupling. We can see from the figure that for small β H 10, the suppression is virtually absent. This is because the coupling is not strong enough for the tachyonic growth of the gauge field to affect the axion dynamics. With moderately large β H 10, a large amount of the gauge field is produced, and the axion abundance starts to be suppressed. As β H increases, the extent of the suppression of the axion abundance becomes more remarkable. For β H ∼ 200, however, one can see that the axion abundance turns upward and continues to increase for a larger coupling. This regime corresponds to the overdamping of the axion coherent oscillation. We have found that it takes a longer time for the system to settle down in equilibrium for larger β H , and that the final energy spectra become relatively hard. We have confirmed that our results are robust against changing lattice resolution, simulated time, and the initial condition. Therefore, we believe that our finding that the axion abundance is enhanced in the overdamping regime is valid, even though the enhancement factor could be slightly modified by more detailed lattice calculations.
The critical value of β H corresponding to the turnaround between the suppression and the overdamping regime can be analytically estimated as follows. Let us focus on the brief period after the axion starts rolling but before it reaches the potential minimum. The back reaction on the axion dynamics becomes important when the energy density of the hidden gauge field becomes comparable to the axion energy density. If this happens before the axion reaches the potential minimum, the system is considered to enter the overdamping regime. From Eq. (7), one obtains the peak wavenumber of the produced hidden photons by tachyonic instability, k peak ∼ β Hφ /2f a . Near the peak wavenumber, one of the circular polarization modes of the hidden gauge field is exponentially amplified as where ∆φ is the axion field excursion during the amplification. Thus, the initial energy density of the hidden gauge field at the onset of the axion oscillation, ρ H,ini ∼ (k peak /a) 4 is amplified by a factor of exp (β H ∆φ/f a ). Equating the energy density of the axion and hidden gauge field when the axion reaches the potential minimum, i.e., ∆φ = φ ini , we obtain Thus we arrive at up to a logarithmic correction, where θ ini = φ ini /f a is the initial misalignment angle. For our reference values, θ ini = π/3 and f a = 10 16 GeV, we obtain β H ∼ 200. This agrees well with our numerical result, as one can see the turnaround around β H = 200 in Fig. 3. Discussion and conclusions -We have studied the evolution of the QCD axion coupled to hidden photons using lattice calculations, and found that the backreaction of the produced hidden photons on nonzero modes of the axion becomes significant in the non-linear regime. In particular, while the axion density is suppressed for a moderately large coupling, it is enhanced for a sufficiently large coupling. The suppression factor of the axion density can be as large as O(10 2 ) for f a = 10 16 GeV and β H = 200. Although the precise value of the maximum suppression factor may be altered depending on f a because it changes the relative time scale of axion oscillation compared to the Hubble time and hence the spectral peak of the axion nonzero modes, we believe that our main result should hold qualitatively for smaller values of f a as well as for axion-like particles coupled to hidden photons. In particular, one may be able to enhance the QCD axion abundance to explain the DM abundance even for f a near the smaller end of the axion window for a sufficiently large β H . More detailed analysis for varying f a and initial misalignment angle is left for the future work.
If U (1) PQ symmetry gets broken after inflation, cosmic strings and domain walls are produced. They may experience an additional frictional force in the presence of a strong coupling to hidden photons, which modifies the resultant axion abundance produced from these defects (see e.g. Ref. [45]). Also, the anharmonic effect as well as axionic isocurvature perturbation may also be affected. We leave these issues for the future work.
Note added: While preparing this Letter, we found Ref. [39] in which its authors studied the same set up as in Eq. (2). They presumed that the backreaction effect was small for their choice of the parameters, and found that the axion abundance can be significantly suppressed by a factor between 10 4 and 10 13 for β H = 20 -30. Although our lattice simulations did not confirm such a large suppression factor, the effect on the axion density becomes important for a similar value of β H 10. Also, for the limited suppression factor found in our calculations, the initial energy density of the QCD axion must be much smaller than the total energy density at the onset of oscillations. Thus, the produced hidden photons cannot give an appreciable contribution to the effective neutrino species.