Disformal Dark Matter

We generalize dark matter production to a two-metric framework whereby the physical metric, which couples to the Standard Model (SM), is conformally and/or disformally related to the metric governing the gravitational dynamics. We show that this setup is naturally present in many Ultra Violet (UV) constructions, from K\"ahler moduli fields to tensor-portal models, and from emergent gravity to supergravity models. In this setting we study dark matter production in the early Universe resulting from both scatterings off the thermal bath and the radiative decay of the inflaton. We also take into account non-instantaneous reheating effects at the end of inflation. In this context, dark matter emerges from the production of the scalar field mediating the conformal/disformal interactions with the SM, i.e. realising a Feebly Interacting Matter Particle (FIMP) scenario where the suppression scale of the interaction between the scalar and the SM can be taken almost as high as the Planck scale in the deep UV.


I. INTRODUCTION
The presence of two geometries, governing the gravitational dynamics and the behaviour of matter fields respectively, is frequent within the landscape of highenergy physics models. For instance, orbifolds in string theories [1] or Kähler metrics in supergravity models [2] are two popular cases where the geometry governing the dynamics of matter are not the same that the one governing the gravitational structure of space-time. This is in fact quite an old idea and was already proposed in Nordstrom gravitational theories [3], Brans-Dicke's [4] or Dirac's [5]. More recently models of emergent gravity [6] modify the metric assuming that gravity springs from vector interactions generated in massive hidden sectors. As a consequence, the dynamical metric in Minkowski space-time can be reduced to [7] where T µν is a linear combination of energy-momentum tensors of hidden sector particles and Λ the scale beyond which the theory breaks down 1 [the Beyond the Standard Model (BSM) scale]. Phenomenological consequences of emergent gravity in dark matter phenomenology and for the early Universe evolution have been studied in [8]. Other possibilities are common in the supergravity framework and involve the presence of moduli 1 Typically the mass scale of the hidden sector.
fields, inducing new couplings to the Standard Model of the type [9] L SM T ⊃ Z H |D µ H| 2 (2) in the Higgs sector with Z H = 1 + 1 Λ t, t being the real part of moduli fields. This setup can also be considered as a modification of the physical geometry, dependent on the moduli fields, especially their stabilized values at the minimum of the Kähler potential. Finally one finds similar examples in high-scale SUSY models [10]. Indeed, the minimal coupling of a gravitino 2 to the Standard Model, whose longitudinal mode is the Goldstino denoted by Ψ 3/2 , is built by first defining a vierbein [11] e α µ = δ α µ − F being related to the SUSY breaking scale 3 . This enters clearly in the category of models where the Standard Model fields interact with the gravitino (Goldstino) through its presence in the physical 4 metric g µν = e µ a e ν b η ab .
As we have just seen, several constructions include two geometries, and the relation between the gravitational metric and the physical one implies modifications of the dynamics and the phenomenology of SM fields by the introduction of new couplings, new interactions and new fields. Some time ago, it was proposed [12] to generalize this in a unique description where all types of dynamical metrics can respect some basic principles. The generalization consists in considering two metrics which are not just conformally related. This more natural relationship involves Finslerian geometry rather than Riemannian geometry. Finslerian geometry is the most general geometry where the squared relativistic interval ds 2 is homogeneous of second degree in the element dx, in other words, ds 2 (x, µ dx) = µ 2 ds 2 (x, dx) .
The expression (8) contains not only the classical conformal transformation induced by C between the two metrics, but also the possibility for a disformal transformation through the coefficient D, disformal in the sense that the space-time structure is stretched differently in each direction proportionally to ∂ i φ in the ith-direction. As expected, if φ is a constant, i.e a homogeneous and isotropic field, both metrics are related by a simple conformal transformation. If not (D = 0), φ is interacting with the matter fields through their kinetic terms. Notice that the metricg µν can also be inferred by requiring general covariance and the absence of derivatives of order larger than two. The latter requirement follows from the generic appearance of ghosts in theories with higher order derivatives. Disformally related metrics have been widely used in the cosmological, gravitational and recently particle physics contexts [14][15][16][17][18][19]. For instance, the authors of [20] have recently given to φ the role of the quintessence field of dark energy and analysed the parameter space defined by (C, D) which is cosmologically allowed. Similarly in the gravitational context, the authors of [21] have constructed disformal versions of the Kerr space-time. We will use disformally related metrics to induce dark matter production. In this setting, we will unravel how the phenomenology of the early Universe and the appearance of dark matter could be modified by the introduction of disformal coefficients. We will also restrict the corresponding parameter space from late-time observables. In particular, we will show that a FIMP scenario for dark matter production can be naturally realised with a disformally coupled scalar field φ to the standard model.
Recently, [22] proposed that φ could play the role of a portal between a WIMP dark sector and the visible sector. However, the WIMP paradigm is nowadays under high scrutiny due to the lack of observed signal, especially in direct detection experiments. In [22] the direct detection constraints were not taken into account, but could drastically reduce the allowed parameter space, especially for dark matter masses below 100 GeV. Indeed, the more recent measurements exclude proton-WIMP cross section σ 10 −46 cm 2 for a 100 GeV dark matter mass [23][24][25], which is more than six orders of magnitude below the cross section for the vanilla models of weakly interacting particles [26,27]. The simplest cosmologically viable extensions of the Standard Model reproducing the relic abundance observed by Planck experiment [28] require to invoke a new physics scale 5 TeV [29] which will be probed in the next generation of experiments [30]. If no signal is seen, this BSM scale will be pushed even further well above 50 TeV.
However, relaxing the requirement of thermal equilibrium between the dark sector and the primordial plasma opens a completely new field of research [33,34]. A Feebly Interacting Massive Particle (or Freeze-In Massive Particle, FIMP) couples too weakly with the Standard Model bath to reach thermal equilibrium in the early Universe. Such seclusion appears naturally in models where the mediator is very heavy, e.g. Z of unified theories [35,36], massive spin-2 particles [37], moduli [9], inflaton-like portals [38] or in the Kaluza-Klein theory framework [39]. Another possibility is to consider theories where the couplings are reduced by a mass parameter of the order of the Planck mass scale M P , as in supergravity 5 [2] or a combination of the supersymmetry breaking scale and the Planck mass in High Scale SUSY scenarios [10]. In all these cases, the temperature dependence of the production rate renders the physics in the earliest stages of the Universe more complex than the vanilla reheating scenarios described in [41]. Noninstantaneous thermalization [42] or non-instantaneous reheating [43,44] modify drastically the distribution function and/or the production rate of particles in the Standard Model plasma, making the dark matter density calculation more complex. Considering quantum effects to the inflaton decay [45] or the possibility of nonstandard inflaton potentials [46,47], show that the study of physics at the end of the coherent oscillation stage at the end of inflation should be treated with care.
In, this work, we propose to consider scenarios where the dark matter is composed of the field φ which defines the physical metric. Indeed, a quick look at Eq. (8) shows that the disformal term is decreasing in magnitude with the BSM scale Λ, above which the dynamical version of the metric breaks down anyway. This scale suppression should seclude φ sufficiently from the Standard Model plasma to suppress the production of φ and making it a perfect FIMP candidate. Moreover, the form of the metric, dictated by the consistency conditions (conservation of the signature and causality) implies a discrete Z 2 symmetry which ensures the stability of φ.
This paper is organized as follows. After a description of our models and the expression of the couplings generated by a disformal metric in section II, we compute the dark matter abundance in section III and its phenomenological consequences before concluding. Throughout this work, we use a natural system of units in which k B = = c = 1. All quantities with dimension of energy are expressed in GeV when units are not specified.

II. THE SETUP The Lagrangian
Supposing that the action is divided into a part driven by the geometrical (gravitational) metric g, whereas matter follows the geodesics of a physical (dynamical) metric g, we can write the matter Lagrangian L m (g, h) being the SM Lagrangian expressed in term of the metricg. For instance, considering one real scalar SM degree of freedom (denoted by h) for simplicity, gives In the literature, the scalar φ is often directly or indirectly related to dark energy, or represents the quintessences field, and a shift symmetry φ → φ + c is then imposed to avoid dangerous mass terms. In our case, as we want to be as generic as possible, we do not impose this symmetry. Concerning the matter fields, we restrict ourselves to a singlet-like scalar h field, representing the Higgs boson, to simplify the equations and explanations. Of course, the complete particle content of the Standard Model is considered for our numerical results. For the same reason, the partial derivatives in Eq. (11) should be understood as covariant derivatives. However, as discussed further on, the covariant part of the derivatives do not contribute significantly to the DM production and therefore are omitted for the sake of simplicity. By expanding the physical metricg µν in terms of the geometrical metric g µν and a small deviation δg µν g µν , justified by the fact that we consider processes occurring at energies much below the BSM scale Λ, the matter action Sg m can be expressed in the Einstein frame as (12) at lowest order in δg µν /g µν 1.
With T m µν the energymomentum tensor of matter fields, from Eq. (8) we have giving where T m µν can be expressed as where the sum is performed over all SM particles of spin i, whose corresponding energy-momentum tensors are given by for scalar (h), fermionic (f ) and vectorial (A µ ) matter fields respectively. F µν = ∂ µ A ν − ∂ ν A µ is the field strength of the spin-1 field whereas V (h) represents the scalar potential. Masses of various SM states are discarded as the typical temperatures involved in early Universe processes are much above the electroweak scale. Terms of the form g µν V (h) can be discarded in the scalar energy-momentum tensor. Indeed, as discussed below, such terms correspond to processes involving a higher number of SM particles and/or suppressed by additional SM couplings, compared to processes relevant for the DM production. Moreover, the term ig µνf / ∂f in the fermionic energy-momentum tensor vanishes for on-shell states. In addition, the trace of the energy-momentum tensor (T i ) µ µ vanishes for i = 1/2 and i = 1 but not for i = 0. This is due to the fact that the energy-momentum tensor for massless states acquires a conformal symmetry in four dimensions for fermions and vectors but only in 2 dimensions for scalars. As an example, the interaction term between our DM candidate φ and one SM real scalar degree of freedom h is given by Notice that we have not yet made explicit the functions C(φ, X) and D(φ, X). The only assumption made was that these functions yield a small δg µν /g µν 1.
The literature is replete with clever propositions, ranging from invoking shift symmetries as in [22] with a quintessence point of view where C and D depend only on X [32], to supposing constant X = ∂ µ φ∂ µ φ in studies of Kerr Black Holes [21]. Other popular examples are Horndeski theories which transform into themselves under special disformal transformations of the metric when C and D depends only on φ and not on X [40]. In this context, the resulting theories form almost the most general class of ghost-free scalar-tensor field theories. C and D can also be considered as dependent on φ only with expression of the type given by In our case, we propose to expand C and D around |φ| 2 (to ensure their positivity) which means Considering processes at energies much below M P , it is reasonable to stop the expansion to the first term as a first approximation. Following (1), notice that such terms emerge from the coupling to a scalar of mass m with where we identify Λ = M P , d = 1,c = 0, c X = −1/2 and c = m √ 2Mp . In the following we will leave these parameters free in a phenomenological way.

Disformal production process
Now that the Lagrangian is defined, one can investigate the DM production processes through scattering off Standard Model particles. As we commented in the previous section, we will focus our analysis on the case of a Standard Model bath composed of a real scalars h, whereas the numerical calculations will be done with the complete set of SM particles. The corresponding Feynman diagram is shown in Fig. 1 The production rate as-< l a t e x i t s h a 1 _ b a s e 6 4 = " s U 8 0 9 g g r S U Y b g j r w u u 8 e S z X n f W 8 sociated to this process 1 + 2 → 3 + 4 where 1, 2 denote particles of the SM and 3, 4 the φ states as represented in Fig. 1 at a temperature T can be written (see appendix for details and notations, as well as [10]) where E i denote the energy of particle i = 1, 2, 3, 4 and represents the (thermal) distribution of the incoming particles 6 . Using the Lagrangian of Eq. (17), the scattering amplitude M can be written where s and t are the Mandelstam variables. We then obtain for the production rate from the scalar scatterers Including the complete spectrum of the thermalized Standard Model species, i.e. production rates from fermions R 1/2 and vectors R 1 , one obtains the total rate as with β d 4 d 2 . The exact expression for the rate is given in Eq. (67) and more details regarding the derivation can be found in the Appendices.
Once we know the production rate R(T ), the relic abundance computation is relatively straightforward. One needs to solve the integrated Boltzmann equation where R(t) denotes the production rate of dark matter (per unit volume per unit time), or in term of temperature supposing an instantaneous thermalization, with T 2 M P , g T being the effective number of relativistic degrees of freedom at the temperature T . Solving the Boltzmann equation when with β a given constant gives for where Y scat φ is constant for T T RH . We assumed a vanishing dark matter density prior to reheating. We have defined the reheating temperature by the condition ρ Φ (T RH ) = ρ R (T RH ), (Φ being the inflaton field) in other words, when radiation and inflaton densities equilibrate. Notice that different definitions of the reheating temperature can lead to slightly different results, but differing never more than by factors of the order of unity as is shown for instance in [46].
The dark matter number density reaches its maximum almost immediately after the reheating process for a temperature of ( 3 10 ) 1/7 T RH and decreases at lower temperature with a constant n φ /T 3 ≡ Y scat . The present relic abundance, at T = T 0 , is given by where ρ 0 c /h 2 = 1.05 × 10 −5 GeV cm −3 is the present critical density and g i is the effective number of degrees of freedom at temperature 7 T i . From Eq. (29) we can compute the relic abundance of the φ field produced by scattering processes implies to focus on heavy dark matter candidates due to its very feeble production in the early stage of the reheating process. We can also extract an upper bound on Λ from the condition m φ T RH for the production to be kinematically allowed. We then obtain d being by definition of the order of unity 8 . This condition reflects the difficulty of producing φ in the earliest stage of the Universe. Planck mass couplings for instance would not be sufficient to produce dark matter with the right abundance, the majority of the reheating models predicting T RH 10 12 GeV [46].
Production from inflaton decay  It was shown in [45] that if dark matter is produced by scattering one cannot avoid the concomitant direct production of dark matter through the loop-induced inflaton decay, as depicted in Fig. 2. The minimal way to couple the Standard Model sector to the inflaton field Φ and to realize the reheating process is via the Higgs SU (2) L doublet H:

w 3 y g y q 3 P J 8 P z P O i c K u Z 4 7 U O x + / H k Y M S i h S y X i O I s 6 4 d B K s 5 z n U g g X W Q k x d E M T 0 i u F E J h w t I I j R m H d y K Q W n U 4 P M + y 1 X w I 5 B y L a d a U y U W v L A M D p g Q K + v B F z H G E 1 P f d Q 2 Q x / S m h S w I G Z m O t R t r g a J E u n 6 d K w t E P 8 g t O 0 C j G E w 5 7 2 2 Q / Y v M U i 6 f g 6 k 5 w r m L T r p u K / 3 H w F h n D m + Y D V W y Y F r d l 2 p B L U A b k K B + 8 k V Y X / f A 8 z w e M Y R p P k n w w H c f 0 o f p 3 E P G I H + y H B 0 U B U Z J / p R t B b 5 z E s i F m o F A u D c S E i S l 6 i u R / Z 0 g m M h k N o H d m g q S B g T I L C o N N 1 H A b F R q K b q O 6 l a 6 t y k 5 q b C v X M 9 z 3 2 6 j T 0 k V X e I E / K Q M d a I 4 j z r I r P g t P v y Q w b y C a J E M 4 W Z V q 0 I o t E I 1 n p H N l a 0 C j q h A 3 w O j 5 q 1 9 + R E f R C i d t O C p U H 2 6 v C K X s s j 2 k k J 4 P V K L B P F X P 5 p I r 8 m i + m h U 5 D A b 5 V 7 g E p N l U 2 B g U u B d M M f i V u m T X j 2 Z s 3 F A Z b i D 1 E S r w O P R A 5 f l q R l G + Y 6 y D 6 j D e g / j Y 7 j o s j + J D Q w + r A u B g x 2 E D G W E + g 4 a w c i k Y F 6 c N 0 A c 1 G X V m H 9 h t k p f T W J A G 1 U A 4 9 E 8 r C u s e Z h w n E 1 s J u H Z D H N Z R 6 S 1 P K K B m O S X C C b 8 H G 3 J 2 6 e Y n O M y z e Y T n h L r Z 2 T k 5 8 2 8 d y Y u m Y 6 K E 1 7 A J J 8 S p o A 1 H d j h d Q h t z r w G X t b E p + x q w r i R v A s 5 Z l j X w R z 1 g e x 6 r y 6 T 2 m L 2 O Z r P V k C 7 q e H a k z v C R J 4 k S I V z y i Q x o r + f L N 7 b C V a x g k h 5 r p z U 1 4 m V j 8 y d l 9 H 3 e l d V l w 6 G 3 n a h 4 u T p P N + x O 8 R J u X 4 Q 3 c q r r Y 5 c x s 1 N Z Z 6 g v q u m C p 9 T Z 3 Z / 2 a Z z M 3 H r b 1 J W b p N d A m U x N M P S S l z E d Z Y L j y y H h 3 O 4 u n c f w D p S l Z r z q I Q X t o x p T 6 o o X s Q R 8 F 5 M k I k g N N H U B U i 1 S 3 z B 1 v y z 3 H q g W Z S R l v 9 I S U G 3 W Z Z M q j K Y F q R U t K t r k g J F V O a H l K u O N s E x / L a m q x 0 j r c t K E n G l G q O a b e b K c Y t V z e q h p a V l p R q i H i T m F 7 G z V M V S b 0 5 K y h x q R 7 q h a J o v C S K o C q c 6 H V / X x Y D Z Y f r F l 1 S T S 8 m q Y G K C e L h / t H 9 0 o j N z 0 D X s H W 1 4 X o k 3 Y T q 4 v J 6 6 V l h F W 9 7 U Z 1 3 3 2 V N J U f Z 0 x k H X B c R h n O 7 u B O 5 S s L B d S t W Z s 0 o 2 k M s g 0 F u O 5 c h B W n t N z U U v L 1 m K E u t G Y
where h i , with i = 1, 2, 3, 4, denote the four real scalar degrees of freedom of the Higgs doublet above the electroweak symmetry breaking scale. The decay width of the inflaton into these fundamental scalars is given by The loop-induced decay width of the inflaton to a DM pair, whose corresponding diagram is depicted in Fig. 2, can be expressed as which gives for the number density of φ where we used ρ Φ (T RH ) = ρ R (T RH ) and the branching ratio is given by where N φ is the number of φ particles produced per decay, which is N φ = 2 in the present case. The DM relic abundance produced from inflaton decays is thus given by (30) It is remarkable that whilst at tree level, one needs to fine tune tiny dark matter couplings to the inflaton sector to ensure a branching ratio B R 10 −9 to avoid overproduction of dark matter, when one considers radiative production, for a BSM scale Λ of the order of 10 14 GeV, the disformal coupling d can easily reach unity without overclosing the Universe.
Moreover, comparing Eqs. (32) and (39), we see that the production has the same order of suppression in Λ, although for a reheating temperature below 10 12 GeV, it is clear that the radiative decay dominates over the scattering processes. To be more precise, we can ask ourselves for which value of T RH the scattering rate will begin to produce more dark matter than the radiative decay. We obtain T RH T eq RH = 2 × 10 12 GeV.
It is remarkable that this temperature does not depend either on m φ or Λ.

Conformal production
It is relevant to compare the disformal production to the one generated by the conformal coupling of Eq. (19). It is easy to understand that the part proportional to X = ∂ µ φ∂ µ φ will not be very different from the disformal part we just discussed. We computed the production rate in appendix, Eq. (69), and obtained a value of R(T ) = β c X T 12 /M 8 P , with β c X 10 c 2 X , i.e. with a numerical prefactor of the same order of magnitude as for β d 4 d 2 . Considering the coupling c X should then give similar phenomenological results as for the coupling d, when c X ∼ d(M P /Λ) 4 . However, the presence of a constant c term in Eq. (19) can affect drastically the dark matter production. The rate will then be given by where β c 1.1 × 10 −2 c 4 . The exact expression is given in Eq. (70). This is computed the in the same manner as β d , i.e taking into account all the Standard Model spectrum in the initial state. From the production rate R c we can deduce the relic abundance after integration on T : for the scattering processes, and for the decaying process, where we used More details regarding the calculations can be found in the appendix. We see then that for lower reheating temperature, T RH 10 11 GeV, the conformal couplings dominate the dark matter production from scattering over the disformal source. That is understandable because the dependence on the production rate is lower for conformal coupling than disformal coupling. The same can be said concerning the decay channel Φ → φφ which dominates for the disformal coupling. The possibility of having both conformal and disformal coupling at the same time will be discussed below.

Instantaneous reheating case
We show in Fig. 3 the parameter space yielding the correct relic abundance for d = 1 and c = 0 in particular the dashed blue curve corresponding to observed dark matter abundance. We recognize clearly the two regimes (scattering and decay) from their different dependence on the reheating temperature, especially the change of regime for T RH = T eq RH 2 × 10 12 GeV, as expected by our approximation (40). While the scattering process gives a mild dependence m φ ∝ T 1/7 RH for T RH T eq RH , the decay processes implies a harder dependence, m φ ∝ T −1 RH . Notice also that for BSM scales above GUT scale, Λ 10 16 GeV, it becomes almost impossible to generate the correct amount of dark matter, neither from scattering nor from the inflaton decay, both processes being too slow to compete with the expansion rate driven by H(T ).
We also show in Fig. 4 the allowed region in the plane (m φ , d) assuming disformal couplings only (c = c X = 0) for different values of Λ and T RH = 10 11 GeV. We observe that fairly natural values of d, of the order of loop factors 1/(4π) 2 , make it possible to obtain dark matter in sufficient quantity while avoiding overabundance. Still, larger values of Λ imposes relatively heavy dark matter, above the TeV-PeV scale to respect the cosmological observations.

Non-instantaneous reheating effects
Until now, we have considered a thermal Universe where reheating took place instantaneously, i.e. the energy available in the last oscillations of the inflaton was instantaneously transferred to the radiative bath at t (Γ Φ ) −1 , where Γ Φ is the decay width of the inflaton. However, we know that non-instantaneous perturbative phenonomena can have serious consequences on the thermal evolution of the Universe [43], and even more on the production of dark matter in its first instants if the modes of production show a large dependence on the energy of the processes involves [46], is the case for the disformal (conformal) scenario with a rate R(T ) ∝ T 12 (T 8 ) respectively.
More generally, the dark matter production during the reheating may not be negligible, especially when R(T ) ∝ T n+6 with n ≥ 6, due to the effects of noninstantaneous reheating [43], non-instantaneous thermalization in [42], and non-quadratic inflaton potentials during the reheating stage [46]. In our case the contributions from the conformal and disformal coupling terms proportional to c X and d, respectively, correspond to n = 6, whereas the reaction rate of the conformal coupling term proportional to c corresponds to n = 2. The interference term corresponds to n = 4. In this section we discuss the non-instantaneous reheating effect, while assuming instantaneous thermalization and no preheating contributions.
As a specific example, we consider the Starobinsky model for inflation [57] where the inflaton oscillation is described by V (Φ) = 1 2 m 2 Φ Φ 2 after the end of inflation.
Then, we may use the result for the enhancement of the DM production discussed in Ref. [8]. Solving the complete set of combined equations for the inflaton density ρ φ , the radiation density ρ R and the dark matter has been carried out and analyzed in [43] for any kind of dark matter production cross section and more recently in [46] for any type of inflationary potential. To summarize these works, we just need to understand that the reheating process being non-instantaneous, the temperature of the primordial plasma evolves from a null temperature to a maximum value T max before decreasing until the radiation density ρ R catches the inflaton density ρ φ , defining the thermal era, happening at the reheating temperature T RH . The evolution between T max and T RH is rather complex, but the main point is that the production of dark matter for cross-sections with a large temperature dependence, of the order T n+6 with n ≥ 6, is largely affected by the maximal temperature as most of the dark matter is produced at this instant. In comparison with an instantaneous treatment, there is a boost factor which is a function of T max /T RH . We summarize the results in the following paragraph.
The maximal temperature T max and T RH are obtained as T max = 45 32 1.6 × 10 13 GeV × 106.75 g * (T max ) 1.9 × 10 11 GeV × 106.75 g RH where again T RH is defined by ρ Φ (T RH ) = ρ R (T RH ), and we assume g * (T max ) = g * (T RH ) = g RH in the following analysis. We have used the inflaton decay width Γ Φ HH ≡ y 2 m Φ /8π with y ≡ µ Φ /m Φ from Eq. (35), where m 2 Φ 24π 2 A S * M 2 P /N 2 * with ln(10 10 A S * ) = 3.044 [28,58] and N * 55 + 0.33 ln y [8]. Then, for n = 6, we obtain the boost factor B scatt ≡ n non−inst φ (T RH )/n φ (T RH ) given by where f 1.2 to match the numerical results. Notice that for n = 2, which is the case of d = c X = 0 with c = 0, we do not have such an enhancement, since the DM production is dominated at T RH . Figure 5 shows the contours of Ω scat φ h 2 + Ω dec φ h 2 = 0.1 where only the disformal coupling contributes, namely c = c X = 0 and d = 0, and we take d = 1, taking into account the effect of non-instantaneous reheating just discussed above. Notice that in the bottom-right corner of the figure, the dark matter mass is in excess of T RH , and thus the scattering contributions get further suppressed, which is however irrelevant for smaller Λ ( 10 15 GeV), since the decay contribution dominates in that parameter space, the domination occuring from Eqs. (40) and (46) for µ Φ 3 × 10 −3 m Φ . It would be interesting, in this framework, to compare, the disformal to the conformal production of dark matter.

Combining conformal and disformal coupling
The conformal and disformal couplings may coexist. For instance, we can take both d and c being nonzero, while c X = 0. Incorporating nonzero contribution from the c coupling, we obtain the scattering reaction rate where β c 1.1 × 10 −2 c 4 is given in Eq. (70) and β d 4 d 2 is given in Eq. (67). The quantity β cd , arising from interferences between conformal and disformal, couplings is given by The radiative inflaton decay is also affected by the conformal coupling contributions, and thus we obtain giving a branching ratio to a DM pair of Combining the relic abundance produced by scattering integrating Eq. (26) 9 with the rate (48), combined with the boost factor due to non-instantaneous thermalization (47) and adding the decay process (50) where we used the following results (see, for instance, Ref. [46]) which is the main result of our work. Eq. (53) gives the total amount of dark matter produced in a model with a combination of disformal (d) and conformal (c) couplings, taking into account production through scattering from the thermal bath and radiative decay of the inflaton, together with instantaneous effects due its the finite width. We illustrate our results in Fig. 6 where we plot the region of the parameter space allowed in the plane (c, d) 9 By taking the non-instantaneous reheating into account, one ob- with H(T ) = (2/5)Γ Φ HH (T /T RH ) 4 , instead of Eq. (27).
for different dark matter masses m φ , fixing T RH = 10 11 GeV and Λ = 10 14 GeV. We clearly distinguish the two regimes, and for which values of c the conformal couplings begin to dominate over the disformal one. For our choice of parameters, the decay rate dominates the production of φ in Eq. (53). It is interesting to notice that for any dark matter mass, there exists a point in the parameter space, with reasonable value of c and d, respecting the cosmological constraint despite the large suppression due to high BSM physics scales.
If one looks into more details at the zone of influence of the disformal coupling versus the conformal one, we find that for both processes gives a similar contribution to the relic abundance, smaller values of Λ favouring of course the disformal production. We illustrate this situation in Fig.  7 where we take c = 100 and d = 1. For Λ below ∼ 10 15 GeV we recognize the characteristic of disformal production observed in Fig.(5) whereas for Λ 10 15 GeV, the production begins to be independent of Λ, which is a clear signature of a conformal production of dark matter.

V. CONCLUSIONS
We have shown that in models where the geometrical metric governing gravitational physics is different from the dynamical metric felt by the standard model particles via a scalar field φ, this scalar φ can play the role of dark matter. In this scenario, dark matter is produced via freeze-in and it is possible to respect cosmological constraints on the relic abundance of dark matter. Moreover, this can happen even though the suppression scale of the coupling between the scalar field and matter particles could be almost as large as the Planck scale.
The number density of scalar dark matter particles would then be produced in the early stages of the Universe, by a freeze-in mechanism due to its very feeble coupling to the Standard Model sector, i.e. due to the large suppression scale. In such a disformal dark matter scenario where the dark matter field would be disformally coupled to the standard model fields, the low production rate of φ would be counterbalanced by a large mass m φ , making disformally coupled scalars perfect candidates for experiments looking for superheavy dark components like ANITA or IceCube [48]. Assuming that DM is predominantly produced by 2 → 2 annihilations of SM particles, the Boltzmann equation for the DM number density can be written as

Acknowledgments
where the quantity on the right-hand-side R(T ) represents the temperature-dependent DM production rate per unit of volume and time. The rate can be expressed as a sum of the contribution of SM species of spin i to where N i is the number of the SM species of spin i. The partial rate R i can be expressed as with p j (E j ) is the 4-momentum (energy) of particles j = 1, 2, 3, 4 for processes 1 + 2 → 3 + 4 with 1, 2 being particles of the SM and 3, 4 dark matter states. f i represent the Bose-Einstein (i = 0, 1) and Fermi-Dirac (i = 1/2) statistics distribution functions. θ 13 and θ 12 are the angle formed by momenta of 1,3 and 1,2 respectively. The differential solid angle can be expressed as dΩ 13 = 2π d cos θ 13 . These kinematics quantities are related to the Mandelstam variables in the ultra-relativistic limit t = (s/2)(cos θ 13 − 1) and s = 2E 1 E 2 (1 − cos θ 12 ). More details can be found in the Appendices of Ref. [8].
B. PRODUCTION RATE: SCATTERING

Rate for a generic amplitude
Assuming an amplitude squared for the process i+i → DM+DM, where i denotes one SM particle of spin i, of the form 10 .
the integrated amplitude squared reads Taking the integral expression of Eq. (58), the contribution of a particle of spin i to the rate is c i nk (−1) k 2 2(n+k) Γ 2 (n + k + 2)ζ 2 (n + k + 2)T 2(n+k+2) 128π 5 (k + 1)(n + k + 1)Λ 2(n+k) × 1 , (i = 0, 1) , and the corresponding contribution to the relic density is given by 10M P m DM 256π 8 g 3/2 * 2 2(n+k) Γ 2 (n + k + 2)ζ 2 (n + k + 2)T 2(n+k)−1 RH (k + 1)(n + k + 1)(2(n + k) − 1)Λ 2(n+k) with m DM being the dark matter mass, where the first and second cases correspond respectively to the Bose-Einstein (i = 0, 1) and the Fermi-Dirac (i = 1/2) statistics for the initial state particles. We used the expression of the Hubble rate in terms of the SM temperature in the radiation domination era H(T ) = (g * π 2 /90) 1/2 T 2 /M P and considered constant relativistic degrees of freedom for simplicity g * = g RH . s 0 and ρ 0 c are the entropy density and critical density of the present time. Assuming that for each SM particle of spin i, the DM production amplitude squared is given by Eq. (59), the total relic density can be expressed as

Rate for disformal couplings
The amplitudes |M i | 2 for the processes i + i → DM+DM, where i denotes one SM particle of spin i, are 10 As s + t + u = 0 in the ultrarelativistic limit, our expression contains all the possible processes.
given by The total rate is given by with β d 4d 2 .

Rate for conformal couplings
As previously mentioned, in this case only the scalar particles contribute to the rate. The amplitude is given by For the case c = 0, the total rate is given by with β c X 9.74c 2 X . For the case c X = 0, the total rate is given by R c (T ) = c 4 π 3 2700 with β c 1.1 × 10 −2 c 4 .