Kaluza-Klein FIMP Dark Matter in Warped Extra-Dimensions

We study for the first time the case in which Dark Matter (DM) is made of Feebly Interacting Massive Particles (FIMP) interacting just gravitationally with the standard model particles in an extra-dimensional Randall-Sundrum scenario. We assume that both the dark matter and the standard model are localized in the IR-brane and only interact via gravitational mediators, namely the graviton, the Kaluza-Klein gravitons and the radion. We found that in the early Universe DM could be generated via two main processes: the direct freeze-in and the sequential freeze-in. The regions where the observed DM relic abundance is produced are largely compatible with cosmological and collider bounds.


Introduction
The nature of Dark Matter (DM) and its interactions remain an open question in our effort to understand the Universe. Up to now, the only evidence about the existence of such dark component is via its gravitational effects. It could well be that DM has no other kind of interaction and, thus, it will be undetectable by current and future particle physics experiments. Moreover, in such a case the reheating temperature needs to be quite high (typically 10 13 GeV for DM mass of 1 MeV) in order to generate the observed DM relic abundance via a purely gravitational interaction [1][2][3][4], given the value of the Planck mass, m P ∼ 10 19 GeV, which determines its strength. This is true, however, only if we live in a four-dimensional space-time: in extradimensional scenarios, the gravitational interaction may be enhanced, either because the fundamental Planck scale in D dimensions is m D m P (as in the case of Large Extra Dimensions (LED) [5][6][7][8][9]), or due to a warping of the space-time which induces an effective Planck scale Λ in the four-dimensional brane such that Λ m P (as in Randall-Sundrum models (RS) [10,11]), or by a mixture of the two mechanisms (as it occurs in the more recent ClockWork/Linear Dilaton (CW/LD) model [12][13][14][15]). As it is well known, this feature of the extra-dimensional scenarios has been advocated as a solution to the so-called hierarchy problem, i.e., the huge hierarchy between the electroweak scale, Λ EW ∼ 250 GeV, and the Planck scale, which would generate corrections of order of the Planck scale to the Higgs mass. These corrections would destabilize the electroweak scale unless either an enormous amount of fine-tuning is present or the Standard Model (SM) is the ultimate theory, which seems unlikely given the questions that are not explained within this framework (for instance, neutrino masses and baryogenesis, besides DM itself). In the extra-dimensional models mentioned above, the large hierarchy between the electroweak scale and the fundamental (or effective) Planck scale is eliminated, since the latter can be as low as O(TeV).
As a consequence of such lower Planck scale in extra-dimensional models (either fundamental or effective), the gravitational interaction is enhanced, and a DM particle with just such interaction could become a WIMP, that is, a stable or cosmologically long-lived weakly interactive massive particle, with mass typically in the range 100 -1000 GeV, and whose relic abundance is set via the freeze-out mechanism. This possibility has been thoroughly studied in the framework of the RS scenario [16][17][18][19][20][21][22][23][24] and in a series of recent papers that study generic spin-2 mediators [25][26][27][28]. It has also been considered in the context of the CW/LD model [29].
In this work we again explore the RS framework for DM, yet analyzing a different scenario in which the relic abundance of DM is set via the so-called DM freeze-in production mechanism [30][31][32][33][34] (for a recent review see Ref. [35]). In this case DM is a feebly interacting massive particle (FIMP), so that it never reaches thermal equilibrium with the SM thermal bath, and as a consequence its abundance remains smaller than the equilibrium one along the history of the Universe. More specifically, here we focus on the sub-case of ultraviolet (UV) freeze-in [36] for which the temperature of the thermal bath is always lower than the scale of new physics, which in our model is the effective Planck scale in the 4-dimensional brane, Λ, at which the gravitons become strongly interacting.
In our setup we assume that both the SM and the DM particles are localized in the same 4-dimensional brane, and by definiteness we consider real scalar DM, only. We relax the request for the RS model to solve the hierarchy problem, and allow Λ to vary in a wide range (Λ ∈ [10 2 , 10 16 ] GeV) to fully explore the parameter space that could lead to the correct DM relic abundance via freeze-in from a purely phenomenological perspective. In order to have a consistent model, we stabilize the size of the extra-dimension by using the Goldberger-Wise mechanism [37], which generates the required potential for the four-dimensional radion field. Then, besides the interaction through Kaluza-Klein (KK) gravitons, we also take into account that the SM and DM species can interact with the radion. We consider both SM particle annihilation into DM through KK-gravitons and the radion (direct freeze-in), as well as production of DM from out-of-equilibrium KK-gravitons and the radion (sequential freeze-in). We solve numerically the relevant Boltzmann equations in all cases and also provide analytical approximations for the final DM relic abundance in different ranges of the temperature, useful to understand our main results. We always work within the sudden decay approximation for the inflaton, and shortly comment on how our findings would be affected by a non-instantaneous inflaton decay.
We vary the DM mass, radion and KK-graviton masses and the scale Λ, determining the reheating temperature T rh which leads to the correct DM relic abundance in each case, within the validity range of our effective four-dimensional theory. We find that in this scenario the observed DM density can be generated even with a reheating temperature lower than the electroweak scale. Recall that the only constraint on T rh is that it has to be higher than the Big Bang Nucleosynthesis temperature of around a few MeV [38][39][40][41][42].
The outline of the paper is as follows: in Sec. 2 we briefly remind the main features of the RS scenario; Sec. 3 is devoted to the analysis of DM production via freeze-in within our model, both via direct and sequential freeze-in; finally, in Sec. 4 we present our conclusions, while the relevant interaction rates used in our calculations are collected in the Appendix A.

Theoretical Framework
In this Section, we shortly review the Warped Extra-Dimension scenario (also called Randall-Sundrum model [10]) and introduce our setup to include DM in the model.

A (very) Short Summary on Warped Extra-Dimensions
The popular Randall-Sundrum scenario (from now on RS or RS1 [10], to be distinguished from the scenario called RS2 [11]) consider a non-factorizable 5-dimensional metric in the form: where σ = k r c |y| and the signature of the metric is (+, −, −, −, −). In this scenario, k is the curvature along the 5 th -dimension and it is O (M P ). The length-scale r c , on the other hand, is related to the size of the extra-dimension: we only consider a slice of the space-time between two branes located conventionally at the two fixed-points of an orbifold, y = 0 (the so-called UV-brane) and y = π (the IR-brane). The 5-dimensional space-time is a slice of AdS 5 and the exponential factor that multiplies the M 4 Minkowski 4-dimensional space-time is called the "warp factor". Notice that, in order to have gravity in 4-dimensions, in general η µν → g µν , with g µν the 4-dimensional induced metric on the brane. The action in 5D is: with M 5 the fundamental gravitational scale, G (5) and R (5) the 5-dimensional metric and Ricci scalar, respectively, and Λ 5 the 5-dimensional cosmological constant. As usual, we consider capital Latin indices M , N to run over the 5 dimensions and Greek indices µ, ν only over 4 dimensions. The reduced Planck mass is related to the fundamental scale M 5 as: where M P = m P / √ 8π 2.435 × 10 18 GeV, being m P the Planck mass. We consider for the two brane actions the following expressions: where f IR , f UV are the brane tensions for the two branes, L SM and L DM the SM and DM Lagrangians densities, respectively. Dots in eq. (2.6) stand for any possible new physics on the UV brane and, thus, decoupled from us.
In RS scenarios, in order to achieve the metric in eq. (2.1) as a classical solution of the Einstein equations, the brane-tension terms in S UV and S IR are chosen such as to cancel the 5-dimensional cosmological constant, f 4 IR = −f 4 UV = −24M 3 5 Λ 5 . Throughout this paper, we consider all the SM and DM fields localized on the IR-brane, whereas on the UV-brane we could have any other physics that is Planck-suppressed. We assume that DM particles only interact with the SM particles gravitationally. 1 Alternative DM spectra (with particles of spin higher than zero or with several particles) will not be studied here. Notice that, in 4-dimensions, the gravitational interactions would be enormously suppressed by powers of the Planck mass. However, in an extradimensional scenario, the gravitational interaction is actually enhanced: on the IR-brane, in fact, the effective gravitational coupling is Λ = M P exp (−k π r c ), due to the rescaling factor √ G (5) / −g (4) . It is easy to see that Λ M P even for moderate choices of σ. In particular, for σ = k r c 10 the RS scenario can address the hierarchy problem. From a purely phenomenological perspective, here we will work with Λ = [10 2 , 10 16 ] GeV, relaxing the requirement that the RS model should provide a solution to the hierarchy problem.
Any 5-dimensional field φ µν can be written as a KK tower of 4-dimensional fields as follows: being χ n (y) the wave-functions of the KK-modes along the extra-dimension. The equation of motion for the n th KK-mode is given by: where m n is its mass. Using the Einstein equations we obtain [43]: from which: being J 2 and Y 2 Bessel functions of order 2 and z n (y) = m n /ke σ(y) . The N n factor is the n th KK-mode wave-function normalization. In the limit m n /k 1 and e kπrc 1, the coefficient α n becomes α n x 2 n exp (−2k π r c ), where x n are the are the roots of the Bessel function, J 1 (x n ) = 0, and the masses of the KK-modes are given by: (2.11) Notice that, for low n, the KK-modes masses are not equally spaced, as they are proportional to the roots of the Bessel function J 1 . At large values of n, on the other hand, the roots of the Bessel function become approximately x n = π n + 1 4 + O n −1 . In this limit, the KK-modes masses are approximately equally spaced (as in LED and the CW/LD scenarios) and proportional to a characteristic length scale R such that: m n k π e −k π rc n = n R , (n 1) (2.12) The normalization factors can be computed imposing that: In the same approximation as above, i.e. for m n /k 1 and e k π rc 1, we get: (2.14) Notice the difference between the n = 0 mode and the n > 0 modes: for n = 0, the wave-function at the IR-brane location y = π takes the form whereas for n > 0: χ n (y = π) = k r c e k π rc = √ r c e k π rc M 3/2 5 Considering the coupling between KK-gravitons and brane matter (being h M N the graviton field), we have: from which is clear that the coupling between KK-graviton modes with n = 0 is suppressed by the effective scale Λ and not by the Planck scale.

Adding the Radion
Stabilizing the size of the extra-dimension to be y = π r c is not easy. Long ago it was shown that bosonic quantum loops have a net effect on the border of the extra-dimension such that the extra-dimension itself should shrink to a point [44][45][46]. This feature, in a flat extra-dimension, can only be compensated by fermionic quantum loops and, usually, some supersymmetric framework is invoked to stabilize the radius of the extra-dimension (see, e.g., Ref. [47]). A popular mechanism implemented in RS models to stabilize the size of the extra-dimension was proposed in Refs. [37,48] and can be summarized as follows: if we add a bulk scalar field S with a scalar potential V (S) and some ad hoc localized potential terms, δ(y = 0)V UV (S) and δ(y = π r c )V IR (S), it is possible to generate an effective potential V (ϕ) for the four-dimensional field ϕ = f exp (−k π T ) (with f = 24M 3 5 /k and T = r c ). The minimum of this potential can yield the desired value of kr c without extreme fine-tuning of the parameters.
The S field will generically mix with the graviscalar G 55 (notice that the KK-tower of the graviscalar is absent from the low-energy spectrum, as they are eaten by the KK-tower of graviphotons to get a mass due to the spontaneous breaking of translational invariance caused by the presence of one or more branes). On the other hand, the KK-tower of the field S is present, but heavy (see Ref. [49]). The only light field present in the spectrum is, then, a combination of the graviscalar zero-mode and the S zero-mode. This field is usually called the radion, r. Its mass can be obtained from the effective potential V (ϕ) and is given by where v v is the value of S at the visible brane and = m 2 /4k 2 (with m the mass of the field S). Quite generally 1 and, therefore, the mass of the radion can be much smaller than the first KK-graviton mass.
The radion, as for the KK-graviton case, interacts with both the DM and SM particles. It couples with matter through the trace of the energy-momentum tensor T [16]. Massless gauge fields do not contribute to the trace of the energy-momentum tensor, but effective couplings are generated from two different sources: quarks and W boson loops and the trace anomaly [50]. Thus the radion Lagrangian takes the following form [49,51]: where F µν , F a µν are the Maxwell and SU (3) c Yang-Mills tensors, respectively. The C 3 and C EM constants encode all information about the massless gauge boson contributions and are given by: with x q = 4m q /m r and x W = 4m w /m r . The explicit form of F 1/2 and the values of the one-loop β-function coefficients b are given by [50]: where n is the number of quarks whose mass is smaller than m r /2.
Possible couplings between KK-modes of the bulk scalar field S, the DM and SM fields are usually allowed, in the absence of some ad hoc bulk symmetry to forbid them. In the rest of the paper we will not include them, since we want to focus on just gravitational mediators (radion and KK-gravitons) between the SM and the dark particles.
Finally, we want to comment about the AdS/CFT correspondence, which suggests a duality between strongly coupled conformal field theories in 4D and weakly coupled gravity in 5D (see, for example, Ref. [52] and refs. therein), also called holography. Within this framework, the extra-dimensional model described above can be interpreted as a strongly interacting theory in which the particles localized at the IR-brane are bound states, while the presence of gravity mediators (KK-gravitons and radion) is a consequence of the conformal symmetry of the composite sector, spontaneously broken by the strong dynamics. The radion is thought to be the Goldstone boson of dilatation symmetry in 4D, i.e., the dilaton, although the dual interpretation of the massive gravitons is not so well understood [16]. The scale Λ in the holographic dual corresponds to the scale of conformal symmetry breaking in 4D.

Dark Matter Production in the Early Universe
In Refs. [23,29] some of us have studied how to reach the observed DM relic abundance in the freeze-out scenario. Freeze-out occurs if the interactions between DM and SM particles are strong enough to bring them into chemical equilibrium. However, if the interaction rates between the visible and the dark sectors were never strong enough, the observed DM relic abundance could still have been produced in the early Universe by non-thermal processes. This is what occurs in the so-called freeze-in mechanism.
The evolution of the DM, radion and KK-gravitons number densities (n, n r and n K respectively) is given by a system of coupled Boltzmann equations: where H corresponds to the Hubble expansion rate, and n eq i are the number densities at equilibrium of the species i. Interactions that only involve bulk particles, namely KKgravitons and radions, both in the initial and final states are subdominant due to a strong suppression of 1/Λ 8 . The quantity γ Φ →SM is the interaction rate density for the 2-to-2 annihilations of a field Φ (either DM, KK-graviton or radion) into SM particles. Similarly, γ d Φ →DM and γ d Φ →SM are the interaction rate densities for the 2-body decay of a field Φ into DM and SM particles, respectively. Let us notice that in this extra-dimensional picture we need a Boltzmann equation like eq. (3.3) for every KK-mode.
In the freeze-in paradigm DM never gets in thermal equilibrium with the rest of the SM particles of the primordial plasma. It is usually assumed that after inflation the abundance of DM was negligible, and slowly produced via interaction between the SM particles. Along the evolution of the Universe, the DM abundance was generated via two main processes: 1. Direct freeze-in. The DM abundance is generated directly by the annihilation of SM particles via an s-channel exchange of KK-gravitons or a radion.
2. Sequential freeze-in or freeze-in from the dark sector. The DM abundance is generated by decays of KK-gravitons or radions, previously produced by annihilations or inverse decays of SM particles via direct freeze-in. This scenario has been doubted "sequential freeze-in" [53].
Another production channel corresponds to the case in which the DM abundance is set entirely in the hidden sector by 4-to-2 interactions [54,55]. However, such a possibility is sub-dominant due to a strong suppression by higher orders of the scale Λ. It has been also shown that, independently of the nature of DM, it is possible to populate the relic abundance through a freeze-in mechanism via the exchange of a massless spin-2 graviton [1][2][3][4]. However, for this mechanism to be dominant, reheating temperatures T rh of the order of 10 13 GeV for a DM mass of 1 MeV are required. We will see in the following that, in this warped extra-dimensional setup (with KK-gravitons and the radion as additional fields playing the freeze-in mechanism) a much wider range of T rh is indeed possible. These two main mechanisms previously mentioned, i.e. the direct and the sequential freeze-in, will be described in detail in the following subsections.

Direct Freeze-in
As it was briefly sketched above, in the case of direct freeze-in the DM abundance n is generated by the annihilation of SM particles via an s-channel exchange of KK-gravitons or a radion. 2 A standard way to rewrite the Boltzmann equations is using the dimensionless yield Y ≡ n/s, with s the SM entropy density (not to be confused with the Mandelstam variable s). The SM entropy density is defined, as a function of the temperature, as s(T ) = 2π 2 45 g s (T ) T 3 (where g s (T ) is the effective number of relativistic degrees of freedom [58]). If the production cross-section is small enough to keep DM out of chemical equilibrium 2 Another possibility corresponds to the interactions mediated by Higgs bosons. However, we focus here on the extra-dimensional portal ignoring the Higgs one. This can be reached by assuming a quartic coupling λ hχ between the Higgs and the DM such as λ hχ 10 −10 [56,57]. with the SM bath, and the evolution of the DM abundance n (or of the yield Y ) is largely dominated by the interaction rate density γ DM→SM , eqs. (3.1) to (3.3) can be simplified to: In a Universe dominated by SM radiation the Hubble expansion rate is , where the SM energy density is ρ SM (T ) = π 2 30 g (T ) T 4 and g (T ) is the effective numbers of relativistic degrees of freedom for the SM radiation [58]. Then, eq. (3.4) becomes: where T rh is the reheating temperature which, in the approximation of a sudden decay of the inflaton, corresponds to the maximal temperature reached by the SM thermal bath. In order to get eq. (3.5) a vanishing initial DM abundance at T = T rh was assumed and the temperature dependence of g (T ) and g s (T ) has been neglected. The asymptotic values g and g s correspond to the SM values for T m t , g = g s = 106.75 (which take into account all SM degrees of freedom). Since this approximation is reliable for temperatures above the QCD phase transition, we explore the range T rh 1 GeV.
The interaction rate density γ DM→SM can be computed from the total DM annihilation cross-section into SM states σ DM→SM which, in the limit where the DM and SM particle masses are negligible, can be expressed as: 3 The details of the individual cross-sections are reported in Appendix A.1.
where the two terms correspond to the exchange of KK-gravitons and the radion, respectively. Left panel of Fig. 1 shows with a solid black line an example of the DM annihilation cross section σ DM→SM for a particular point in the parameter space, m r = 10 3 GeV, m 1 = 10 7 GeV and Λ = 10 9 GeV. Notice that this cross-section is largely independent of the DM mass, m χ , as long as m χ s. Then, also the interaction rate density becomes independent of m χ provided m χ T . Therefore, in the following we will consider as a benchmark point m χ = 1 MeV to illustrate our results, but keeping in mind that they can be extended to a wide range of DM masses, typically between the keV and PeV scale. The first peak at s = m 2 r corresponds to the resonant exchange of a radion, whereas the following well-separated peaks correspond to the lightest KK-graviton modes. The non-trivial behavior for s m 2 1 is due to the sum over poles and interferences of many different KK mediators. For very large values of the KK-number n, the widths of the KK-graviton resonances become comparable to their mass gap, Γ n ( √ s) ∆m. This happens approximately for: as at large n the KK-modes separation is a constant, ∆m m 1 /x 1 , see eq. (2.12). In this regime the resonances overlap and become individually indistinguishable. They eventually merge into one single contribution to the cross-section, as it can be seen in the rightmost region of Fig. 1 (left). Finally, the red-shaded region corresponding to s > Λ 2 is beyond our EFT approach, being the center-of-mass energy of the process larger than the effective scale of the theory. In order to solve eq. (3.5), we need to compute the interaction rate density γ DM→SM as a function of the temperature. In general, for the process where two particles (i, j) annihilate into two states (k, l), the interaction rate density i + j → k + l is defined as: where s min ≡ max (m i + m j ) 2 , (m k + m l ) 2 , σ R is the reduced cross-section summed over all the degrees of freedom of the initial and final states, and K 1 is the modified Bessel function. σ R can be written in terms of the cross-section σ(s) as: Several useful approximations can be implemented for different ranges of T , such that the interaction rate density γ DM→SM for the DM annihilation into SM states becomes: (3.10) The right panel of Fig. 1 shows the DM interaction rate density for m r = 10 3 GeV, m 1 = 10 7 GeV and Λ = 10 9 GeV with a black solid line, whose behavior as a function of the temperature can be easily understood using the approximations in eq. (3.10). At low temperatures (T m r /2) all the mediators are very heavy and decouple from the lowenergy theory; the rate presents a strong temperature dependence, γ ∝ T 12 , represented by a red-dotted straight line in the plot. When T m r /2, the resonant exchange of a radion dominates and γ ∝ T K 1 (m r /T ). This corresponds to the first bump in the plot, again coinciding with a red-dotted (curved) line. In the intermediate regime, m r /2 T m 1 /2, the temperature is higher than the radion mass but still smaller than all KK states. The interaction is, thus, driven by the exchange of the light radion, with γ ∝ T 8 . This is shown by the second straight red-dotted line in the plot, with a slope smaller than the first one (as it is proportional to T 8 , compared to T 12 in the first region). We reach then the region in which the KK-gravitons dominance takes over: first, at the peak of the first KK-graviton mode (T m 1 /2) for which γ ∝ T K 1 (m 1 /T ), corresponding to the second bump in the plot; eventually, when the increase of the temperature makes heavier KK-graviton states to have a sizable contributions to the rate, with a constructive interference giving a γ ∝ T 7 behavior. We can see that all the different regimes in T follow extremely well the curved and straight blue-and red-dotted lines, corresponding to the approximate behaviors depicted in eq. (3.10). As for the left panel, the red-shaded region corresponding to T > Λ is beyond our EFT approach.
Notice that a big hierarchy between m r and m 1 was chosen in order to avoid an overlap between the two bumps, such that the five regimes in eq. (3.10) can be clearly seen in the plot. For generic choices in the parameter space, overlap between regions may occur.
Using the approximated expressions of γ DM→SM from eq. (3.10), the Boltzmann equa-tion (3.5) can be analytically solved, finding for the different regions in T : 11) where Y 0 corresponds to the asymptotic value of Y (T ) for T T rh . The final DM yield in eq. (3.11) has a strong dependence on T rh , characteristic of the UV freeze-in production mechanism.
Finally, let us emphasize that for the previous analysis to be valid, the DM has to be out of chemical equilibrium with the SM bath. One needs to guarantee, therefore, that the interaction rate density be γ DM→SM n eq H, which translates into:  where W −1 [x] corresponds to the −1 branch of the Lambert W function computed at x. Fig. 2 shows the reheating temperature T rh required to reproduce the experimentally observed DM abundance, Ω χ h 2 , for a fixed value of the DM mass, m χ = 1 MeV. In the left panel, we show T rh as a function of the first KK-graviton mass, m 1 , for fixed Λ = 10 11 GeV; in the right panel, we show T rh as a function of Λ for fixed m 1 = 10 5 GeV. The radion mass has been chosen as m r = m 1 /10 3 (therefore, it is a variable parameter in the left panel, whereas it is a fixed one in the right panel). In order to compute T rh , the DM yield has been held fixed so that m χ Y 0 = Ω χ h 2 1 s 0 ρc h 2 4.3 × 10 −10 GeV, where ρ c 1.1 × 10 −5 h 2 GeV/cm 3 is the critical energy density, s 0 2.9 × 10 3 cm −3 is the entropy density at present and Ω χ h 2 0.12 [59]. The gray-shaded areas are the regions where chemical equilibrium with the SM is reached and, therefore, where the freeze-in cannot occur and the analysis performed here is not valid. The black-dotted lines, representing T rh = m 1 /2 and T rh = m r /2, have been added for reference. Eventually, the red-shaded areas (m 1 > Λ) represent the regions for which the EFT approach breaks down. For the sake of completeness, notice that the s-channel exchange of a (massless) graviton gives an irreducible contribution to the total DM relic abundance [1][2][3][4]. The corresponding interaction rate density is given by: and, therefore, its contribution to the DM yield is: (3.14) We stress that this expression is a function of T rh , only, being naturally independent from Λ and the masses of the KK-gravitons and the radion. This contribution, indeed, comes from 4-dimensional gravitons or, in the case considered here, from the long distance (low-energy) limit of 5-dimensional gravitons (corresponding to the KK-graviton zero-mode). Due to the large hierarchy Λ M P , the contribution of the massless graviton is typically subdominant and can be disregarded. For example, for a DM mass m χ = 1 MeV it would only be relevant for reheating temperatures T rh ≥ 10 13 GeV, i.e. well above the range of T rh depicted in Fig. 2.

Sequential Freeze-in
In this case the DM abundance comes from decays of KK-gravitons or radions, previously produced via the freeze-in mechanism. Such states are mainly generated by 2-to-2 annihilations or inverse decays (2-to-1) of SM particles. We will now review one by one the two possibilities.

Via Annihilations
KK-gravitons and radions with masses below the reheating temperature can be created onshell in the early Universe via annihilations of two SM particles by the freeze-in mechanism. Once created, their decay products may contribute to the total DM relic abundance. In this scenario, eqs. (3.1) to (3.3) can be simplified as: where the rates are: Notice that the m −8 n factor in γ KK → SM (and, hence, the strong temperature dependence) comes from the polarization tensor of the KK-gravitons (as it was shown in Refs. [16,23] for spin-2 massive particles). Such a suppression is not present in the case of radions (that have spin 0). The branching ratios into DM particles are: Notice that in eq. (3.20) only the lightest KK-graviton is taken into account. This is a consequence of the strong suppression with the KK-graviton mass m n in eq. (3.16). Even if all of the KK-gravitons do contribute to the total DM density, the only relevant contribution is given by the lightest state. For the previous analysis to be valid, the KK-gravitons and the radion must be out of chemical equilibrium with the SM bath, which corresponds to the conditions γ KK→SM n eq K H and γ r→SM n eq r H. The reheating temperature in this limit satisfies the tightest of the following conditions (depending on the mass of the lightest KK-graviton, m 1 ):   For T rh < m r , on-shell KK gravitons and radions are not produced in the early Universe, and therefore this mechanism can not account for the DM relic abundance. If m r < T rh < m 1 , only radions are created. In this region T rh is independent on m 1 (and therefore on m r ) due to the fact that the interaction rate in eq. (3.17) does not depend on m r , as it can be seen in the left panel of Fig. 3. Now, if T rh > m 1 the KK-gravitons are also produced and their decay dominate the DM production. The reheating temperature needed to reproduce the observed value of Ω χ h 2 in this region is very near to the border of the gray-shaded area for which the DM is in equilibrium with SM particles and freeze-in does not occurs (remember, though, the log-log scale of the plots).

Via Inverse Decays
Alternatively, frozen-in KK-gravitons and radions are also created on-shell via inverse decays of SM particles (a 2-to-1 process), and subsequently they can decay into DM particles. Similarly to the analysis done in the case of freeze-in via annihilation, in this scenario eqs. (3.1) to (3.3) can be simplified to: where the interaction rate densities for decays are defined by: with Γ the decay width obtained by summing (rather than averaging) over the degrees of freedom of the decaying particle. Using eqs. (A. 19) and (A.21) we get, then: Eq. (3.22) admits the following approximate analytical solution: (3.26) In this case, most of the DM production happens at T m n /2.5 and T m r /2.5 for KKgravitons and radions, respectively. However, the sum over KK-modes should be performed up to KK-graviton states with mass below the reheating temperature, m n < T rh . For this reason, the total contribution due to the decay of KK-gravitons explicitly depends on T rh (whereas the second term in eq.  The gray-and red-shaded areas are the regions where chemical equilibrium with the SM is reached, 4 and where the EFT approach breaks down (as m 1 > Λ), respectively.
As in the case of sequential freeze-in via annihilation, for T rh < m r on-shell KKgravitons and radions are not produced in the early Universe and, therefore, this mechanism can not account for the DM relic abundance below the T rh = m r black-dotted line. In the region m r < T rh < m 1 , only radions are created and, in this case, the DM yield is independent on T rh (as the second term in eq. (3.27) does not depend on T rh ). This can be clearly seen in Fig. 4. For T rh > m 1 , the KK-graviton states are also produced. Their decay eventually dominate the DM production and the reheating temperature is proportional to √ m 1 (left panel) or Λ (right panel).
So far, each individual production channel has been studied separately. Fig. 5 depicts the parameter space favored by the observed DM abundance for m χ = 1 MeV and Λ = 10 14 GeV as a function of m 1 (upper left panel), or m 1 = 10 3 GeV as a function of Λ (upper right panel), taking into account all of the three DM production mechanisms described previously (i.e. direct production, sequential production via annihilation and sequential production via inverse decay). The thin blue lines correspond to the partial contributions of each of the mechanisms, whereas the black thick line to the total abundance. Eventually, in the lower panel we show the correlation between Λ and m 1 at a fixed value of the reheating temperature required to achieve the observed DM abundance, T rh = 10 5 GeV. In all panels, the radion mass is related to the first KK-graviton mass as m r = m 1 /10 3 . As always, the gray-and red-shaded areas represent the regions where chemical equilibrium between DM and the SM particles is reached and where the EFT breaks down since m 1 > Λ, respectively. Finally, in Fig. 6 we show the correlation between Λ and m 1 required to reproduce the observed DM abundance for m χ = 1 MeV and m r = m 1 /10 3 for several representative values of the reheating temperature, T rh = 1, 10, 10 5 , 10 8 and 10 10 GeV (notice that the range of Λ plotted in Fig. 6 differs from that in the lower panel of Fig. 5). Let us note that the lines corresponding to T rh = 1 and 10 GeV overlap when m 1 10 3 GeV and Λ 10 9 GeV. This can be understood seeing that in that region, the DM relic abundance is mainly generated by sequential freeze-in via inverse decays of the radion and is therefore independent of T rh , see eq. (3.26). The red-shaded area, as always, represents the region where the EFT approach breaks down. On the other hand, the upper left green corner corresponds to radion lifetimes higher than 1 s, potentially problematic for BBN (all the KK-graviton states are heavier than the radion and therefore will naturally have shorter lifetimes). Eventually, the blue-shaded regions depict present and future experimental bound coming from resonance searches at the LHC. The proton-proton collision can generate resonant KK-gravitons that later decay into SM particles. ATLAS and CMS put bounds over these processes in γγ and lepton-lepton channels as a function of the mass of the resonance (the lightest KKgraviton). These bounds can be translated into limits over Λ as a function of the mass of the first graviton m 1 . The present bounds (dark blue) come from the resonant searches at LHC with 36 fb −1 [60] and [61], whereas future bounds are estimated assuming 300 fb −1 (medium blue) and 3000 fb −1 (light blue) for the LHC Run-III and High-Luminosity LHC, respectively. Notice that in this plot we do not show the gray-shaded region for which DM is in equilibrium with SM particles (where freeze-in does not occur), as we should draw a different region for each value of T rh .

Beyond the Sudden Decay Approximation of the Inflaton
While reheating is commonly approximated as an instantaneous event, the decay of the inflaton into SM radiation is a continuous process [62]. During reheating in which the SM thermal bath is produced and the Universe makes a transition to SM radiation domination, the bath temperature may rise to a value T max which exceeds T rh [63]. Away from the instantaneous approximation for reheating, the maximum temperature of the thermal bath may reach T max T rh prior to cooling. It is plausible that the DM relic density may be established during this reheating period, in which case its abundance will significantly differ from freeze-in calculations assuming radiation domination. In particular, it has been observed that if the DM is produced during the transition from matter to radiation dom-ination via an interaction rate that scales like γ(T ) ∝ T n , for n > 12 the DM abundance is enhanced by a boost factor proportional to (T max /T rh ) n−12 [64], whereas for n ≤ 12 the difference between the standard UV freeze-in calculation differ only by an O(1) factor from calculations taking into account non-instantaneous reheating. More recently, it has been highlighted that the critical mass dimension of the operator at which the instantaneous decay approximation breaks down depend on the equation of state ω, or equivalently, to the shape of the inflationary potential at the reheating epoch [65,66]. Therefore, the exponent of the boost factor becomes (T max /T rh ) n−nc with n c ≡ 6 + 2 3−ω 1+ω , showing a strong dependence on the equation of state [65]. Subsequent papers have explored the impact of this boost factor in specific models [4,[67][68][69][70][71][72][73][74][75][76][77]. Finally, another way for enhancing the DM abundance occurs in cosmologies where inflation is followed by an epoch dominated by a fluid stiffer than radiation. In such scenarios, even a small radiation abundance, produced for instance by instantaneous preheating effects, will eventually dominate the total energy density of the Universe without the need for a complete inflaton decay. In particular, a strong enhancement if DM production happens via interaction rates with temperature dependence higher that n c = 6 [78].
The present model of KK FIMP DM in warped extra-dimensions features processes where the interaction rate has a particularly strong temperature dependence, the most relevant ones being: i) the DM annihilation into SM states for reheating temperatures much lower than the radion mass γ DM → SM (T ) ∝ T 12 ; ii) the same process near the resonances T rh m r /2 and T rh m 1 /2, where γ DM → SM (T ) ∝ T K 1 (m i /T ) (with m i being the radion or the lightest KK-graviton mass, respectively); and iii) the KK-graviton annihilation into SM particles γ KK → SM (T ) ∝ T 16 . In these regimes, the non-instantaneous decay of the inflaton is expected to generate a strong boost factor to the DM yield, which translates into a reduction of the reheating temperature required to match the observed DM relic abundance. As the precise determination of such boost factors depends on the details of the inflationary model (in particular on the energy density carried by the inflaton and its equation of state parameter previous to its decay), it is beyond the scope of this study.

Conclusions
Dark Matter (DM) is typically assumed to be made of weakly interacting massive particles (WIMPs), produced in the early Universe via the freeze-out mechanism. Freeze-out occurs if the interactions between DM and SM particles are strong enough to bring them into chemical equilibrium. However, if these interaction rates were never strong enough, the observed DM relic abundance could still have been produced by non-thermal processes, like the freeze-in mechanism. In that case, DM is called a feebly interacting massive particle (FIMP).
In a warped extra-dimensional scenario, DM could naturally be a FIMP, if the effective gravitational scale Λ is much higher than the electroweak scale. In this case, DM is produced in two main ways: i) promptly by annihilations of SM states via the s-channel exchange of KK-gravitons and radions, i.e. the so called direct freeze-in, and ii) by decays of KK-gravitons or radions, previously produced by annihilations or inverse decays of SM particles via direct freeze-in. This scenario has been doubted sequential freeze-in.
In this paper we have systematically studied the different regions of the parameter space that generate the observed DM abundance in the early Universe, within a warped extra-dimensional model . We assume that both the SM and the DM particles are localized in the IR-brane, where the effective four-dimensional Planck scale is given by Λ, which is allowed to vary in a wide phenomenological range, [10 2 , 10 16 ] GeV, relaxing the requirement for the RS model to solve the hierarchy problem. We also include the radion, using the Goldberger-Wise mechanism [37] to generate the required potential to stabilize the size of the extra dimension. For definiteness, we consider scalar DM and focus on its interactions with gravitational mediators, i.e., the radion, the graviton and the KK-gravitons.
As the interaction rates between the visible and the dark sectors have a strong temperature dependence, the bulk of the DM density is typically produced at the highest temperatures reached by the SM thermal bath, which in the approximation of a sudden decay of the inflaton corresponds to the reheating temperature, T rh . This is a characteristic of the so-called UV freeze-in. We found however a case where the DM abundance was mainly produced at much lower temperatures, corresponding to the sequential freeze-in where the radion was generated via inverse decays. In that case the peak of the production happens when the temperature approaches the radion mass, T m r /2.5.
The possibility of generating the DM relic density within the RS scenario via the usual freeze-out mechanism was analyzed in Refs. [16][17][18][19][20][21][22]24]. After including the DM annihilation channel into KK-gravitons previously disregarded, it was found that even when both SM and DM particles live in the IR-brane there is a region compatible with the experimental and theoretical constraints where it is possible to reach the correct DM relic abundance [23]. The allowed region corresponds to m χ ∈ [1,15] TeV and Λ ∈ [10, 10 4 ] TeV. The upper limit on the DM mass comes from unitarity, while the lower limit is an indirect one, derived from searches at LHC of KK-graviton resonant production, which constrains the scale Λ as a function of the first KK-graviton mass. This bound is very relevant, since it determines the minimum value of the DM mass for which the annihilation channel into the first KK-graviton mode is kinematically open, leading to the observed DM relic density. In the freeze-out scenario, the LHC prospects for the near future exclude most part of the allowed region.
In the present work we find that it is also possible to obtain the correct DM relic abundance in the same RS model via the freeze-in mechanism, for DM masses in a much wider range spanning typically from the keV to the PeV scale, and larger values of the scale Λ than in the freeze-out scenario. This implies that the LHC bounds on the parameter space of the model are weaker than in the freeze-out case. This can be seen in Fig. 6, where we summarize our results in the (m 1 , Λ) plane for the benchmark DM mass m χ = 1 MeV, finding that only the lower-left corner will be probed by HL-LHC. On the other hand, other constraints are relevant, such as the life-time of the radion, which we require to be larger than 1 s to avoid problems with BBN, and excludes the upper-left corner. The results are not strongly dependent of the radion mass: for this reason we fix m r = m 1 /10 3 , in agreement with the expectation within the Goldberger-Wise mechanism. We find that the observed DM relic density can be obtained in a wide range of reheating temperatures, T rh ∈ [1, 10 10 ] GeV. Notice that we find some region of the parameter space for which the observed DM relic abundance is achieved with Λ as low as a few TeV (with lower values excluded by LHC data). In this region, the hierarchy problem is mostly solved, leaving only a remnant little hierarchy to be explained.
Finally, we argued that a more detailed analysis of the present model will require to go beyond the usual approximation where the inflaton decays instantaneously, and therefore the reheat temperature is the maximal temperature reached by the SM thermal bath. This is due to the strong temperature dependence of some interaction rate densities that enter in the determination of the DM relic abundance. A complete analysis must take into account the details of the inflationary model (in particular on the energy density carried by the inflaton and its equation of state parameter previous to its decay), and is therefore beyond the scope of this study. where P ≡ (s − m 2 r ) 2 + m 2 r Γ 2 r −1 is the radion propagator. For the SM massless vectors the vertex is generated by the trace anomaly and, therefore, the cross sections are proportional to α EM and C EM for the photon case, and to α 3 and C 3 for the gluon case, as given in eqs. (2.19) and (2.20).

A.5 Radion Decays
Eventually, the decay widths of radions into SM and DM particles are: where again all SM masses were neglected for simplicity.