Quantum Field Theory of Interacting Dark Matter/Dark Energy: Dark Monodromies

We discuss how to formulate a quantum field theory of dark energy interacting with dark matter. We show that the proposals based on the assumption that dark matter is made up of heavy particles with masses which are very sensitive to the value of dark energy are strongly constrained. Quintessence-generated long range forces and radiative stability of the quintessence potential require that such dark matter and dark energy are completely decoupled. However, if dark energy and a fraction of dark matter are very light axions, they can have significant mixings which are radiatively stable and perfectly consistent with quantum field theory. Such models can naturally occur in multi-axion realizations of monodromies. The mixings yield interesting signatures which are observable and are within current cosmological limits but could be constrained further by future observations.


Introduction
Close to 95% of the Universe is invisible. About a quarter of it is compressible, behaving as dark matter (DM). The remaining three quarters have negative pressure, are called dark energy (DE), and approximate the cosmological constant. To leading order these two components seem to be mutually non-interacting in the absence of gravity. The origin of their scales and their relative normalization remains mysterious, particularly in light of the cosmological constant problem [1,2,3,4]. While we expect to observe and study systematically the dark matter, DE is extremely difficult to explore beyond simply establishing its existence.
Whatever we think about the cosmological constant problem, nature solves it in some way. So one might wonder if the residual components of the dark sector might be some relics of the mechanism which stabilizes the vacuum energy. This might suggest that the dark degrees of freedom are not totally decoupled from each other. They might have interactions which transfer energy from one dark sector to another. This question is particularly appealing since such energy rebalances might be probed with direct observations of large scale structure and CMB in the forthcoming experiments.
It is easy to imagine DM/DE mixing [5], especially since one might resort to a fluid approximation to describe their influence on astronomical scales [6]- [10]. The question is whether such setups make sense in a realistic quantum field theory. One might be tempted to ignore this issue, arguing that we don't really know much about dark sectors. For example, some question whether QFT even works as usual in the dark sector, since after all it suffers from the cosmological constant problem. However, QFT has proven time and again to be the most reliable tool for the description of consistent interacting models in nature. So our starting point here will be to seek for a QFT description of interacting DM/DE.
We will outline the key obstructions to formulating an interacting DM/DE theory in QFT both conceptually and phenomenologically. The main conceptual obstacle stems from decoupling: one expects that DM should be heavy fields, and DE -if dynamical -light. Because of decoupling, in a normal perturbative QFT very heavy and very light fields do not interact very efficiently. If one introduces strong couplings by hand, these lead to the renormalization of the dark energy scales by the exchange of dark matter quanta. That would render the DE potential very steep, and so DE could not be a field in slow roll any more. The only way for such a potential to behave as DE is if the field is in the minimum with a nonzero cosmological constant. Further, if for any reason the dark energy quanta remain light, their exchange between dark matter particles yields extra long range forces between them, which can easily violate observational bounds on DM interactions.
The combined implications of these problems for interacting DM/DE models are quite severe for a broad class of proposals 1 . Combining the bounds from the weakness of long range forces and flatness of the light DE field potential after DM-induced radiative corrections are included prohibits significant DM/DE mixing when DM are heavy particles with masses m DM m DE . However, we will show that significant DM/DE mixings are quite likely if a fraction of DM is composed of ultralight boson condensates, with mass not too much greater than the current Hubble scale H 0 . Thus it is possible to have the mass of this part of the DM close enough to the mass of DE to evade decoupling. Using scaled down axion monodromies, which can arise, for example, in the 'axiverse' framework [11], we will see that a system of several axions can easily yield such a model with dark energy and a fraction of dark matter, with significant DM/DE interactions. We will focus on a set of three axions, of which two are very light, one being DE at the present and another, with a mass m > ∼ H 0 , being a component of DM. Since the amount of ultralight dark matter is bounded [12,13], in this case most of DM is provided by a heavier axion 2 , with a mass m H 0 . As a result, systems of coupled axions are completely realistic frameworks for interacting DM/DE, and so probing for DM/DE interactions becomes quite sensible. In addition, such tests may be a new portal for exploring the 'axiverse. ' 2 Problems with DM/DE Interactions So imagine a model of interacting DM/DE within QFT. The observations suggest that at cosmologically large distances any such model should be described by a weakly coupled effective field theory, involving an interaction Lagrangian L I (φ, ψ). Otherwise it would not even be possible to identify dark matter and dark energy as separate constituents of the universe. Let the field φ be the DE field (a 'quintessence' [6,14,15,16], which is 'fundamental' up to some cutoff scale µ H 0 ) and ψ a DM field. Depending on the mass of the DM field, dark matter is either DM particles (if the mass is large enough so that the particles are nonrelativistic from matter/radiation equality onwards) or an almost uniform, time dependent, condensate of a DM field zero mode (especially if the mass is small; but note that even in this case, m DM m DE ). Clearly, DM could also be a combination of many separate components.
Suppose first that DM is made up of sufficiently heavy particles. For simplicity 3 , we could take them to be fermions, with a mass term m ψψ ψ. To model DM/DE interactions, let the DM mass depend on the DE field expectation value, m ψ = m ψ (φ). Now, the DE field must not change too quickly, in order to act as dark energy. Further, since we are postulating that the underlying theory of DE is an effective field theory, to ensure its validity we require that φ/M Pl is at most of order unity. So we can expand Here g = c m 0 ψ /M Pl , and c is a dimensionless number. So unless c is exactly zero, the DM/DE interactions would involve a Yukawa term Clearly, prohibiting such a term would require imposing a symmetry in the DE sector.
In the absence of a symmetry, however, this terms yields new interaction channels. This coupling gives rise to DM/DM interactions depicted by the Feynman diagram of Fig. (1). The exchange of the virtual DE quanta generates an attractive long range force between a pair of DM particles, with an effective potential where we have used the fact that 1/M 2 Pl = G N . This expression, of course, is just the 3D Fourier transform of the Euclidean DE propagator, ( p 2 + m 2 DE ) −1 . Note that while the Yukawa suppression cuts the DE-mediated force off at distances r 1/m DE , at shorter distances the resultant force follows the inverse square law. Another important consequence of this term arises from the DE self-energy diagram of Fig. (2), which yields a 'dressing' of the DE potential and specifically the DE Lagrangian mass term renormalization due to the emission and reabsorption of the virtual DM quanta. The Lagrangian mass term of the DE field, controlling the steepness of the DE potential, is corrected additively by the contributions from the virtual DM quanta. The corrections can be calculated using the standard Feynman rules for the diagram (2), which yield the truncated scalar two-point function where / k = γ µ k µ , and k and p are the external and internal momenta, respectively, and m 0 ψ is the dark matter mass. Rationalizing the integrand by using the Feynman trick to combine the factors in the denominator, performing the trace and analytically continuing to d = 4 − dimensions leads to wreak havoc on the cancellation scheme, requiring a completely different prescription for the finite counterterms to maintain the cancellation. This is beside any concern one might have about the UV sensitivity of the theory, which we have completely set aside here. The point is that the sensitivity of the DE Lagrangian mass term to DM scales is a real thing, as long as the DM particles are much heavier than DE, we insist that the standard rules of QFT apply, and there are no additional symmetries between the scales m DE and m DM to cancel these corrections 4 .
What are the implications of these considerations? First of all, the DE field should be very light to simulate dynamical dark energy -a slowly varying scalar field. Typically this is realized by taking the mass m DE to be smaller than the current value of the Hubble parameter H 0 for the field to yield at least an e-fold of cosmic acceleration now. Hence, at all subhorizon scales, r < 1/H 0 , the DE field is effectively massless, giving rise to an additional long range force among dark matter particles, V DM −DM −c 2 G N m 0 ψ 1 m 0 ψ 2 /r. In the linearized limit, this behaves as a scalar correction to gravity, and its strength relative to the Newtonian potential is This quantity is constrained by the bounds on the weak equivalence principle violation between dark and ordinary matter [17]- [22] to be less than about 0.1. Therefore, Note that if we rewrote (10) in terms of the Yukawa coupling g = cm 0 ψ /M Pl , we would get 2 , which is essentially the parameter β of Eq. (18) of [5], introduced by the comparison of long range dark energy and gravitational fields sourced by dark matter agglomerates, if we ignore the variation of the DE field between the interior and the exterior of the DM agglomerate 5 .
Let us now turn to the quantum corrections to DE mass (9). As we noted above, DM contributes to the DE Lagrangian mass term regardless of any additional aspects of the UV sensitivity and hierarchy problems. To keep the DE field potential flat despite these corrections, one must require ∆m DE m DE . Then (9) and where in the last equation we substituted numerical values of the Planck scale and the Hubble scale. The point of this equation is to note that as the mass hierarchy between DE and DM increases, the cross coupling, as parameterized by c, becomes extremely small quickly. Again, this is perfectly reasonable since it is but a manifestation of decoupling of heavy and light modes, which is universally valid in QFT. Moreover, with the actual data reflecting the real world included, equation (12) shows that it is generically extremely difficult to couple any DM heavier that a milli-eV to DE without completely destroying the DE sector. We summarize (11) and (12)   Here DE is an ultralight quintessence field with a sub-Hubble mass.
So, for standard quintessence fields with masses H 0 , the DM/DE interactions are extremely suppressed. This renders simple interacting DM/DE models, with large m DM and tiny m DE , such as those described in [7], completely unrealistic. The DM/DE interactions must be irrelevant when m DM is large. The problem with this, however, is that the phenomenological bounds on DM masses coming from large scale structure [23] require that most of DM is considerably heavier than 10 −3 eV. Hence, if there are very light 'dark' degrees of freedom, which can interact significantly with quintessence, they cannot constitute much of the cosmic contents. As a result, they cannot affect DE evolution very much. Note that while we have used fermions as an example of DM, we would get exactly the same bounds if DM were bosonic particles with trilinear couplings to DE. Note also that we have taken here DE to be a scalar, as is clear from the fermion mass term ∝ψψ in (2). If φ had been pseudoscalar, coupling to the fermion bilinear ∝ψγ 5 ψ, the bounds from long range forces would essentially disappear since they would be spin-dependent and drop with a higher power of distance [24]. The bounds from the radiative corrections would however remain essentially the same as discussed here, still prohibiting the couplings of heavy DM to DE.
One might try to alleviate this by using a heavier DE field. That could help both with the bounds on the corrections to the DE potential and the bounds on the extra long range forces, since the mass corrections could be larger and the force Yukawa-suppressed at distances longer than 1/m DE . An example is provided by the so-called MaVaN models [25]. The original proposal assumed that DE is a single scalar with the mass parameter much larger than the Hubble scale, mixing with a visible and sterile neutrino. Although the scalar was too heavy to be in slow roll on its own, the coupling to neutrinos could allow the relic neutrino background to slow down its cosmological evolution. However, this relied on neutrino-DE couplings which were too strong and led to additional attractive long range forces at distances < 1/m DE , which, while shorter than the scale of the universe, were still too long range. This forced non-relativistic neutrinos to clump, forming nuggets of size ∼ 1/m DE , which would not suppress the DE rolling [26]. Indeed, this is in agreement with our calculations of |c| and V DM −DM , which show that with m DE H 0 one would get a DE-mediated attractive force much stronger than gravity at distances 1/m DE .
Evading this requires either lowering DE mass down to H 0 , or using much more involved DE models, as in the latter of [25] (see also [8]). These models utilize a hybrid inflation dynamics at extremely low scales, which drives current cosmological acceleration by an effective cosmological constant term (and assumes that all other contributions to the cosmological constant are somehow cancelled). It also includes a coupling of DE to a neutrino species which is always relativistic in order to erase the effects of strong DEmediated couplings. Yet, in all these cases, most of DM is left completely decoupled from DE. So while the interactions might lead to imprints in the DE evolution, the evolution of DM is not directly affected.

Axion Monodromies to the Rescue
The considerations above show that the main obstruction to DM/DE interactions comes from decoupling in QFT: light and heavy fields do not mix very well. So, this naturally suggests that the strongest interactions occur if DM and DE mass scales are not very different. Now, since DM must be non-relativistic, to accommodate this -i.e. have DM with ultrasmall masses and yet ensure it is not relativistic -one needs DM to be a condensate like DE. In other words, DM is a spatially homogeneous expectation value of an ultralight scalar field very much like quintessence, but with its mass somewhat larger than the current Hubble scale H 0 . In this case, instead of being stuck in slow roll, the DM vev is oscillating about its minimum. To leading order, its energy density behaves just like CDM.
To build models with such fields consistently in QFT, one needs to protect the theory from large radiative corrections, both from the ultralight modes in the dark sector and from any other heavier degrees of freedom. A simple way is to provide the low energy EFT with continuous shift symmetries, which ensure that the ultralight dark modes are only derivatively coupled in perturbation theory. Hence their potentials are automatically very flat and radiatively stable. In fact, the only fully self-consistent models of dark energy rely precisely on constructions involving pseudo-Nambu-Goldstone bosons as DE fields (pNGBs, or 'axions') [27]- [31], which utilize shift symmetry as the protection mechanism of the effective theory. So in practice all one needs is to add an extra axion, over and above the one playing the role of DE, and pick its scales so that it can be DM today [12,32,33]. However, if it is very light, this DM axion is constrained by cosmological data and cannot be all of DM today [12] (see also [13] for the updated bounds). Even so, it could comprise as much as ∼ 5% of total amount of DM, which is a far larger fraction of critical energy density than could be in, for example, sterile neutrinos and other light fermions.
This reasoning points directly to a possible framework for realizing such a scenario 6 : the 'axiverse' [11]. Multiple axions, with varying scales which control their couplings and masses, arise from string compacifications. The vast diversity of the multiverse and the proliferation of the pNGB fields suggests that such hierarchical axion models are readily realized. In such models the axions are naturally mixed, which means that in low energy EFT they couple to each other [34,35]. In fact, the simplest mechanisms to generate small axion masses based on axionic monodromies [36,37,38] automatically ensure that the interactions between the lightest axions are so strong that the axion mixing energy is comparable to the energies in the normal modes. Hence such models can leave quite strong cosmological signatures, affecting directly both DE and DM evolution, by rebalancing the energy contents of individual channels in the course of cosmic evolution.
In what follows we will analyze a simple model, which is essentially a combination of the ideas of [12] and [34,35,36] to demonstrate this. Ultimately, we will be led to a setup with three axions where one is heavy, another is lighter than the current Hubble scale H 0 , and the third has a mass in between but relatively close to H 0 . By decoupling, the interactions between the heaviest and the two light axions are small, and so we will ignore them completely. On the other hand, we will devise a simple setup where two light axions mix strongly, so as to ensure that the axion decay constants are not super-Planckian, while one still happens to be lighter than H 0 . This sector is basically a scaled down version of misaligned axion inflation of [36,39]. The underlying monodromy which generates ultralight axions in slow roll with sub-Planckian axion decay constants implies that the two lightest axions automatically have strong interactions. We will analyze the cosmological evolution of the system, both on the cosmological background and with perturbations, to identify the possible signatures and show that the model is consistent with the current bounds.

Decoupled Limit
To set the stage for the perturbative description of the multi-axion dark sector, we start with the case where the interactions are completely turned off [12]. So imagine three axions, each with a potential The unity here merely cancels the bare cosmological constant (by hand). The potentials arise from nonperturbative gauge dynamics and are radiatively stable: once flat, they are always flat [27]- [31], [38]. The radiative stability of the cosine potentials is ensured by the underlying shift symmetry φ i → φ i + C i , broken down to the discrete subgroup φ i → φ i + 2πN i f i by the non-perturbative effects, where N i are integers. The unbroken discrete symmetries ensure that the full nonperturbative potentials are harmonic, periodic functions, where the higher harmonics are suppressed by exponents of (M Pl /f i ) 2 . As long as f i 's are not very large the corrections can be ignored [40]. The homogeneous modes' dynamics are described by the equations of motion The fields for which ∂ φ i V i 3Hφ i are DM, with virialized kinetic and potential energy. Their homogenous mode oscillates rapidly about the minimum of its potential, with frequency m 2 i = µ 4 i /f 2 i , and, by the virial theorem, with an amplitude φ i (t) ∝ 1/a 3/2 . The fields for which ∂ φ i V i 3Hφ i are DE, being in slow roll. The standard approach to properly normalize the amplitude is to take each light axion to have an initial value f i set during early inflation by the inflationary thermal drift affecting all fields with masses smaller than H inf lation . Then, for decoupled axions, f i is the typical distance from the nearest axion minimum, taken as the origin in the field space. This axion stays in slow roll until H drops to ∼ m i , and after that starts to oscillate around its nearest minimum, with the amplitude diluting as ∼ 1/a 3/2 . If all (or most) of DM today is one such axion, its amplitude at present is given by 4 . If f i is fixed, the mass m i should be picked correctly to match this condition, accounting for the dilution during the evolution of the Hubble scale from H ∼ m i to H 0 .
On the other hand, the DE axion is still in slow-roll, frozen at its initial value φ i O(f i ). In this case one needs to have f i M Pl , µ 4 3M 2 Pl H 2 0 , in order to match the scale of DE today. Finally, a very light axion which has a mass above H 0 , but not by much, has been in slow roll until relatively recently. It started to roll around the minimum in a very late universe. This implies that at scales H > m intermediate , there was more dark energy than now, and is one of the reasons behind the bounds of [12]. Such an axion didn't dilute too much since it started rolling, but it cannot be more than about 5% of DM now [12,13]. A priori, its presence seems to require a fair bit of tuning. Yet, as we will see shortly, such tuning can be accommodated rather straightforwardly in monodromies used to generate ultralight mass scales without super-Planckian axion decay constants f i .

Mixing and Interactions
Let us now turn on the interactions between the two lightest axions. We use the effective potential generated by the axions' couplings to dark gauge fields as in a simple monodromy model [36,39], where we again subtract the cosmological constant term by hand, picking V = 0 at global minima. We now expand (14) to quadratic order around a minimum at φ 1 = φ 2 = 0, where M ij is the mass matrix We can easily transition to the system of normal modes, using a field space rotation by an angle tan 2θ = , which yields the light and heavy axion eigenmodes, respectively, with eigenvalues λ ± = [m 2 1 + m 2 2 ± (m 2 2 − m 2 1 ) 2 + 4m 4 12 ]/2. The light field l is DE, and the heavy one h is a DM component. This requires λ − 10 −33 eV. If λ + 10 −24 eV, then its abundance is limited to no more than 5% of DM, and depending on the mass, possibly even less. In this case, most of DM has to be something else, e.g. a third axion. On the other hand, if λ + 10 −24 eV, the heavier light axion is not constrained strongly and could be all of DM. For that case, one needs a large separation between the eigenvalues, λ + /λ − 10 18 . In any case, some hierarchy between mass eigenvalues is always needed.

Large Hierarchy, Planckian Decay Constants, Decoupling
To achieve a large hierarchy between the axion masses, we need det M (tr M ) 2 , or in terms of the Lagrangian parameters, In this limit, the eigenvalues become approximately λ − det M tr M , λ + tr M , and the mixing angle is tiny, so that l φ 1 , h φ 2 .
To simplify the initial analysis, we consider the case when all dimensional parameters are degenerate: µ 1 = µ 2 = µ 3 = µ, f 1 = f 2 = f 3 = f , and use the integer n to control the hierarchy. If one light axion is to be quintessence and the other heavy one is to be all of DM, then we need n ≥ 10 9 . Clearly, introducing one more axion will relax this. We will consider such cases later on.
The full potential reads In terms of the normal modes, V (2) = µ 4 2f 2 [l 2 + (2 + n 2 )h 2 ], and the quartic is, to leading order in n, Moving away from the degenerate limit, one finds other contributions, and in particular, a term ∼ h 3 l with a coefficient starting at O(n).
The bottom line, however, is that in this case the interactions are strongly suppressed by the hierarchy parameter n. The cosmological evolution is very close to the one of noninteracting dark matter and dark energy components. This is expected: this is just decoupling in action.
In this case, the only dimensionless number that controls the interaction between the light and heavy fields is given by the ratio of the masses, and so by the very hierarchy it must be small. So the dark matter is h, which will oscillate rapidly, and for sufficiently large n can be all of DM in the universe today.
The quintessence field l, on the other hand, stays in slow-roll. Since the dimensional coefficients are degenerate, the numerical values will be as in the single-axion case: µ 10 −3 eV, f O(M Pl ), and the amplitude of h is of order f /n. However, this shows that this case does not utilize field space monodromies since the axion decay constant f is O(M Pl ). All the dynamics unravels in the attractive basin of a single minimum in the theory.

Less Degeneracy, Sub-Planckian Decay Constants and Monodromies
More interesting cases can be realized with less degeneracy in the parameter space and by meeting the conditions for monodromy. To explore this, let us still take the axion decay constants f 1 = f 2 = f for simplicity. Since we imagine that they are determined by a symmetry breaking in some UV complete framework such as string theory, this is probably realistic anyway. Their magnitude is sub-Planckian, of order M Pl /S, with S 1 an action of the breaking sector. We will take f 0.01 − 0.1M Pl . Since the scales µ i are generated by non-perturbative physics, they are generically exponentially sensitive to S, µ i ∼ e −S i . To illustrate some interesting dynamics, we will focus on µ 1 µ 2 µ 3 without loss of generality. For calculational simplicity, let us further introduce a single parameter 1 and define µ 1 = µ, µ 2 = µ, µ 3 = µ/ .
The potential now reads (with V 0 chosen to cancel the vacuum energy in the minima) Expanding in around the vacuum φ 1 = φ 2 = 0, the eigenvalues are λ − In these variables, Note that when = 1 this reduces to the previous case. However, having introduced an extra dimensionless parameter , we can now explore the richer parameter space (n, 1/ ) to search for physically interesting regions. First off, we see that if we again require m h = µ 2 /( 2 f ) 10 −23 eV, then the heavy field oscillations can be all of dark matter with mass m h . To ensure that the other mass is ≤ H 0 , we need 10 −5 to get the requisite eigenvalue separation. This is not all: if the heavy field vacuum is h = 0, the two fields are essentially decoupled, but the light field l behaves as quintessence only if f M Pl , which we have been trying to avoid in order to guarantee full control over EFT.
However: if we pick a generic heavy field vacuum h = 2kπf / √ 1 + n 2 , with k 1 integer, keep f < M Pl , and then integrate h out, the effective potential for the light field becomes The frequency of the second oscillatory term is very high, and so we do not expand it. It represents merely a small modulation on top of the flat quintessence potential given by the first term, yielding small bumps in dark energy density. Note that the effective light axion decay constant is f l,eff nf , and so it is O(M Pl ) for n ∼ 10 − 100 if f 0.01 − 0.1M Pl . Further, the light field vacuum is not at zero anymore but at l vac = 2knπf / √ 1 + n 2 2kπf , which will be super-Planckian for k > 10 or so. This is how monodromy sets up large field excursions from respective vacua. Essentially, the heavy field pulls the light one's minimum far from the trivial one. Ergo, the slow roll flat plateaus are set up by entirely sub-Planckian local physics. The effective potential remains tiny, of the order k 2 µ 4 /n 2 , even when the field space distance to be traversed is super-Planckian.
For completeness: although DE and DM in this limit remain decoupled, we can easily get an approximate ΛCDM evolution from entirely sub-Planckian microphysics. To start, we need the equation of state to transition from w 0 at early times towards −1 in the future. The oscillations of h around its minimum provide all of the required DM. Writing h(t) = 2kπf / √ 1 + n 2 + χ(t), the field χ(t) oscillates around zero with frequency m h and amplitude decreasing as a −3/2 . From the potential above, it is clear that l is in slow roll today as long as k is big enough. The resulting energy density and pressure are ρ m 2 h χ 2 + 2µ 4 = µ 4 ( 1+n 2 4 χ 2 f 2 ) Clearly, this realizes the required limiting behavior: before χ dilutes enough, so that its energy density dominates, w total 0. As time goes on, w total converges to −1 for as long as l remains in slow roll. Specifically, under the assumption that the heavier axion is all of DM today, this sets the value of the heavy field displacement from the minimum at the present epoch More generally, this is the upper bound on the value of χ at the current time.

Mixing and Monodromies
The most interesting examples involve a very small hierarchy between the two lightest axions, while still realizing monodromy in the field space. This means that we need the third axion to be most of DM. To realize this case, we retain the parametrization of scales in the full potential of the previous section, but allow to be larger. As we noted previously, the masses of normal modes are m 2 l µ 4 Hence taking 1 and n ∼ O(few) will easily generate a small hierarchy between them. Further, the potential (22), which we repeat here, Since most of DM is another, third axion, much heavier than the normal modes l, h, we now do not require 1. We do require, however, that m l H 0 in order for the light mode to be quintessence. Further, as above, we demand that f < M Pl such that the EFT is under control, including the corrections from quantum gravity. We see that n ∼ O(few) will help bring the mass m l down by about an order of magnitude, and combining this with µ 10 −3 eV will ensure that m l ∼ H 0 , while m h is a few orders of magnitude heavier.
Next, we imagine that the heavy field h resides in an attraction basin of a generic vacuum h vac = 2kπf / √ 1 + n 2 . Shifting the heavy field to h = 2kπf / √ 1 + n 2 + χ, we rewrite The first term in (25) is the potential of the heavy field h -or χ, after the field redefinition.
In writing it, we have dropped the phase shift 2kπ. At times when m h > H, this potential forces χ to oscillate around the minimum, contributing to DM. Clearly, the oscillations are bounded, as shown by [12]. This simply means that, at present, the amplitude of χ cannot be larger than a fraction of (24) -i.e., if the amplitude of χ is less than a fifth of χ 0 in (24), the energy density contribution to DM from the heavier ultralight axion cannot be more than few percent. This would fit the bounds of [12,13]. The second term gives the leading order contribution to the light field potential. It features a monodromy. To see this we normalize the argument of the cosine by taking out the light field frequency as a prefactor, so this term is ∝ cos The phase shift 2knπ/(1 + n 2 ) generically cannot be removed by periodicity of the cosine, and for 1 < k < n, it is ∼ 2kπ/n < ∼ 2π. Yet, this automatically induces a displacement of l from its vacuum by 2knπf / √ 1 + n 2 2kπf . Even for relatively moderate values of k, this makes the effective field excursions of l super-Planckian although f M Pl . The third term is an additional modulation of the light field potential -as before -producing bumps on the quintessence potential. The reason we are adding it is that, in principle, it helps generate the hierarchy between l and h.
It is important to notice, however, that when n, 1/ are not extremely large, all the terms in the potential are normalized approximately the same. This means that if we expand (25) beyond the quadratic order, we will find nonlinear terms, describing interactions between DE (l-field) and a component of DM (h-field), with energy densities that are comparable in magnitude to the energy density contributions from the normal modes at times before the heavier normal mode dilutes away. In other words, this guarantees that the DM/DE interactions are significant during at least a brief period in cosmic history. In fact, we will show that the interaction can be strong enough to facilitate a classical transition of the heavier normal mode from one vacuum to a neighboring one, creating a domain wall, whose tensions are fortuitously small enough to meet the observational limits [41]. In those cases, we find the largest deviations from the ΛCDM cosmology. In the next section we resort to numerical evolution of the model in this regime to determine the signatures of the interactions. We also consider the effects on cosmological perturbations.

Cosmological Evolution
To investigate the cosmological signatures of the dynamics, we now turn to numerical analysis. First, we integrate the homogeneous equations to explore the evolution of the background and find out the behavior of the leading order cosmological observables as a function of the redshift z. After that, we consider the subleading effects using cosmological perturbation theory. Our focus is on the imprints of DM/DE interactions. We exhibit the similaritiesand crucial and interesting differences -of the interacting models with monodromy when compared to ΛCDM and models with ultralight decoupled axions as discussed by [12]. The interacting models may leave nontrivial imprints when compared to these benchmarks, especially when the interactions are strong enough to force a vacuum transitions of the heavier ultralight mode. After that we turn to perturbative analysis of the dynamics and find that interacting models meet the current observational limits. For the background evolution, we start our integration from high redshift, fixing the parameters to their ΛCDM values Ω c0 = 0.26, Ω b0 = 0.05, Ω m = Ω c + Ω b , T CM B = 2.275K, assuming 3 flavors of massless neutrinos. We normalize the Hubble parameter today to H 0 = 68 km/s/Mpc. We choose initial conditions for the axion potential such that the Universe is flat, and we set the initial field velocity to zero, as the axions are negligible at the initial epoch. A more complete analysis, including data and a broader allowed range of values of cosmological parameters, would clearly be interesting but is beyond the scope of the present work.

Background
To maximize the mixing between the two lightest axions, we choose parameters that ensure that one of the cosine terms goes over maxima before settling in a minimum 7 . The vacuum transition is most easily understood using the form of the potential given in Eq. (10): Let us now take µ 2 µ 1 , µ 3 (reordering the scales relative to the previous case) and f 1 = f 2 . With these parameters, φ 2 begins to oscillate while φ 1 is still in slow roll. If the oscillations of φ 2 are large enough that n|∆φ 2 |/f 2 > 2π, then the third term crosses over a maximum of the cosine. This means that, for given n and f 2 , the initial value of φ 2 must be a sufficient distance from the minimum of the full potential. This is easily accomplished, even with the constraint that the heavier light axion field is not more than ∼ 5% of DM today. Further, one should, in principle, also pick the initial conditions for the slowly rolling field so that it rests in a convex section of the cosine, to avoid excessive tachyonic perturbations 8 [31]. These conditions ensure that there will be mixing, with some effect on cosmic expansion. The parameters we choose are listed in table 1. For comparison, we also consider a decoupled model where the mixing between the two ultralight axions is completely turned off. The parameters describing it are listed in table 2.    Table 2: Parameters for the decoupled model. Figures (4-7) show the evolution of these models, normalized to ΛCDM. Both the decoupled and interacting models differ from ΛCDM by a temporary increase in H(z) at low 7 One can consider the interactions between fields that remain in a single vacuum, but in this case the effect of mixing is diminished and short-lived. This is because once the heavy field settles in its minimum the cosine is well approximated by a quadratic potential, and therefore the mixing term can be set to zero by a field redefinition. 8 For a short lived stage of late acceleration these bounds are relatively weak.
redshift (and the corresponding variation of d A (z)), Fig. 4. Initially, we normalize H/H ΛCDM to unity fixing the amount of CDM in the early universe to reflect the data. We add DE at early times in the form of the heavier ultralight axion, which will decay before today. This means that early on there was more DE. The equilibration of the heavier ultralight axion is delayed by its oscillations and also by the interactions. Nevertheless, asymptotically it removes some fraction of DE and converts it into DM as time goes on, which redshifts away under the influence of DE. In the asymptotic future, this implies that H/H ΛCDM will dip below unity. In both models, the oscillations of the intermediate-mass axion field produce oscillations in the total equation of state, Fig. (5), with corresponding effects in H and d A .
Interestingly, these oscillations can mimic a DM/DE interaction even in the decoupled case, when the actual interaction is zero. As the equation of state of the intermediate-mass field oscillates, it goes from being more DE-like to DM-like to DE-like, and so on. For a field that has m ∼ O(10)H 0 these oscillations are slow enough to be observed and do not effectually average out as they do for the main component of dark matter. Essentially, since the mass of the heavier light axion is so close to the Hubble scale, the kinetic and potential energies in the field oscillations do not have time to virialize. Their energy rebalance, as well as the fact that the field just fell out of slow roll, simulate the energy transfer between DE and DM, as if the interactions were present. Constraints on axions in this mass regime are discussed in [16], assuming that DE is a single axion.
A scrutiny of the evolution of the equation of state parameter with redshift reveals a striking difference between the coupled and decoupled models, shown in Fig. (5). This happens when the coupled model has interactions strong enough to force a vacuum transition of the heavier light axion. This jump and the subsequent field ringdown is responsible for differences in the shapes of the plots of H(z) and d A (z). While the vacuum transitions are clearly very interesting for inducing the largest allowed deviations from the ΛCDM evolution, at the same time, they are worrisome since they lead to the formation of domain walls inside the horizon. These domain walls could in principle lead to the large perturbations in the late universe, and are strongly constrained by observations [41]. To evaluate the constraints, we need to get a reliable estimate of the domain wall tension, which controls the scale of distortions. Since the transition is between nearest neighbor minima in the cosine potential governing the heavier light axion, we can approximate the cosine by the quartic potential, obtained by expanding the cosine about the maximum between the two minima. This yields where f ef f f 2 /n is the effective period of this cosine, as discussed after eq. (25) above. This means that the field vev in the core of the wall is ∼ f ef f , but crucially, the self-coupling is λ 4 (µ 3 /f ef f ) 4 . Its scaling with four inverse powers of f ef f is the reason the domain walls pass the observational bounds. Indeed, the tension of a domain wall separating two nearest neighbor minima is On average there would be one such domain wall per Hubble volume, and so the fraction of the domain wall energy to the total energy inside a Hubble patch would be of order Since µ 2 3 < M Pl H 0 -as the heavier light axion is just a component of the total energy density of the universe now, and the interactions are a subleading contribution to its energy density -we find This inequality 9 can be satisfied with the monodromy configurations by taking f 2 M Pl , which ensures that the domain wall-induced distortions of the cosmological geometry are small enough. With further expansion of the universe, the energy density in domain walls decreases with respect to dark energy, so its contribution remains subdominant. However, there could be interesting small effects if the bound (28) is close to saturation. We note that the energy density in the heavier field in the coupled model dissipates away faster than that of simple DM, as displayed in the evolution of the decoupled model, Fig. (6,4). This is due to the sharp increase in kinetic energy of the intermediate-mass field as it goes over the extrema, as shown in Fig. (7). During this stage, the kinetic energy dissipates much faster, approaching the 1/a 6 law. This constitutes an important signature of the coupled model: the coupled model may have a comparatively large impact on cosmic 9 We note that we are somewhat sloppy here, since we are ignoring the numerical prefactor in (28). This prefactor is easily O(10 −2 ) for realistic choices of parameters, and for n ∼ O(10) − O(100) it allows f 2 M Pl /10 to satisfy the bound (28). We have also neglected to scale the energy density in the domain walls up by a factor of redshift to the epoch when the walls were produced and compare that to the critical energy density then. Since the relevant factor goes as √ z, and the redshift when the domain walls are made is z < ∼ 100 this is readily compensated by other uncertainties.
expansion at a particular redshift and a comparatively small effect at later times. 10 Finally, we note that the fluctuations about the background satisfy a redshift-dependent dispersion relation. This is because the effective mass of the intermediate-mass axion is background dependent due to the interaction with the DE axion.

Perturbations
Finally, we turn to the linearized cosmological perturbations of the models with two light axions and cold dark matter. Our purpose is not a detailed study of the perturbation effects, i.e., a full analysis that would be needed for the detailed comparison with data. We merely want to demonstrate the consistency of the axion DM/DE models with the current bounds, by comparing their late time evolution to the standard ΛCDM model. This is the key consistency test, since at early times the matter contents and behavior is approximately the same in both models. We stress that some aspects of this have already been done. Specifically, an extensive analysis of the effects of a single ultra-light axion, which is a limit of the decoupled model with two ultralight axions we have been using here, has been done in [12,13]. A similar analysis of the coupled dynamics with nontrivial mixing would therefore seem very warranted but is left for future work.
Here we will work in synchronous gauge, comoving with cold dark matter, defined by the line element where the gauge variable η is given by and the additional condition δu m = 0 [43]. Here we decomposed the axion fields as φ I =φ I + ϕ I (I = 1, 2), splitting them into the background and the order perturbations, respectively. Then, the perturbation field equations are [44]: The conservation equation for matter isδ To complete the system of independent equations, we usë We set initial conditions deep in the matter dominated era [45], choosing, without any loss of generality, z = 20 as the initial time for the numerical integration of the perturbation equations. Given an initial δ m (t in , k), the last equation translates intoδ m (t in , k) = H(t in )δ m (t in , k). The initial field perturbations can be set to zero, as they will quickly reach the attractor solution. Numerical results are given in Fig. 8, 9. We find that in our model there is a small rescaling of the matter power spectrum with respect to the ΛCDM cosmology, with a suppression at higher wavenumbers, as illustrated in Fig. 8. By matter perturbations here we mean those in the main component of the dark matter. Since the heavier axion is still very light, it will not cluster on the scales probed by the present galaxy surveys.
The evolution of the field perturbations is shown in Fig. 9, where we plot the fraction δρ ax /ρ ax as a function of redshift and wavenumber. An interesting observational signature of the spectra are the enhanced oscillations around a particular scale, particularly in the interacting model. While the detailed implications for cosmological data are yet to be understood, we see that at least at present it appears that the multiaxion models do provide viable alternatives to ΛCDM.

Summary
The dark sectors of the universe could be totally separate components, referred to as dark matter, behaving as completely neutral non-relativistic billiard balls, and dark energy, modeled by a single number, the 'cosmological constant'. Such a simple picture, with carefully tuned normalizations, fits the observations very well at the present time. Yet, possibilities for other options remain open.
The 'cosmic coincidences' problem, namely the uncanny similarity between the fractional contributions of very different cosmic constituents, remains an open question. Why would dark matter and dark energy influence the bending of cosmic geometry comparably at present, even though the time evolution of their effects is dramatically different today? What sets their normalizations? Could the conundrum of cosmic coincidences be resolved by some interaction between the hidden sectors of the universe, and could those also influence the visible contents? There are many very simplistic schemes attempting to address these questions by proposing to relate various sectors with interactions that could equilibrate between cosmic constituents at the level of GR+fluid sources.
At the microscopic level, the situation is not so simple. The problem rests on the fact that the evolution of different cosmic components is controlled by very different scales. If these components are allowed to interact, quantum mechanical effects communicate the presence of these scales from one sector to another. This basically precludes any significant interactions between heavy dark matter and light dark energy. Specifically, the dark energy-generated long range forces between dark matter particles and the quantum radiative corrections to the dark energy mass induced by virtual dark matter particles constrain the cross-couplings to essentially zero when the masses of DM and DE are very different.
Yet, there is a possible way out of this impasse: it involves a somewhat more complicated, but consistent, set of models of multliple axions. Some of them are light enough to be DE, some are heavy enough to be (mostly decoupled) DM, and some are too heavy to be DE today, but nevertheless can couple strongly to DE. Setups like this emerge naturally in the construction of radiatively and non-perturbatively consistent models of field-driven cosmic acceleration, using monodromy to explain the origin of super-Planckian field displacements in effective field theory. To be clear, this does not solve the cosmic coincidence problem. However, if DE and DM are axions with comparable decay constants, and if their masses have correct values, then during inflation they will get the right initial conditions to be DM and DE today. Also, in such cases we easily find significant interactions between DE and a component of DM, with interesting observational implications.
From the low energy point of view, this is a self-consistent procedure for designing a field theory which could have a viable UV completion. But it might be more: the constructions involving many light, mixed axions could be a signal of the presence of the 'axiverse', a string-theoretic realization of a low energy theory with many light axion-like fields. Such constructions make the presence of many mixed axions more natural.
We find that at the present time such models are basically degenerate with ΛCDM but have significant deviations which could be looked for in future observations. The lightest axions can mix significantly to generate such signatures. Among the specific predictions, we note oscillations in the equation of state, sharp transitions in the Hubble parameter around a particular redshift, and possible signatures of domain walls produced by the classicallyinduced phase transitions of dark matter fields. These are consistent with current data but could lead to interesting small signatures at very large scales. We believe that such specific predictions warrant a dedicated analysis with a detailed comparison to data.