Gravitational Waves and Dark Photon Dark Matter from Axion Rotations

An axion rotating in field space can produce dark photons in the early universe via tachyonic instability. This explosive particle production creates a background of stochastic gravitational waves that may be visible at pulsar timing arrays or other gravitational wave detectors. This scenario provides a novel history for dark photon dark matter. The dark photons may be warm at a level detectable in future 21-cm line surveys. For a consistent cosmology, the radial direction of the complex field containing the axion must be thermalized. We explore a concrete thermalization mechanism in detail and also demonstrate how this setup can be responsible for the generation of the observed baryon asymmetry.


Concrete
Spontaneously broken global symmetries are essential ingredients in solutions to a variety of problems in the Standard Model. Examples include a Peccei-Quinn symmetry [1] (a solution to the strong CP problem), lepton symmetry [2] (the origin of neutrino masses), and flavor symmetries [3] (the origin of the pattern of the fermion masses and mixing). Associated with the breaking of such symmetries are light degrees of freedom-Nambu-Goldstone bosons.
The resulting particles go by various names: a QCD axion [4,5] for the Peccei-Quinn symmetry, a Majoron for lepton number [2], and a familon [6,7] in the case of a flavor symmetry. Here we refer to these Nambu-Goldstone bosons generically as "axions." It is commonly assumed that an axion field is initially static and begins oscillation when its mass becomes comparable to the Hubble expansion rate. However, the axion field may instead be initially rotating in field space, with important consequences. Such rotation may be initiated by an explicit breaking of the global symmetry coming from a higher dimensional operator containing the symmetry breaking field. When the radial direction of the symmetry breaking field takes on large values, these operators are of increased importance and can give the axion a kick. These ingredients are generic and are found, e.g., in Affleck-Dine baryogenesis [8]. This setup has been studied in the context of baryogenesis: axiogenesis [9], ALP-genesis [10], and lepto-axiogenesis [11]. Its impact on the axion dark matter abundance has been explored under the name of the kinetic misalignment mechanism [12,13]. In the above examples, axions directly couple to Standard Model particles.
However, if the coupling of the axion to particles in our sector is too weak, it can prove challenging to probe. In this case, one may have to rely on signals from gravitational forces, e.g., superradiance [14] and mini-clusters [15], or perhaps signals from dark sectors. Here, we consider a rotating axion that couples to a dark photon and explore two consequences: a gravitational wave (GW) signal and the production of dark photon dark matter.
When the axion undergoes motion, it will modify the dispersion relation of a dark photon to which it couples. If the coupling is sufficiently large, the dark photons become tachyonic for a range of dark photon wave numbers. In the conventional case where the axion is assumed to oscillate, the production of dark photons via tachyonic instability is ineffective for natural coupling strengths. Nevertheless, if the coupling is much larger, dark photons are explosively produced [16]. The produced dark photons are highly inhomogeneous and source gravitational waves [17,18]. This scenario was recently considered to explain the North American Nanohertz Observatory for Gravitational Waves (NANOGrav) signal [19], but it was found that the axion oscillation leads to an excessive amount of dark matter [20] unless the axion mass decreases after the beginning of the oscillations [21,22]. The produced dark photons may be dark matter [23][24][25].
As we will explore in this work, an axion that rotates in field space can also produce dark photons via tachyonic instability. In contrast to the oscillation case, the sign of the velocity of the axion field does not alternate. A consequence is that dark photon production may occur without any need to enhance the axion-dark photon coupling. The tachyonic instability from a rotating field was invoked as way to produce the hypercharge gauge field in the context of magnetogenesis [26]. Furthermore, the rotational velocity, which determines the strength of the GW signal, is not directly tied to the axion mass. In this case, the axion may be nearly massless around the minimum of the potential, so the NANOGrav potential signal can be explained without overproduction of axion dark matter. Throughout this work, we will assume such a nearly massless axion, although the mechanism may also work for a massive axion such as the QCD axion. This mechanism can produce GW signals in a wide range of frequencies of interest to current and future GW detectors such as NANOGrav [27], the Square Kilometer Array (SKA) [28], the Laser Interferometer Space Antenna (LISA) [29], the DECi-hertz Interferometer Gravitational wave Observatory (DECIGO) [30], and the Big Bang Observatory (BBO) [31].
We study under what assumptions the dark photon can comprise the dark matter and when an alternate candidate is needed. In the latter case, the dark photon comprises dark radiation.
As mentioned above, the axion field's rotation may be induced by the combination of a large field value in the radial direction and an explicit breaking of the symmetry. In this case, radial motion is also expected. To avoid cosmological difficulties such as a moduli problem, radial motion should be dissipated via thermalization. The necessity for thermalization places non-trivial constraints on the theory.
In the concrete thermalization model we consider, couplings between the axion and Standard Model fermions are generated. The global charge asymmetry in the axion rotation is partially translated into the particle-antiparticle asymmetry in the thermal bath, including the baryon asymmetry [9,10]. We show that the observed baryon asymmetry can be simultaneously explained in some of the parameter space.
The structure of the paper is as follows. In Sec. 2, we discuss the dynamics of the axion, including its rotation and how its coupling to a dark photon can lead to tachyonic instability.
In Sec. 3, we discuss how the strength and frequency of produced GWs are related to the underlying axion parameters. In Sec. 4, we turn to a discussion of when and how the dark photon may comprise the dark matter in this model. The warmness of the dark photon gives an important constraint as well as a possible future signal. In Sec. 5, we investigate a specific realization of this scenario, paying particular attention to the thermalization of the radial component of the symmetry breaking field and the baryon asymmetry produced from the axion rotation. Technical details on the cosmology are postponed to Appendix A. Sec. 6 is devoted to the summary of the results and discussion.

Axion rotations
We consider a field theoretical axion, a pseudo-Nambu-Goldstone boson associated with the spontaneous breaking of a global U (1) P symmetry. The axion φ is the angular direction of a complex scalar field P , where S is the radial direction that we refer to as the saxion. Its value at the minimum of the potential determines the decay constant f φ = S , and θ = φ/f φ is the angular direction.
If, as conventionally assumed, the axion field has zero initial velocity, it begins oscillations when its mass is comparable to the Hubble expansion rate, m φ ∼ H. This picture, however, may be too simplistic. In the early universe, the saxion field value may be much larger than f φ . For instance, the saxion may be dynamically driven to a large field value due to a negative Hubble-induced mass during inflation [32]. The field value may be stabilized by a higher dimensional operator. When the saxion takes on such large field values, it is not guaranteed that the potential of P remains nearly U (1) P symmetric. A higher dimensional potential term that explicitly breaks the global symmetry, ∆V ⊃ P n M n−4 + h.c., (2.2) where M is a cut-off scale, is enhanced at large field values. Such terms are in fact expected when the U (1) P symmetry accidentally arises as a result of another exact symmetry [33][34][35][36].
This explicit breaking can initiate an angular motion for P . Cosmic expansion will cause the saxion field value to decrease, with an accompanying improvement in the quality of the global symmetry. Once the global symmetry is approximately conserved, the field P continues to rotate so as to preserve the angular momentum in field space. This is because the angular momentum may be identified with the charge of the global symmetry associated with the rotation, n θ . The rotation of a complex field induced by a higher dimensional operator was originally considered in the context of Affleck-Dine baryogenesis [8]. As we will discuss below, this rotation has important cosmological consequences, including the possible production of gravitational waves and dark matter.
In general, the initial rotation is not perfectly circular-it contains both angular and radial motion. If the radial motion is not damped, it can come to dominate the energy density of the universe. Here we assume that P couples to the thermal bath in a way that eventually allows for its thermalization and the removal of this possible moduli problem. We will discuss an explicit realization of this thermalization mechanism in Sec. 5. Upon thermalization, the rotation becomes circular (up to thermal fluctuations). Even if the thermalization occurs via scattering with U (1) P -charged particles so that the rotation could lose its charge, most of the charge is maintained in the rotation, since that is the state with the least free energy [9].
In this paper, we consider the scenario where thermalization occurs sufficiently early so that the motion is circular during the era of tachyonic instability.
The circular motion then evolves as follows [9]. The angular velocity is determined by the equation of motion. Long after the beginning of the rotations, when V (S)/S H 2 with V (S) the potential energy of S,θ 3) The evolution of S and the equation of state of the rotation may be determined from charge conservation, where R is the scale factor of the universe. For example, if V (S) is quartic and S f φ , where ρ θ ≡θ 2 S 2 is the energy density associated with the rotation. If V (S) is nearly quadratic and S f φ , Once the saxion reaches its minimum, S f φ , charge conservation requires Then the energy density of the rotation ρ θ scales as so-called kination [37,38].
In the following subsection, we discuss the consequences of axion rotation when the axion couples to a dark gauge field.

Dark photon production
If the global symmetry has a mixed quantum anomaly with a gauge interaction, the axion couples to the corresponding gauge field. We consider the case where an anomaly is present for a dark U (1) D gauge interaction, resulting in a coupling where e D is the dark gauge coupling and F is the field strength of U (1) D . The mixed quantum anomaly arises from the coupling of P to fermions ψψ simultaneously charged under U (1) D and global U (1) P , The equations of motion for the two transverse modes of the dark photon, A ± , are where m A is a possible dark photon mass and H is the Hubble expansion rate. The dark photon has a tachyonic instability centered around the momentum if m A < k TI . Inside the instability band, the growth rate of the dark photon fluctuation is as large as k TI . The global charge n θ is transferred to the dark photon helicity density.
If the axion field has negligible initial velocity in field space, it begins oscillating when H ∼ m φ . The sign ofθ then alternates with a frequency m φ , so efficient tachyonic instability requires k TI > m φ ∼θ. By the above equation, it can be seen that unless e D ∼ 4π, the tachyonic instability is ineffective. This conclusion may be avoided if the axion-dark photon coupling is much larger than Eq. (2.8), so that k TI is comparable to m φ at the beginning of the oscillations. Such an enhancement is possible by utilizing the Kim-Nilles-Peloso mechanism [39] and its clockwork generalization [40][41][42][43][44], as considered in [45].
On the other hand, if the axion field is initially rotating, the sign ofθ remains the same.
Also, for S > f φ , the angular velocity evolves as in Eq. (2.3), and it may decrease more slowly than the Hubble expansion rate. As a result, the growth rate via tachyonic instability, as set by the k TI of Eq. (2.11), may eventually outpace the Hubble expansion. The consequence is explosive dark photon production-without resorting to an enhanced axion-dark photon coupling via baroque model building. Note that once S = f φ ,θ decreases in proportion to R −3 as in Eq. (2.7), and the tachyonic instability will not become effective if it has not already done so.
Once k TI > H, the dark photon density is exponentially amplified, seeded by quantum fluctuations. This continues until back-reactions stop the tachyonic instability. This can occur either when the angular velocity decreases substantially via the loss of n θ stored in P or when the production of axion fluctuations from dark photons changes the angular velocity.
The former can occur when an O(1) fraction of the energy of the rotation is transferred into the dark photon. Computation of the scattering rate between A and the rotating field, taking care to include Bose enhancement for A , shows that the latter is relevant at a similar time. We assume that the production of dark photons from the axion rotation is saturated where r p is an O(10) factor that allows for the exponential amplification. The subscript p indicates evaluation at the time of dark photon production so that S p ,θ p , and T p denote the saxion field value, the angular velocity, and temperature at production. Determining the precise spectrum of the dark photons requires a numerical lattice simulation beyond the scope of this paper. Motivated by the detailed study of tachyonic instability in other contexts [20,22], we assume that when H = H p , the momentum distribution of dark photons is sharply peaked at k TI and the dark photon energy density is as large as the energy stored in the rotations at the onset of the tachyonic instability, ρ θ .

Axion rotations and fluctuations
We now discuss the fate of the axion field after the dark photon production. In this paper, we assume that the axion is massless but comment on the possibility of a massive axion and axion dark matter in Sec. 6.
Even after the production of dark photons via tachyonic instability is saturated, non-zero angular momentum n θ remains in the axion rotation. The rotation should not dominate the energy density of the universe during Big-Bang Nucleosynthesis (BBN) or after recombination. As we will see in Sec. 5, this can impose a non-trivial constraint. We parameterize the residual rotation by where Y θ,i is the yield of the global charge before tachyonic instability becomes effective.
We expect that r θ = O(1) if the dark photon production is saturated when S is close to f φ , since as argued above Eq. (2.12), tachyonic instability rapidly becomes ineffective after S reaches its minimum. Even for S > f φ at the time of production, it is plausible that after the saturation of the production of dark photons, the tachyonic instability ceases when the scattering between axions and dark photons makes the rotation incoherent so that . We take r θ as a free parameter, leaving the rigorous determination of r θ by a numerical simulation to future work.
The produced dark photons can create axion fluctuations via the coupling in Eq. (2.8). As long as the axions are light, as assumed here, the axions produced in this way at most have the same energy as dark photons and contribute to dark radiation, which will be discussed around Eq. (3.7).

Effects of U (1) D -charged fields
In the above discussion, we have integrated out the U (1) D -charged fermion ψ and ignored its dynamics. However, it could conceivably be produced in several ways, as enumerated below. Should this occur prior to the era of tachyonic instability, electric conductivity from the fermions [46][47][48] acts as a friction term in the equation of motion for the gauge field and may prevent efficient tachyonic instability from occurring. It is important to ensure none of these production modes are present.
1. The fermions can be produced from scattering between the particles in the thermal bath. This can be avoided if either y ψ S T so that the production of ψ is kinematically forbidden, or the couplings of ψ with the bath are sufficiently small that the production is inefficient.
2. When the motion of P is not completely circular, the mass of ψ oscillates and ψ may be produced by parametric resonance. However, as long as y ψ S θ , the oscillation of the mass is adiabatic, and production by parametric resonance is exponentially suppressed.
3. The fermions could be produced by the Schwinger effect from the dark electric field E D [49,50], or by the chiral anomaly from the dark helical electromagnetic field [48,51,52]. This production is exponentially suppressed if the fermion mass y ψ S is larger than (e D E D ) 1/2 at the time of gauge field production.
In the concrete realization discussed in Sec. 5, we impose the conditions listed above but find that only those in 1) give relevant constraints because T θ , E 1/2 D . When we identify the produced dark photons as dark matter, we also impose these conditions even after the dark photons are produced. This ensures that energy is not removed from the dark photons by electric conductivity or via the fermion production discussed in 3).
When the condition in 3) is satisfied, this also ensures that the mass of ψ is larger than k TI , justifying integrating out the fermion. It is therefore consistent to use the dimensionfive operator in Eq. (2.8). After ψ is integrated out, the dark gauge field interacts via a dimension-eight operator suppressed by the mass of ψ called the Euler-Heisenberg term [49].
It is natural to ask whether these self-interactions cause any difficulties. However, the energy density from the Euler-Heisenberg term is much smaller than the energy density from the kinetic term as long as y ψ S < (e D E D ) 1/2 . Therefore, it is not expected that this term causes any difficulty for efficient tachyonic instability. Furthermore, while self-scattering via the Euler-Heisenberg term may affect the spectrum of the dark photons, in the concrete realization we discuss in Sec. 5, we find that the scattering rate is smaller than the Hubble expansion rate, so it is irrelevant.

Deviation from coherent circular motion
In the above analysis, we have assumed that the rotation is coherent and circular. Indeed, in this paper, we require that thermalization of the rotation occur at a temperature above T p so that the motion is nearly circular by the time of tachyonic instability.
First, we comment on a potential subtlety. After the beginning of the rotation and before the completion of thermalization, the rotation is not yet circular. During this phase, it is possible that parametric resonance [53][54][55][56][57] may produce fluctuations of P with wave numbers ∼ V (S) [11,13]. If there is an epoch where parametric resonance becomes efficient before the completion of thermalization, the motion of P is no longer coherent and has fluctuating components with a field magnitude ∼ S. However, thanks to charge conservation, the field P still rotates as a whole. Once thermalization occurs, the fluctuations are dissipated, and the field value of P is nearly homogeneous within a horizon (up to thermal fluctuations) and continues to rotate. Tachyonic instability will then proceed as described above.
We now make a few brief comments on what might happen in the case where thermalization has not occurred by the era of tachyonic instability. As we now discuss, it is plausible that tachyonic instability still occurs, which might allow a relaxation of our requirement that thermalization occur prior to tachyonic instability. Of course, thermalization would in any case ultimately be necessary to avoid a moduli problem. First we consider the case that parametric resonance has not occurred before the epoch of tachyonic instability. In this case, the rotation is coherent, but owing to non-thermalized motion in the radial direction, θ oscillates with a frequency ∼ V (S). Despite this, modes with k V (S) only feel the time-averaged angular velocity θ ∼ V (S) so that tachyonic instability still occurs.
We have confirmed this by numerically solving the equation of motion of the dark photon.
We next consider the case where the parametric resonance becomes effective before the time of tachyonic instability, but thermalization has not yet occurred. Although the field value of P fluctuates, we expect that the dark photons with wave numbers V (S) feel the spatial-and time-averagedθ, which is O( r V (S)) because of charge conservation, where r is the ratio between the U (1) P charge density and the number density of the radial direction. We expect that tachyonic instability still occurs.
Despite these plausibility arguments for an effective tachyonic instability even in the presence of delayed thermalization, for simplicity, in our analysis, we will require thermalization to occur before the dark photon production by tachyonic instability.

GRAVITATIONAL WAVES FROM ROTATIONS
Gravitational waves are created at the epoch of explosive vector field production [17,20,21]. The GW equation of motion is given bÿ where M Pl = 2.4 × 10 18 GeV and Π T T ij represents the transverse-traceless anisotropic stress tensor. Importantly, Π T T ij receives a contribution from the quadrupole of the energy of the dark photon, which sources the GWs.
The peak frequency of the GW f GW is set by that of the dark photons at the time of creation. In Eq. (2.12), we defined r p as the constant of proportionality that relates the dark photon wave number to the horizon. We expect that the GW wave number at emission is given by k GW ≈ k TI ≡ r p H p . This redshifts to the present day, resulting in a peak frequency f GW = 2 nHz r p 20 1 2 k TI 3 × 10 −13 eV 1 2 10 g * (T p ) 1 12 . (3.2) Here, we have normalized to a frequency relevant for a signal at NANOGrav [19] or SKA [28].
The GW frequency determines the temperature at production Gravitational waves are initially produced with energy density ρ GW = M 2 Pl ḣ ijḣij /4. This energy density redshifts as that of radiation ρ R ; the ratio ρ GW /ρ R is invariant. Examination of Eq. (3.1) allows us to evaluate this quantity at the time of production as with k TI the wave number at production as defined above. In the second equality, we have noted that the energy density in dark photons at production ρ A is expected to be of the same order of magnitude as the energy in the rotation prior to the tachyonic instability, Here the subscript onθ p emphasizes that this quantity should be evaluated at the time of the tachyonic instability production. This results in a present day GW energy density at the peak frequency f GW .
In our explicit realization in Sec. 5, we show a case in whichθ p m S with m S the saxion mass. Here in the second (third) line, we have normalized the equation to a possible signal at NANOGrav (and at DECIGO/BBO), but this scenario can instead give rise to signals visible at other pulsar timing arrays such as SKA or space-based GW detectors such as LISA. To observe a signal in a similar frequency range but with a smaller amplitude than that suggested by NANOGrav, it may be necessary to efficiently subtract "background" events from super massive black hole mergers [58,59].
For the dark photon spectrum sharply peaked at k TI , the GW spectrum is peaked around the frequency f GW in Eq. (3.2) with a magnitude in Eq. (3.5). It decreases in proportion to f 3 for frequency f < f GW [60], and is sharply cut off above f GW . The scattering of the dark photon and the axion will make the spectrum broader [20,22]. The scattering is expected to be more efficient than the case of axion oscillations as we will discuss in Sec. 4.
Also, the polarization of the GWs could differ from the oscillation case. Here, a linear combination of the U (1) P charge in the axion rotation and the helicity of the dark photons is conserved, so strong polarization is expected to persist even in the presence of scattering and backreaction. Finally, should non-negligible production of dark photons from the axion rotation occur even after the transfer of O(1) fraction of the axion energy, contrary to our assumption, the spectrum will become broadened towards higher frequencies.
From the relationship between k TI and the angular velocity given in Eq. (2.11), the dependence of the dark gauge coupling on the GW signal is as follows: For a fixed magnitude of GWs and smallerθ p , the field value S p is required to be larger according to Eq. (3.4). The coupling of P with the thermal bath then becomes weaker, and henceθ p is bounded from below for the thermalization of the radial mode to be successful.
We will see this explicitly in Sec. 5. For frequencies of GWs relevant for pulsar timing arrays, a small dark gauge coupling constant is required.
Depending on the mass of the dark photon, it may act as either dark radiation or dark matter. The possibility of dark photon dark radiation from tachyonic instability is also discussed in Refs. [16,17,22]. In the (nearly) massless case, the contribution to N eff is intimately tied to the size of the GW signal. This is because the abundance of GW is related to the dark photon production. We can write Axion fluctuations produced from the dark photons can also contribute to dark radiation, but their abundance does not exceed that in Eq. (3.7). Without a detailed lattice simulation dedicated to the case of an axion with rotation, such as the one done in the non-rotating case [22], one cannot more precisely determine the relative energy transferred from the rotation to the GWs and the dark photon. For larger GW signals, a more precise computation could be of interest to establish a firm prediction for CMB-S4 [61].

DARK PHOTON DARK MATTER FROM ROTATIONS
If the dark photon is massive, 1 it is an attractive dark matter candidate. Since the GWs are associated with dark photon production, the present day dark photon number density is related to the size of the GW signal.
It is important to check that dark photon dark matter does not thermalize with the Standard Model bath following its production via tachyonic instability. Otherwise, the would-be dark matter abundance is depleted down to a thermal abundance, and while such a population would generally contribute an acceptably small contribution to N eff , it could not be cold dark matter. This constraint depends on details of the implementation, including, e.g., the masses of the U (1) P -charged fermions integrated out to generate the coupling between the axion and the dark photons. In our specific realization in Sec. 5, we discuss this in more detail.
The spectrum of dark photons are initially sharply peaked at the momentum k TI . The spectrum can, however, evolve via the scattering by the axion-dark photon coupling in Eq. (2.8). The effect of the scattering is investigated with lattice simulations in Refs. [20,22] for the case of axion oscillations, and it is found that the dark photon momentum is approximately conserved up to the redshift. In our case, there are two important differences.
First, there exist massless excitation modes around a background of a nearly circularly rotating P even when S > f φ ; these would correspond to axions at the minimum of the potential S = f φ . In the oscillation case, the dark photon scattering rate at low momenta is suppressed by the axion mass, but no such suppression occurs here. Second, the effective decay constant of the axion ∼ S continues to decrease after the production as shown by Eqs. (2.5) and (2.6) so that the scattering rate may increase, even after accounting for the dilution of the dark photon by Hubble expansion. For these reasons, we expect that the scattering can be efficient and may change the spectrum appreciably.
The evolution of the spectrum of the dark photons may affect the predictions of the mass of the dark photons as well as their warmness. We first discuss the predictions neglecting the possible scattering as a baseline. We then qualitatively discuss the effects of the scattering and how the predictions may be altered, leaving a rigorous discussion via numerical simulations to future work.

The case without scattering
We first estimate the warmness and mass of dark photons neglecting the scattering of dark photons following their production by tachyonic instability. The amount of the GWs fixes the fractional energy density of dark photons before they become non-relativistic at T = T NR , after which the fractional energy density increases in proportion to T −1 . Requiring dark photons to be dark matter, we can predict T NR using Eq. (3.4), with ρ DM /s = 0.44 eV the observed dark matter co-moving energy density.
Notably, the warmness exclusively depends on (the square root of) the magnitude of the GW signal. It is this quantity Ω GW h 2 that sets the initial density in dark photon dark matter. If this initial density is larger (smaller), the dark photon must spend a smaller (larger) amount of time redshifting as matter to ensure that the right dark matter density ρ DM /s will be obtained. This, in turn, is directly related to the warmness of the dark matter.
The current bound on the warmness arising from the Lyman-α measurements [63] indicates that the dark matter should become non-relativistic by a temperature Future measurements relying on observations of 21-cm radiation should be sensitive up to T NR = 100 keV [64]. This warmness constraint is shown in Fig. 1, where it shows up as the red solid curve-as expected from Eq. (4.1). The small deviation from horizontal at low frequency is due to changes in g * . This figure indicates the dark photon dark matter scenario is inconsistent with a putative NANOGrav signal. However, a future GW signal from the SKA-or a satellite-based detector such as DECIGO-could be consistent with this framework.
Notably, a signal observed at SKA or DECIGO would result in dark photon dark matter with a finite residual velocity that could leave an imprint on structure formation that might be observed by future 21-cm line observations. However, we emphasize that these conclusions may be modified by scattering in some cases, as discussed in the following subsection.
We obtain from Eqs.  10 -6 ↑ dark matter too warm without particle cascade to IR ↑ warmness probed by 21-cm without particle cascade to IR Correlations between GW signals and dark photon dark matter. In the limit that scatterings subsequent to production are neglected, the region above the red solid curve is excluded in the dark photon dark matter scenario because the produced dark photon would be too warm.
Scatterings among dark photons and axions can relax this constraint-see Sec. 4.2. Above the red dashed curve, the warmness can be probed via future 21-cm measurements [64]. The predictions of the dark photon mass are given by the black dashed lines. If scatterings increase the dark photon number, smaller dark photon masses will be predicted, and these contours should be regarded as the maximum possible values m max A . The Standard Model g * (T ) is used, and the quantity r p defined in Eq. (2.12) is fixed at r p = 20. The warmness constraint and the mass contours change with r p according to Eq. (4.1) and Eq. (4.3). The sensitivity curves of the gravitational waves are provided by Ref. [58] for EPTA [65], NANOGrav [27], SKA [28], LISA [29], (B-)DECIGO [30,66], BBO [31], and ET [67].

Effects of scattering
Although the coupling between the axion and the dark photon in Eq. (2.8) is suppressed by a loop factor, the dark gauge coupling, and a decay constant, dark photons and axions can still interact efficiently because of the enormous number densities involved. Using kinetic theory [68][69][70], which is essentially a perturbative computation that incorporates the relevant Bose enhancement, the scattering rate involving N i initial and N f final state particles with a typical momentum k is is the number of axion-dark photon vertices involved in the process.
f k is the occupation number of particles at the momentum k. Just after the production, scattering is fact expected to be efficient for any N V , and it is unclear which N V dominates.
A linear combination of the U (1) P charge in the axion rotation and the helical density of the dark photon is conserved. This conserved charge is analogous to the conserved particle number in a cold atom system, and it allows the scattering to evolve the system towards the formation of a Bose-Einstein Condensate (BEC) via particle cascade [71,72]. It is conceivable that the BEC-and its associated charge-could be in the form of the axion rotation, in which case dark photon dark matter would not be viable. However, for a fixed U (1) P charge, the free-energy density of the system is smaller if the charge is stored in small momentum modes of helical dark photons rather than the rotating axion. Thus, we expect that the charge is predominantly stored in the form of helical dark photons. As the system evolves towards a BEC, these dark photons cool, and the dark matter warmness constraint may be relaxed relative to the case without scattering presented in Sec. 4.1. Since the helicity density is conserved, the prediction on the dark photon mass remains the same.
A numerical lattice computation is necessary to follow the evolution of the dark photons and axions. As the system evolves towards low momenta, the dark photon occupation number increases, and the rate in Eq. (4.4) becomes larger for larger N V ; a perturbative computation based on the kinetic theory fails. If the scattering rate is suppressed at low momenta by a non-perturbative effect, which happens in the case of real scalar fields with a quartic potential [72,73] order interactions. However, given the importance of non-perturbative effects, as referenced above, we cannot rule out this possibility. In fact, this behavior has been observed in a lattice simulation of relativistic real scalar fields with a quartic potential [72]. To conclusively settle the fate of this non-helical component, a lattice simulation is necessary.
We emphasize that the number density of the non-helical component is not conserved as a whole. This opens the possibility that the number density of the non-helical component could in fact dominate the helical one. For this to occur while satisfying energy conservation, number-increasing processes would need to be effective either at the initial stage (before the evolution towards a BEC occurs) or in the high momentum regime while the evolution towards a BEC occurs in a low momentum regime. If the dark matter is primarily provided by these non-helical dark photons, the prediction for its mass becomes smaller.
The mass in Eq. (4.3), which was predicted by the helical density, can then be understood as an upper bound on the dark photon mass. Because of the cooling of the non-helical dark photons during the evolution towards the BEC by number-conserving processes, the warmness constraint can be relaxed compared to that from Eq. (4.1).
If instead the dark photon is massless, it contributes to dark radiation. Despite the possible scattering, the estimation in Eq. (3.7) is generic. It still applies as long as the axion does not continue to produce dark photons after T p .

CONCRETE REALIZATION
We have thus far discussed the production of dark photons and gravitational waves in as model-independent a fashion as possible. However, it is of interest to construct a concrete model and discuss a complete cosmology as an existence proof.
We explore the parameter space motivated by the potential GW signal reported by NANOGrav [19] and also for a smaller signal strength within the reach of SKA [28]. We then discuss a possible future GW signal at much higher frequencies, e.g., relevant for LISA and DECIGO/BBO, and discuss implications for dark matter in that part of the parameter space. In all cases, there are two particularly important questions: the approximate form of the saxion potential and the mechanism by which the saxion is thermalized.

Scalar potential
To produce an observable amount of GWs, the energy density contained in the axion rotation must be large at the time of tachyonic instability. The rotation of the complex field P begins when the Hubble expansion rate is comparable to the saxion mass. Since the energy density does not redshift until the onset of the rotation, a smaller saxion mass, namely, a flatter saxion potential, is advantageous for enhancing GW signals.
We first discuss the possibility that the saxion potential is approximately quartic at large field values, λ 2 |P | 4 , as might be expected in the case where the theory does not possess supersymmetry. We will show that a successful tachyonic instability is difficult in this case.
For an initial field value S i , the initial saxion mass in the quartic potential, ∼ λS i , is bounded by the Hubble scale during inflation H I . The bound on the tensor to scalar ratio from the CMB then imposes λS i < 6 × 10 13 GeV [74]. For a given strength of a GW signal, this translates to a bound on the quartic coupling as λ 10 −3 (10 −12 /Ω GW h 2 ) 1/4 (20/r p ) 1/2 .
However, for a given value of the quartic coupling λ, the coupling y ψ between the saxion and U (1) P -charged fermions y ψ Pψψ is bounded by the requirement that there be no excessive quantum correction to λ. Therefore, we expect y 4 ψ 16π 2 λ 2 . This bound on y ψ in turn imposes an upper bound on the mass of the U (1) P -charged fermions at the time of the tachyonic instability. This bound is sufficiently strong that U (1) P -charged fermions are unavoidably kinematically accessible at the time of production, m ψ < T p , and there is a danger that there is a thermal population of these fields.
A thermal population of such fermions is problematic because their presence leads to a non-trivial dark electric conductivity of the universe [46][47][48], which acts as a damping force in the dark photon equation of motion. This impedes the amplification needed for an effective tachyonic instability. One might have thought that tiny y ψ might allow only a sufficiently small abundance of the fermions to be produced, and hence a sufficiently small conductivity. However, y ψ cannot be taken vanishingly small. Most importantly, a sufficiently large fermion mass is required to suppress the production of fermions via the analog of the Schwinger effect [75]. We find the lower bound of on y ψ coming from this consideration leads to an excessive freeze-in abundance of the fermions, and hence an excessive electric conductivity. This consideration is robust, even in the potential presence of a large self-interaction amongst the U (1) P -charged fermions that reduces the conductivity.
For this reason, we now concentrate on a supersymmetric scenario, where the flatness of the saxion potential is radiatively stable. We focus on the case where the saxion has an approximately quadratic potential with a typical mass m S . This can be realized in a model with a global U (1) P symmetry spontaneously broken by dimensional transmutation from the running of the soft mass [76], and can also be realized in a supersymmetric two-field model with soft masses, where X is a chiral multiplet whose F -term fixes the global symmetry breaking fields P and P along the moduli space PP = V 2 P . The fermionic superpartner of the axion, the axino, obtains a mass from supersymmetry and/or R symmetry breaking. In the model described by Eq. (5.2), the tadpole term for X generated via supergravity effects, λXV 2 P m 3/2 , induces a vacuum expectation value X ∼ m 3/2 /λ, and the axino obtains a Dirac mass ∼ m 3/2 . When the saxion is thermalized, axinos are also thermalized; this leads to production of axinos. For the relevant range of gravitino masses discussed later, the thermalized axinos would have to be diluted by entropy production, e.g., from the decay of long-lived particles, to avoid overclosure. In this case, the analysis on the thermalization of the rotation in Sec. 5.2.2 should be modified by taking into account the dilution.
For sufficiently large m 3/2 , we may instead introduce R-Parity violation so that the axino with a mass ∼ m 3/2 can decay into Standard Model particles. And for m 3/2 > O(1) TeV, there is yet another possibility: the axino may decay into the Lightest Observable-sector Supersymmetric Particle (LOSP). Notably, this case is compatible with gravity mediation for the visible sector soft masses. As we will see in Sec. 5.2.2, such a large m 3/2 is consistent with a high-frequency region motivated by LISA, DECIGO, and BBO, but not with a lowfrequency region motivated by NANOGrav and SKA. We comment on these possibilities again in Sec. 5.3, where we discuss baryogenesis from the axion rotation.
However, for simplicity, we will focus on the model in Eq. (5.1), where the axino may remain sufficiently light that even its thermal abundance may not cause any conflict with observations, and no dilution is required. In this model, the axino mass is generated by a one-loop quantum correction from the coupling of P withψψ in Eq. (2.9) and is given by where A ψ is the scalar trilinear term associated with the coupling, V ⊃ y ψ A ψ P ψ ψ. In this expression, we have assumed the mass ofψψ is dominated by the supersymmetric mass term y ψ P . At minimum, the A term receives a contribution from anomaly mediation [77,78], eV [79]. To achieve this light axino, y ψ and/or m 3/2 must be sufficiently small.
However, the coupling y ψ cannot be arbitrarily small. The potential in Eq. (5.1) in terms of the parameters of theory is

5)
2 A coupling between the supersymmetry breaking field Z and P , ψ, andψ in the Kähler potential, e.g., P P † ZZ † , also generates an A term after Z obtains a vacuum expectation value, as required in gauge mediation, but we find this contribution is typically negligible.
where m P and m ψ are the soft masses of P and ψ/ψ at a UV scale µ. The second term is the quantum correction that arises from loops of the U (1) P field ψ. To successfully break the U (1) P global symmetry by dimensional transmutation, the soft mass must be smaller than the supersymmetric mass, i.e., m ψ < y ψ f φ . Furthermore, the saxion mass is dominated by the quantum correction so m S y ψ m ψ /4π. Taken together, these conditions imply For simplicity, we approximate the saxion potential at S > f φ by an exactly quadratic one with a mass m S , ignoring the logarithmic correction (except for the discussion on the Q-ball formation by a thermal potential around Eq. (5.15)), but we will impose the condition of Eq. (5.6).
Successful dimensional transmutation also imposes a constraint on the underlying soft supersymmetry breaking parameters of the ultraviolet theory. Since y ψ is bounded from above to suppress the quantum correction to the axino mass in Eq. (5.4), and m 2 P ∼ −y 2 ψ m 2 ψ log(y 2 ψ f 2 φ /µ 2 )/(8π 2 ), the soft mass m 2 ψ must be much larger than m 2 P for a large gravitino mass. For example, for m 3/2 ∼ 1 GeV, m ψ > 10m P is required. This could be achieved by a modest hierarchy in the suppression scales of the couplings of the supersymmetry breaking field with ψ and P .
To give the axion a rotational kick, we consider an explicit U (1) P symmetry breaking superpotential term, where Λ is a cut off scale. The scalar potential of P is then given by where H is the Hubble scale, c H is a dimensionless parameter, and A is a dimensionful parameter. In the case that supersymmetry breaking is mediated to P via supergravity, we expect scalar masses and A-terms to both be of order the gravitino mass m 3/2 .
We assume that c H during inflation is negative and O(1), or |c H | 1. For the former case, the field value of the saxion S ≡ √ 2|P | during inflation is fixed at a large value with V (S) = 0 [32,80], where the second and the third terms of Eq. (5.8) balance with each other. For the latter case, assuming that the initial field value of the saxion is large, e.g., around the fundamental scale, the saxion field value follows a classical attractor solution with V (S) ∼ H 2 [81]. 3 The angular field value may be initially random but is homogenized by inflation. We assume that c H after inflation also has the property described above.
After inflation, as the Hubble scale decreases, the saxion field value tracks a point with [32,80,86]. Once the Hubble scale becomes smaller than m S , the P field begins rotations with a kick induced by the A-term. One can show that the ratio between the potential gradient in the angular direction and that in the radial direction is O(A/m S ) (e.g., Ref. [11]). If A ∼ m S , as expected if the supersymmetry breaking soft terms of P are mediated by gravity, the rotations have O(1) ellipticity. In what follows, we assume m 3/2 m S to achieve O(1) ellipticity. Since the field value of P is homogenized by inflation, the rotation is nearly coherent.

Allowed parameter space
The main result of the paper is shown in Fig. 2. In the white regions, the axion rotation, In these cases, we may assume the dark photon mass is significantly smaller than that in Eq. (4.3) so that the dark photon contributes to dark radiation at a level given by Eq. (3.7). Then another dark matter candidate is needed; one example is the QCD axion. 4 Furthermore, due to small f φ in the top panels, the U (1) Pcharged fermions ψ are necessarily light, and they may have millicharge, Q ψ , via a kinetic mixing χ. The supernova constraint [87] requires Q ψ 10 −9 for m ψ 10 MeV. However, given that the millicharge Q ψ χe D /e EM is already suppressed due to small e D in Eq. (3.6), the resultant bound χ 0.05(nHz/f GW )(m S /300 keV) 1/2 is very weak.
On the other hand, in the bottom right panel, dark photon dark matter is consistent with the warmness constraint shown in Fig. 1 for a GW strength within the reach of DECIGO [30] and BBO [31]. The remainder of this section is devoted to explaining the origin of the exclusion regions. TeV for DECIGO/BBO. Since the saxion field value at production is always greater than or equal to its value at the minimum of the potential, S p ≥ f φ , andθ m S , the expression Eq. (5.9) taken with ρ θ (T p ) =θ 2 p S 2 p bounds f φ as a function of m S . This is shown by the blue regions of Fig. 2.
While m S ≤ f φ is required for the perturbativity of the saxion potential, a stronger condition is provided by Eq. (5.6). The precise relation between f φ and m S depends on y ψ . This in turn will depend upon the requirement of avoiding fine-tuning in the axino mass, as in Eq. (5.4). This requirement excludes the orange regions of Fig. 2. Taken together, these constraints show that the saxion must be light.
After the production of dark photons by tachyonic instability, the axion field may continue rotating. Its energy density behaves as matter for S > f φ and kination for S f φ according to Eqs. (2.6) and (2.7), respectively. The energy density of the axion rotation should be subdominant throughout BBN to prevent modification of the expansion rate of the universe.
Below the red dashed lines in Fig. 2, the energy density of the axion rotation is larger than the equivalent of ∆N eff = 0.4 for some temperature between T = 1 MeV and 0.1 MeV, and we expect light elements abundances may be modified to an unacceptable level. Here we show this constraint for r θ = 1 and 0.1 in Fig. 2, where r θ ≤ 1 is defined in Eq. (2.13) and parameterizes the residual amount of U (1) P charges in axion rotations after the dark photon production by tachyonic instability. We note that the upper bound on N eff during BBN dominantly comes from its effect on the neutron-proton conversion around T ∼ MeV and its impact on the helium abundance. When the axion rotation behaves as matter between T = 1 MeV and 0.1 MeV our requirement that the axion energy density be subdominant during this entire period comes from the lowest temperature, T = 0.1 MeV, by which time the neutronproton conversion is largely complete. In this case, our constraint is likely somewhat too strong; a more precise determination of the constraint is beyond the scope of this paper.
In deriving such a bound, the contribution to the energy density from thermalized axions, saxions, and axinos (∆N eff 0.1) and dark photons produced by tachyonic instability in Eq. (3.7) should be also included.
The thermalized saxions will eventually decouple and then decay to axions with a rate The contribution to dark radiation is negligible if the decay occurs before the decoupled saxions become non-relativistic. The gray region in the bottom panels of Fig. 2 is excluded because the decay happens only after the saxions become non-relativistic and the dark radiation abundance is enhanced.
The tachyonic instability can be obstructed if there exists a dark electric conductivity for the dark photons provided by dark-charged particles. While the precise dark content is model-dependent, at minimum, there exist particles ψ that appear in the Lagrangian with coupling y ψ Pψψ and are integrated out at low energies to realize the interaction in Eq. (2.8).
If these particles are sufficiently heavy as not to be thermally populated at T p , then this dangerous dark electric conductivity for the dark photon is avoided. In the entire parameter space of interest, ψ particle's mass, m ψ = y ψ S(T ), can be made sufficiently large by using a Yukawa coupling y ψ less than unity and consistent with the upper bound from the quantum correction to the axino mass in Eq. (5.3).

Thermalization constraints
We consider saxion thermalization via the scattering of P with new vector-like leptons.
We add a pair of vector-like Standard Model charged leptons L andL as follows: The coupling y L will be responsible for the thermalization of the saxion field. When the vector-like fermions are in the thermal bath, the thermalization rate of P is Γ th α 2 y 2 L T, (5.12) which decreases more slowly than the Hubble expansion rate does. We require thermalization to occur before L falls out of thermal equilibrium, i.e., T = m L , and also before the g * (T th ) 100 1 4 , (5.13) with the thermalization temperature T th = max(m L , T p ). 5 Production temperatures for the different benchmarks are given below Eq. (5.9). The lower bound on m L from collider search at the LHC is 790 GeV when the dominant decay mode is to the third generation leptons [90]. For concreteness, we fix m L = 1 TeV in our analysis.
While this constraint favors large y L , there are a number of upper bounds of y L that we will discuss below. Tension between these upper bounds and the lower bound from the thermalization requirement will lead to the green regions in Fig. 2.
The fermions responsible for thermalization must be thermally populated themselves.
This means that they must have a mass lower than the temperature at thermalization, If this condition is violated, the mass of L is always larger than the temperature for T > T th since S ∝ R −3/2 and T ∝ R −1 .
When the vector-like fermions are in the thermal bath, S obtains a thermal mass ∼ y L T .
If the thermal mass is larger than the vacuum mass m S , P rotates in the presence of a thermal potential withθ ∼ y L T > m S . The thermal potential is not exactly quadratic, but is flatter than a quadratic one. The consequence of such a potential is that Q-balls are formed [91][92][93][94][95]. The Q-balls evaporate into an inhomogeneous configuration of P once the thermal mass becomes subdominant. For evaporation to occur, the logarithmic correction to the nearly quadratic zero-temperature potential in Eq. (5.1), which makes the potential steeper than a quadratic one, is essential; otherwise, even a subdominant thermal potential makes the potential slightly shallower than a purely quadratic potential. The evaporation occurs when either y L T < m S or T < m L . For the latter case, however, the inhomogeneous configuration is not subsequently thermalized since L decouples from the thermal bath.
For successful thermalization of these would-be inhomogeneities, we require that the saxion thermal mass become smaller than its vacuum mass before L decouples. We further require that the Q-balls evaporate before dark photon production occurs so that tachyonic instability 5 Another lower bound on y L could arise from avoiding parametric resonance production of the scalar partner of L during the non-circular rotation of P . This can be avoided if y L S > m S , so that the mass of the scalar partner oscillates adiabatically. However, the scalar partner obtains a thermal mass ∼ T , implying that y L S < T is sufficient to achieve adiabaticity. Since T th > m S in the relevant parameter region, whenever y L S > m S is violated, y L S < T is satisfied.
is not disturbed. 6 In summary, for the Q-ball constraint, we require that the thermal mass be smaller than m S when T = max(m L , T p ) = T th , The coupling y L introduces quantum corrections to the axino and saxion masses. The axino receives a mass contribution from a loop of the vector-like fermions, wherem L is the soft mass of L, and the requirement that the saxion mass not be fine-tuned, ∆m 2 S < m 2 S , also leads to an upper bound on y L . Two-loop corrections from the vector-like fermions and weak or hypercharge gauge bosons necessarily introduce the following flavor-universal coupling, The current bound on the axion-electron coupling from observations of white dwarf stars [96] leads to an upper bound on y L . Two-loop corrections involving gluons (present if vector-like quarks are introduced) only induce vector couplings with quarks, which may be removed by baryon number rotations. 6 We expect that it may be possible to relax this condition. Since the thermal potential is of the form where F is some function, the radius of the Q-balls increases in proportion to the scale factor R after their production, and the Q-balls are not localized. In this case, we expect that the average angular velocity θ is given by V /S ∼ y L T and the tachyonic instability could still occur with this modified angular velocity.
After integrating out L andL, the axion couples to , In general, the coupling may violate the lepton flavor symmetry and induce rare processes such as µ → ea [7,97]. Also, the axion-electron coupling is constrained by the observations of white dwarf stars [96]. A muon-electron-photon coupling is also generated by the quantum correction from L and weak gauge bosons, inducing µ → eγ. These bounds can also be evaded if L dominantly couples to the tau leptons. For this reason, we do not include these constraints in the figure but comment below on how they would impact the parameter space.
In summary, the green regions in Fig. 2

A baryogenesis connection
The dynamics proposed in this paper may have intriguing connections with baryogenesis.
The charge in the rotation of the axion field can be transferred into the particle-antiparticle asymmetry in the thermal bath, allowing an explanation of the baryon asymmetry [9,10].
In the concrete model studied in Sec. 5.2.2, the charge in the rotation is transferred into the chiral asymmetry of Standard Model particles via the coupling in Eq. (5.11), which is transferred into the baryon asymmetry via the electroweak sphaleron processes. Using the result in Ref. [10], the baryon yield is given by where T EW is the temperature when the electroweak sphaleron processes fall out of thermal equilibrium. The Standard Model predicts T EW = 130 GeV [98].
For the parameter choices relevant to Fig. 2, we find that the Yukawa coupling needed to reproduce the observed baryon asymmetry, Y obs B = 8.7 × 10 −11 [99], is given by In the regions to the left of/above the cyan lines of the bottom panels in Fig. 2, the baryon asymmetry is underproduced because the required y L in Eq. We note that relaxing the constraint shown by the orange shading expands the parameter space compatible with the baryogenesis scenario. In fact, for m S > TeV, we may take the axino mass ∼ m S above the LOSP mass, so that the axino decays into the LOSP via the coupling y L as discussed in Sec. 5.1. In this case, only a milder constraint m S < f φ is applicable, actually opening up a parameter space with m S > TeV in the lower two panels of Fig. 2. This is not possible in the upper two panels.
Finally, in the absence of flavor structure in y L , the cyan lines may be impacted. This is because the flavor dependent constraints of the last section impose stronger upper bounds on y L , which in turn impacts the baryon asymmetry via Eq. (5.23). Once these bounds are imposed, it is no longer possible to produce the full baryon asymmetry at all in the bottom right panel. In the bottom left panel, values of large f φ × m S ( 10 4 GeV 2 ) that previously were able to reproduce the baryon asymmetry can no longer do so without running afoul of these constraints.

SUMMARY AND DISCUSSION
Axion fields frequently appear in extensions of the Standard Model that involve spontaneous breaking of global symmetries. Due to their small masses, axion fields in general do not rest at the minimum of the potential in the early universe. Rather, they are expected to be dynamical and, as a result, may play meaningful roles in the cosmological evolution of the universe. In this paper, we investigated consequences of an axion rotating in field space.
Such a rotation is expected when the radial direction of the global symmetry breaking field takes on a large initial field value and the global symmetry is explicitly broken by higher dimensional operators.
When a rotating axion couples to a dark photon, it modifies the dispersion relation of the dark photon, and the dark photon can become tachyonic for a particular range of momenta.
Dark photons are then efficiently produced by tachyonic instability. The produced dark photons generate gravitational waves whose frequencies may be of interest to current and future GW detectors such as NANOGrav, SKA, LISA, DECIGO, and BBO.
The produced dark photon may be dark matter. The strength of the gravitational waves is correlated with the warmness of dark photon dark matter that can be probed by galaxy surveys and observations of 21-cm lines. The scattering among dark photons and axions can change the warmness, which should be investigated via a numerical lattice simulation.
Since this mechanism most naturally occurs in cases where the radial direction has a large initial field value, it is important to also consider the dynamics of this saxion. The radial motion should be dissipated since otherwise it eventually dominates the energy density of the universe and leads to a moduli problem or dilution of gravitational waves. We investigated the thermalization process in detail for a model where the global symmetry breaking field couples to new particles charged under the Standard Model gauge symmetry. In this case, thermalization is successful for a variety of GW signals, spanning a broad range of frequencies and strengths, including those that can explain the NANOGrav signal or be probed by future detectors such as SKA, LISA, DECIGO, and BBO.
The angular momentum of the rotation is transferred into particle-antiparticle asymmetry in the thermal bath through the sector responsible for the thermalization, which is subsequently transferred to a baryon asymmetry via electroweak sphaleron processes. The parameter space that explains the observed baryon asymmetry is also consistent with the thermalization requirement for a parameter region with high-frequency gravitational waves.
In this paper, we assumed that the axion is massless and becomes dark radiation. With a non-zero mass, the axion that induces tachyonic instability may itself be dark matter. One possibility is that the axion is produced via the conventional misalignment mechanism [100][101][102]. However, peculiar to the mechanism we consider, there are other possible axion dark matter production channels, where the axions are created with non-zero momenta and may be warm in contrast to the conventional misalignment mechanism. We leave detailed discussion of this possibility for future work, but here we outline three interesting possibilities: 1. The kinetic misalignment mechanism [12,13], where the axion energy originates from the kinetic energy of the axion rotation. In Refs. [12,13], the axion rotation was assumed to remain coherent and oscillates whenθ ∼ m φ . However, as pointed out in Refs. [103,104] and Ref. [105] in the context of the axion monodromy dark matter and relaxion respectively, parametric resonance arising from the axion self-interactions is important. The coherent axion motion turns into axion fluctuations at a rate ∼ m 4 φ /θ 3 with a momentum θ /2 [105]. This occurs before the oscillations begin. Since the energy of the axion per quantum θ /2, the yield of the axions Y φ ρ θ /(θs/2) = Y θ = r θ Y θ,i . 7 The axions may be too warm for a small decay constant and/or a small axion mass, where the growth of fluctuations occurs at a late time, which gives the axion momentum an inadequate opportunity to redshift. Estimation of the axion abundance from the kinetic misalignment mechanism requires a more detailed estimation of the fraction of the charge that remains in rotation after tachyonic instability, denoted by r θ in Eq. (2.13).
2. Parametric resonance from the radial oscillation mode [11,13,106], where the axion energy originates from the potential and kinetic energy of the radial mode. Since the energy of the axions per quantum is determined by the frequency of the zero mode motion, m S ∼θ, for a motion with O(1) ellipticity, the yield of axions Y φ ∼ ρ θ,i /(θs) ∼ Y θ,i . If r θ 1, this contribution dominates over that of the kinetic misalignment mechanism.
3. The production of axion fluctuations from the scattering of the dark photons produced by tachyonic instability. Since the momentum of the produced axion ∼ k TI θ , if O(1) fraction of the energy of the axion rotation is eventually transferred into the axion fluctuations in this way, the number density of axions produced in this channel, Y φ ∼ ρ θ,i /(ω(k TI )s), is much larger than those in 1) and 2). Here ω(k TI ) is the energy of an axion quantum with a momentum k TI that depends on the shape of the saxion potential.
Moreover, the scattering among axions and dark photons, as argued in Sec. 4.2, can change the spectrum of axions.
This work motivates future numerical studies. In this mechanism, the energy of the axion rotation is transferred into that of dark photons. We expect that the backreaction to the axion rotation stops the transfer, and some amount of energy remains in the axion rotation. Here we simply treated the residual amount of energy as a free parameter, but this should eventually be fixed by numerical lattice computations. These computations will also reveal the abundance of axions produced in the channels 1) and 3) discussed in the previous paragraph. This will allow the investigation of the interesting possibility of identifying the axion we discuss in this paper with the QCD axion. Also, we estimated the amount and peak frequency of the produced dark photons and gravitational waves analytically. However, a numerical lattice computation could precisely determine the spectral shape. In principle, this could then be contrasted with those of other production mechanisms. Such a simulation would also allow a more precise prediction on the correlation between GW signals and the warmness of dark photon dark matter or the amount of dark photon dark radiation. In this Appendix, we describe the cosmological fate of new particles beyond the Standard Model. We find that all new particles will either be dark matter, comprise harmless dark radiation, or decay before BBN or after it invisibly.
First, we consider gravitinos. Requiring their abundance not to be so large as to create a gravitino-dominated era bounds the reheat temperature. This temperature sets a starting point for understanding the rest of the cosmology. In particular, we will require this reheat temperature be larger than other scales relevant for the cosmology, including the temperature at the time of the tachyonic instability, T p .
We expect that the saxion mass and gravitino mass are comparable, as we imagine that supersymmetry breaking is mediated to the P sector via Planck-suppressed operators. As a reminder, for the case of a NANOGrav or SKA signal in the top panels of Fig. 2 1. The gravitinos are produced during inflationary reheating and later decay to the axinos.
The gravitinos are thermally produced from the scattering of Minimal Supersymmetric Standard Model (MSSM) particles such as gauge bosons and gauginos [107]. When m 3/2 is much smaller than the gaugino masses, the yield of thermal production computed in Refs. [108,109] can be approximated as where T R is the reheat temperature after inflation. Here, m g is the gluino mass. The A-term associated with that top Yukawa coupling is assumed to be of the same order as the gaugino masses, for which case the contribution from the A-term is negligible.
The expression for T R < m g is correct for a sufficiently heavy inflaton, where the nonthermal production is subdominant [110,111]. We also caution that Eq. (A.1) breaks down when T R m g ; gravitino production is typically dominated by processes at lower temperatures during inflationary reheating, and therefore the contribution from lighter binos can dominate over that of the gluinos despite a smaller gauge coupling constant.
If the gravitino decays to the axino and the axion before dominating the energy density, the contribution to dark radiation is negligible. For m 3/2 10 GeV, the decay occurs after matter-radiation equality and hence we require that the gravitino abundance be less than that of dark matter. For m 3/2 10 GeV, the decay occurs before matter-radiation equality and the gravitino energy density at this time should be below that of radiation.
We obtain an upper bound on the reheat temperature from Eq. (A.1), for Noting that m 3/2 ∼ m S , this equation indicates that T R can be chosen consistently higher than T p in Eq. the assumption of a radiation-dominated universe during thermalization is satisfied.
2. The axions contribute to dark radiation. In this paper, we assume that the axion is massless (although the axion may be dark matter via different production channels described in Sec. 6). In the case of massless axions, the remaining energy density in the coherent rotation ρ θ scales as matter, ρ θ ∝ R −3 , when S > f φ and as kination, ρ θ ∝ R −6 , after S f φ . This implies that ρ θ eventually redshifts away and becomes harmless as long as it does not change the energy budget around BBN and the decoupling of the CMB. In addition to this zero-mode contribution, thermalized axions and the axions from saxion decays both redshift as radiation and contribute to dark radiation; relevant discussions are around Eqs. (3.7) and (5.10).
3. The axino is assumed to be the lightest supersymmetric partner (LSP) and constitutes a small contribution to dark radiation. If the axino is lighter than O(eV) [79], even the thermal contribution of the axino is harmless. We impose this constraint, noting that axinos decouple from the bath when the fields responsible for thermalization leave the bath. Since we assume the axino is the LSP, other R-parity odd particles will eventually decay to them, as discussed above for the gravitino. For other R-parity odd particles, we will verify below that their decay occurs when their energy densities are subdominant compared to radiation, and thus do not produce excessive dark radiation. 4. The dark photino mass can be as light as m A e 2 D m 3/2 /16π 2 from the quantum correction in anomaly mediation [77,78] (see also [112][113][114]). 8 From Eq. (3.6) and taking m 3/2 m S , we note that the photino may be quite light, Even for the highest frequency we consider (the DECIGO/BBO benchmark of f GW = 0.1 Hz), this modest mass implies that as long as the dark photino number density is smaller than one-tenth of the thermal one, the dark photino energy density is negligible. Dark photinos are produced by annihilation of axions, saxions, axinos, or vector-like fermions, but since the production rate is suppressed by the decay constant, a loop factor, and the dark gauge coupling, we find that the dark photino number density is much smaller than the thermal number density in the entire viable parameter space.  8 The threshold correction from the dark fermion-sfermion loops is negligible and does not cancel the UV anomaly mediation. This is because in the model with dimensional transmutation, the Bµ term of ψψ is only given by anomaly mediation and is much smaller than y ψ f φ m 3/2 .
On the other hand, thermalization would be disastrous for dark photon dark matter-it will necessarily heat dark photons while reducing its number density. This can occur through the scattering with ψ, saxions/axions, or L in the thermal bath. We find that we may always take y ψ large enough so that the mass of ψ ∼ y ψ S is much larger than the temperature and ψ is never in the thermal bath. The scattering rate with thermal saxions/axions is ∼ (α D /4π) 2 k A T 4 /S 4 , where k A is the momentum of dark photons.
The rate turns out be always smaller than the Hubble expansion rate. The scattering rate with L is∼ (α D /4π) 2 y 2 L k A T 2 /S 2 . After requiring y L S < T , the rate is smaller than the rate with saxions/axions. 6. The fermion ψ charged under U (1) P and the dark gauge symmetry (and their superpartners) is kinematically inaccessible and thus never populated. In a small fraction of the parameter space, it is also possible that ψ are populated after tachyonic instability and the freeze-out abundance of it explains the observed dark matter abundance.
We recall that the coupling to the dark photons was induced via the presence of the fermions ψ charged under the dark gauge symmetry, with an anomaly under the global symmetry, see Eq. (2.9). At T p , when the tachyonic instability takes place, the U (1) Pcharged particles must be absent from the thermal bath to avoid excessive electric conductivity that prevents efficient tachyonic production. As discussed in Sec. 5.2.1, we find that a proper choice of y ψ < 1 allows the fermion mass to be higher than the temperature at T p , thereby suppressing their production during this crucial era.
However, since their mass depends on the S field value and redshifts, these U (1) P -charged particles could, at least in principle, come into thermal equilibrium at a later time. We have found that there are choices of y ψ , consistent with the bound on the axino mass Eq. (5.4) such that this does not occur so that these particles are never produced.
Nevertheless, in the case where y ψ is small enough that production is not always forbidden, it is an interesting question whether they could be the dark matter. The U (1) P -charged We first discuss the case of the NANOGrav or SKA signal where the dark photon is not a viable dark matter candidate due to warmness-unless a cascade towards IR occurs due to scattering-and an alternate dark matter candidate would be especially welcome. To have sufficiently cold dark matter, m ψ 10 keV is required [115], which requires f φ 0.1 GeV to achieve the correct freeze-out abundance. Here we borrow the warmness constraint on sterile neutrino dark matter produced from oscillations of active neutrinos, since the momentum distribution for that case is expected to be similar to the freeze-out case. As can be seen from Fig. 2, in the region with the thermalization constraint, f φ ∼ 0.1−1 GeV.
Then, achieving the right freeze-out density would require m ψ ∼ 10 keV. Since S p ∼ f φ in that parameter region, m ψ = y ψ S is necessarily smaller than T p 5.5 MeV(27 MeV) for NANOGrav (SKA) and is already produced from the thermal bath at the time of GW production, leading to excessive electric conductivity and hence spoiling the tachyonic instability. Therefore, ψ cannot be the dark matter in this case.
While in most of the parameter space in the bottom panels of Fig. 2, ψ is too abundant if its production is kinematically allowed, around the region 10 MeV < f φ 10 GeV of the bottom left panel that is excluded by BBN if r θ > O(0.1), the freeze-out abundance of ψ can explain the observed amount of dark matter while ψ is cold enough.
When the motion of P is not completely circular, the mass of ψ oscillates and ψ may be produced by parametric resonance. However, as long as y ψ S is much larger than m S , which is the case in the viable parameter space, the oscillation of the mass is adiabatic and the production by parametric resonance is exponentially suppressed.
7. The vector-like fermions decay into a Standard Model doublet lepton and a saxion/axion via the Yukawa coupling y L with a rate ∼ y 2 L m L /(8π) and the decay occurs much before BBN. The vector-like fermions mix with Standard Model doublet leptons with an angle ∼ y L f φ /m L . This mixing causes them to decay into a Higgs boson and a singlet lepton, which we assume is a tau, with a rate ∼ y 2 L y 2 τ f 2 φ /m L . These decays occur before BBN for the parameter region in Fig. 2 For m 3/2 < GeV(m LOSP /TeV) 5/2 , the decay of the LOSP into a gravitino occurs before BBN and the resulting gravitino abundance may be below the dark matter abundance [116,117] so that the subsequent decay of the gravitino into an axino does not overproduce dark radiation however slow the decay is. In the benchmark points in Fig. 2, y L 10 −7 is required by Eq. (5.13), so the decay occurs well before BBN. The vector-like leptons then decay promptly as discussed above.
If m L > m B , decay occurs through the Kähler potential y 2 L P P † † /m 2 L . The decay rate is suppressed by the mass of the fermion component of , which we assume to be τ , where we assume m B ∼ m . To ensure decay before BBN requires large enough y L and/or f φ that are inconsistent with the low-frequency parameter regions, but the decay into a gravitino anyway occurs before BBN because of the small m S ∼ m 3/2 .
If the LOSP is the superpartner of the new fermion, L, it decays into + a with a rate ∼ y 2 L m L /(8π), which is much larger than in Eq. (A.5).
If the LOSP is the right-handed stau τ and m τ > m L , the stau decays into L + a with a The decay before BBN requires large enough y L and/or f φ .