Detectable gravitational waves from preheating probes nonthermal dark matter

We describe the challenges and pathways when probing inflaton as dark matter with the stochastic gravitational waves (GWs) signal generated during the (p)reheating. Such scenarios are of utmost interest when no other interaction between the visible and dark sectors is present, therefore having no other detectability prospects. We consider the remnant energy in the coherently oscillating inflaton’s zeroth mode to contribute to the observed relic dark matter density in the Universe. To fully capture the nonlinear dynamics and the effects of back-reactions during the oscillation, we resort to full nonlinear lattice simulation with pseudo-spectral methods to eliminate the differencing noises. We investigate for models whose behavior during the reheating era is of m 2Φ Φ 2 type and find the typical primordial stochastic GWs backgrounds spectrum from scatterings among highly populated inflaton modes behaving like matter. We comment on the challenges of constructing such viable inflationary models such that the inflaton will account for the total dark matter of the universe while the produced GWs are within the future GWs detectors such as BBO, DECIGO, PTA, AION-MAGIS, and CE. We also describe the necessary modifications to the standard perturbative reheating scenario to prevent the depletion of residual inflaton energy via perturbative decay

Interestingly, irrespective of the production mechanisms, the expected signals and constraints on the dark matter candidates are dictated by their interactions with the states of the standard model (SM).In this paper, we consider the scenario in which the DM relic abundance is set via the inflaton oscillations and decay, which is similar to the out-of-equilibrium freeze-in production mechanism.The cosmic inflation resolves the flatness and horizon problem and seeds the initial density fluctuations for large-scale structure formation [33][34][35][36][37], and is considered to be driven by a neutral scalar field -stable on cosmological scales -can also be a natural dark matter candidate.Several such possibilities have been considered in the literature so far .Generation of dark matter in the primordial Universe can be due to the incomplete decay of the inflaton (the remnant inflaton as the DM) or from the decaying inflaton -producing DM along with other SM degrees of freedom.The fact that only one scalar field can minimally explain both inflation (explaining the temperature fluctuations as seen in the cosmic microwave background radiation (CMB) [59]) and make up the currently observed DM abundance in the Universe, is needless to say, an economical idea worth our attention.
However, such a nonthermal production of dark matter is impossible to detect 1 in any conventional astrophysical or laboratory-based dark matter searches, as they do not communicate with the visible sector (SM) and are typically considered to be beyond the paradigm of testability via any other interactions.Thankfully, the gravitational waves generated during the nonlinear preheating phase after inflation can provide an alternative channel for DM searches.
The initial phase of inflationary (p)-reheating (see, e.g., Refs.[61][62][63]) is marked by resonant growth of the inflaton fluctuations and other fields that are coupled to the inflaton.The preheating soon gives rise to large inhomogeneties in the field fluctuations.These inhomogeneities act as a source of classical stochastic background gravitational waves (SGWB) [64,65].The production of such 1 Such non-thermal particle production was proposed to be tested via dark radiation or N eff measurements during Big Bang Nucleosynthesis (BBN) & CMB era, in context to neutrino anomalies [60].
Inspired by the considerations mentioned above, this study addresses the issue of inflaton as dark matter, with a focus on the possibility of detection of the GW signals generated during the preheating phase.Our findings reveal that, with specific parameter choices, these GWs produced during preheating can be observable through GW detectors such as DECIGO, BBO, CE, and ET.Additionally, when exploring the inflaton as a candidate for dark matter, it places constraints on the parameter space of the inflaton coupling to fields.The significance of these alternative probes cannot be overstated, especially in contrast to traditional dark matter detectors, as they enable the exploration of non-thermal dark matter production mechanisms that involve limited or no interaction between the dark and visible sectors.
The paper is organized as follows: In Section II, we outline the basic prerequisites in the scalar field potential to ensure that the frequencies of the resulting GWs signals fall within the detectable spectrum.The challenges and prospects for achieving such potentials, which also satisfy the CMB observations on larger scales, are discussed in Appendix A. We study the generation of classical sub-Hubble gravity waves from large time-dependent inhomogeneities due to (p)-reheating in Sec III.In Section IV, we present the specifics regarding the relic inflaton's role as dark matter, and we describe how to reheat the universe, keeping a stable inflaton relic while also making way for the impending radiation-dominated universe in Sec V. Finally, we show the predictions of the GW spectrum in current and future detectors in Sec VI.This culminates in final discussions and conclusions in the concluding section.We have used M Pl = 1/ √ 8πG = 2.435 × 10 18 GeV to denote the reduced Planck mass.

II. THE SCALAR FIELD POTENTIAL DURING REHEATING
Our preheating analysis will mostly be agnostic about the details of inflationary models during the inflationary phase.Our only requirement is that the field Φ responsible for generating the gravitational waves behaves as the archetypal m 2 Φ Φ 2 during the reheating phase.Most of the current observationally favored inflationary models such as the original Starobinsky model [37] in the Einstein frame or its recent generalization as the α-attractor Emodel [79] or the T-model [80] with the exponent n = 2 behaves as quadratic chaotic inflationary potentials.Other typical cases which can give rise to such behavior are: 1 Hybrid inflationary models: In this scenario, the inflaton responsible for generating the primordial pertur-bations can be of any shape [81,82].The waterfall field Φ dictates the end of inflation and is relevant for (p)-reheating.So, in such a scenario, we can identify the waterfall field ϕ and its inhomogeneities which give rise to the gravitational waves.
2 Curvaton models: For such models, the field Φ is inflaton while the curvaton field seeds the primordial perturbations [83].
These theoretically well-motivated scenarios can alleviate the CMB constraints on the inflaton sector when looking for gravitational wave signatures from the primordial universe within the detectable range of various GWs detectors that are otherwise not accessible to CMB-based observations.The typical frequency of the observed GWs generated during preheating with a canonical inflationary model falls within ultrahigh-frequency (gigahertz) bands.However, as pointed out in [66], we can shift the frequencies towards the detectable ranges if we consider the effective mass m Φ of the scalar field Φ as a free parameter.
In this regard, it is worth mentioning that the production of GWs driven by a system of coherent oscillating scalars taking the effective mass as a free parameter is also studied in the context of axions in [84].However, since we are considering that the scalar field is also responsible for the inflationary universe, we have more constraints on model building.In the vanilla m 2 Φ Φ 2 model or other single field models such as α-attractor models (where the effective mass around the minimum is determined by the amplitude of the plateau Λ and the width M of the potential), m Φ is constrained from power spectrum normalization.For the remainder of the paper, we will use the m 2 Φ Φ 2 simplified model to exemplify the characteristics that emerge when treating the inflaton as dark matter and exploring the related GW signals from preheating.The details of preheating with the simple m Φ Φ 2 model will follow in the next section.We will outline the challenges in model building in Appendix A.

III. GRAVITATIONAL WAVES GENERATED DURING PREHEATING
The preheating phase is generally episodic; the growth of the quantum fluctuations of the fields in the first phase is mainly due to the resonance phenomenon, while most of the system's energy is still in the homogeneous inflaton condensate, and the backreaction of the produced fields is not significant.The resonance, however, soon promptly populates all the modes of the quanta, and back-reaction becomes relevant.The second stage is where we observe the effects of back-reaction: the highly populated modes behave like matter waves and start colliding in a phase dominated by processes of nonlinear dynamics.These scatterings of waves happen in a highly relativistic and asymmetric manner.These processes generate a SGWB.Unlike the inflationary gravitational waves (tensor perturbations) generated from the quantum vacuum and get stretched to classical stochastic waves, GW production during preheating is a consequence of the scattering of the classical inhomogeneities, thus providing us with a complementary channel to test the physics of the early Universe.To fully capture the nonlinear dynamics of preheating and treat the backreaction carefully, we must study all the phases in full numerical simulation.Among several publicly available numerical codes to study preheating dynamics with 3+1 dimensional lattice simulation [85][86][87][88][89] we work with PSpectRe [87] while taking the initial inhomogeneous seeds of the scalar fields to be generated similarly to Defrost [86].Being a pseudospectral code, PSpectRe is free of differencing noise arising due to approximating the derivatives.In our simulation, the scalar field evolution is computed using a second-order-in-time symplectic integration scheme with a fixed time step as default in PSpectRe.We followed [67,90] to compute the transverse-traceless part of the first-order metric perturbation.In synchronous gauge, the metric is written as: with The equation of the motion of the metric fluctuations are: where T T (with ϕ a = {ϕ, χ}).We then evolve the Fourier modes of the metric perturbations (and their sources in the stress-energy tensor) in Fourier space using a fourth-order-in-time Runge-Kutta scheme.The background evolution is self-consistently computed from the scalar field, only neglecting the backreaction of the metric fluctuations on the background.The GW energy density is given by where ⟨• • • ⟩ V denotes a spatial average over V. Finally, the spectrum of the energy density of GWs (per logarithmic momentum interval) observable today is [65][66][67][68]] where, respectively, the subscript 0 indicates quantities evaluated today, e at the end of simulation and * at the end of reheating.The equation of state (EoS) of the universe between the end of inflation to the end of reheating is denoted by w.The critical density of the Universe is Although the EOS at the end of reheating should ideally be radiationlike.However, this may not be the case for inflaton with a massive component as considered here but may reach values very close to 1/3 [91][92][93][94].Usually, perturbative reheating is considered in such cases to drain the inflaton energy density and bring the radiation-dominated universe.However, to identify the inflaton as a stable dark matter candidate, we require some mechanism to initiate the radiation-dominated universe without depleting the inflaton energy density.The details of such reheating mechanisms are described in Sec V. Hence, for our calculation, we will take w = 1/3.As a consequence, the redshifting of a physical wave vector, k, induced by preheating after inflation with corresponding energy scale H e (evaluated at a e ), corresponds to an observed frequency today as Now, for the simulation we will work with where |H| is produced due to inflaton decay.The dimensionless coupling g is generally taken to be small to avoid large radiative corrections during inflation.After the end of inflation, when the inflaton oscillates about the minimum of the potential, we can identify the effective mass parameters for any model that behaves as V (Φ) ∝ Φ 2 around the minimum.As mentioned earlier, we will take m as a free parameter; additionally, to have a sufficient amplitude of GWs, we also need to take the initial amplitude of the inflaton to be large enough (∼ M Pl ).This makes the model building a bit challenging in this case.For instance, in α-attractor models, to decrease the effective mass, we will also need to decrease the value of the α parameter.It, however, will decrease the initial amplitude of the inflaton field (Φ 0 ) at the start of the simulation.We mention some directions on the inflationary model-building aspects in the Appendix.The essential requirement for unifying the inflaton and dark matter sector is that of incomplete decay of the inflation, which is characteristic of preheating with massive inflaton.

IV. INFLATON AS THE DARK MATTER
We describe the post-inflationary dynamics of the inflaton field in this section: the end of inflation is marked by the oscillation of the inflaton field around the minima of the potential minimum, and during this time, due to oscillations of the inflaton field, the energy density of the inflaton get transferred to other fields that are coupled to the inflaton.Due to the inherent Z 2 parity, under which the inflaton transforms as odd, the inflaton cannot decay to any other particles, so perturbative reheating is not possible, following the preheating particle production [61] which in this case happens via, the g 2 Φ 2 |H| 2 coupling, we expect the inflaton to violently produce Φ particles which are also known as the parametric resonance, due to Bose-Einstein scalar field statistics.When the inflaton amplitude becomes smaller than m Φ /g, the preheating stops [61,62].
The oscillating mode of the stable inflaton stops transferring the energy density into any other sector at the end of the preheating era and thus may act as the dark matter, with the ratio of the average energy density of this oscillation mode to the number density of photons being given by [41,95] where Φ * is the oscillation amplitude at time t * when m Φ = H * , the Hubble parameter.We have also used the typical values of the total number of relativistic particle degrees of freedom g * ≃ 100 and preset entropic degrees of freedom g S,0 = 3.9.The averaged inflaton energy density for t > t * scales as TeV and g2 = 0.1.To satisfy the observed relic density of the universe, ξ dm,0 ≃ 5.5 × 10 −27 M Pl , and therefore we conclude that the oscillating mode has a negligible contribution to the energy density of the present Universe 2 .
One crucial issue we must address is how to reheat the universe while keeping the remnant inflaton relic as a stable dark matter.In the upcoming section, we will address the matter of perturbative reheating specifically tailored for scenarios involving the inflaton as dark matter.

V. INCOMPLETE DECAY OF INFLATON AND REHEATING
Preheating after inflation with m 2 Φ Φ 2 type potential around the minima with four-legged interactions, as we have considered in this work, is known to produce an incomplete decay of inflaton [91][92][93][94].Although including a three-legged decay term can lead to a complete decay [96], we bar such interactions by presuming that the Z 2 is (almost) exact symmetry for the inflaton.Thus, under these circumstances, we cannot expect the radiation-dominated era to commence; rather, a phase of a matter-dominated era reemerges after a few e-folds [91][92][93][94].In such cases, once preheating becomes inefficient (once the inflaton amplitude drops below m Φ /g or Φ 0 / √ q), the perturbative reheating with explicit inflaton decay terms are invoked for a complete decay of the inflaton and to bring in the radiation-dominated universe.But if the inflaton is to remain as the stable dark matter relic, we must forbid such perturbative decay.Several proposals are made in the literature in this direction: • In [39,40], the authors utilized the idea of plasma mass effects [97] that kinematically forbids the in-flaton decay in some part of the reheating phase.However, inflaton can decay in such cases once it is subdominant with respect to the radiation fluid.Consequently, such scenarios will not be helpful for the case of a unified description of inflaton as dark matter.
• Authors in [41] propounded the idea of tripleunification by proposing a thermal inflation period following preheating.This phase of second inflation reduces the amount of inflaton/dark matter to the observed level.However, such a phase inflation of thermal inflation will not be suitable when considering the production of GWs from preheating.
• In [48], the authors noticed that as the mass of the inflaton decay products is field dependent, employing it is possible to construct models where inflaton decay only during the initial stage of coherent oscillation while it is kinematically forbidden at late times.By employing a suitable discrete symmetry that facilitates such decay away from minima of the potential, they can show that the scenario can bring in radiation-dominated era while the inflaton particles thermalize in the process and eventually decouple and freeze out, yielding correct DM relic abundance.This scenario has been further considered for identifying the inflaton field as dark energy at late times in the context of warm inflation [57].
Below, we will outline, following [48], the appropriate reheating scenario to bring in the radiation-dominated universe with a stable inflaton acting as the dark matter relic.Here, we consider the inflaton is also coupled to fermionic fields Ψ + and Ψ − via the standard Yukawa interactions, and they satisfy the following discrete symmetry C 2 ⊂ Z 2 ×S 2 such that when the inflaton transform as Φ ↔ −Φ the fermions are simultaneously transformed under the permutation symmetry as Ψ + ↔ Ψ − .The resulting Lagrangian density is Now, because of the above discrete symmetry, the two fermions have the opposite Yukawa coupling to the inflation while having the same tree-level mass m f .The discrete symmetry is preserved at the minimum of the potential at Φ = 0; consequently, perturbative decay to the fermionic channel is only allowed.When m f > m Φ /2, these decays are kinematically blocked, and the inflaton is completely stable, which accounts for the present DM abundance.Away from minimum, the discrete symmetry is broken, causing mass splitting of the fermions m ± = |m f ± hΦ|.Perturbative decay is allowed when field values satisfy |Φ| ≳ m f /h.We also note that due to Fermi blocking, we could neglect the resonance effects of the fermions during preheating.Additionally, the nature of the interaction makes it possible to neglect the fermionic interaction at the early stages of the preheating when m ± ∼ h|Φ| ≫ m Φ , i.e., at the time of the production of GWs due to the parametric resonance via four-legged interactions.
To sum up the reheating epoch, we have the following scenario for the successive phases of inflaton decay: (i) preheating commences after the end of inflation when H ∼ m ϕ with an initial amplitude of the field Φ 0 .Until the amplitude drops below Φ 0 / √ q, preheating is efficient and serves as the primary channel for energy transfer from the inflaton to the decay products.(ii) When m Φ /g ≳ |Φ| ≳ m f /h, we will have perturbative inflaton decay via fermionic channels.(iii) Perturbative decay is prohibited, thanks to the discrete symmetry mentioned above, when inflaton amplitude drops below m f /h, allowing a stable relic inflaton component.Now, the total decay width of the processes is Γ Φ = Γ + + Γ − with the partial decay widths of the fermionic decay channels given by Then we solve the following Boltzmann equations for the perturbative inflaton decay with the expansion given by the following Friedmann equation In Fig. (1), We show the evolution of the densities of inflaton ρ Φ and the radiation bath ρ R at temperature T [ρ R = (π 2 /30)g * T 4 with g * being the effective number of relativistic degrees of freedom forming radiation bath] The reheating temperature is given by GeV.
Thus, the typical amplitude when perturbative inflaton decay is forbidden can be tuned by choosing the fermion mass and the Yukawa interaction.For a typical value of m Φ for detectable GWs, we see from ( 12) that the reheating temperature will respect the BBN bound [98][99][100].With an understanding of how to reheat the universe, keeping a stable inflaton relic, we will describe the details of our numerical simulation in the following section.and radiation (red-dashed) during the perturbative reheating process after preheating.The radiation-dominated universe commences with a stable inflaton relic working as the dark matter candidate.

VI. NUMERICAL ANALYSIS
We performed our simulations with 512 3 lattice points.We kept the dimensionless resonance parameter q(= g 2 Φ 2 0 /m 2 Φ ) fixed at q = 8 × 10 4 for all the simulations.Where the initial values of inflaton oscillation Φ 0 are taken to be 1 M Pl .We show the produced GWs spectra in Fig (2).As mentioned earlier, we are interested in the scenarios where the CMB observables do not constrain the quadratic potential during preheating; consequently, we choose the inflaton mass m Φ as a free parameter.We have taken its value for each run as {5 × 10 −6 M Pl , 5×10 −9 M Pl , 5×10 −12 M Pl , 5×10 −15 M Pl , 5× 10 −18 M Pl , 5 × 10 −20 M Pl }.With these values, the ratio of the average energy density of the oscillation modes to the number density of photons defined in (6), are found to be {8.8 × 10 −9 M Pl , 2.8 × 10 −10 M Pl , 8.8 × 10 −12 M Pl , 2.8 × 10 −13 M Pl , 8.8 × 10 −15 M Pl , 8.8 × 10 −16 M Pl }.Now to satisfy the measured value of the current dark matter per photon ξ dm,0 = 1.1 × 10 −27 M Pl , we require the value of coupling g to have values {1.4 × 10 6 , 2.5 × 10 2 , 4.5 × 10 −2 , 8 × 10 −6 , 1.4 × 10 −9 , 4.5 × 10 −12 }.We have plotted the results in (3).The solid red line in the m Φ -g 2 parameter space satisfies the observed dark matter to photon ratio.Notice that for the cases considered here, when we have inflaton mass below 3 × 10 −11 M Pl or 7.34 × 10 4 GeV, the couplings are well within the perturbation limit.However, we have restricted our simulations for those values of q, providing us with sufficient resonance without any numerical instabilities.Simulations for all g 2 values along the red line in Fig. (3) will result in a large resonance parameter, rendering the lattice simulation afflicted by numerical instabilities.Nevertheless, provided we have sufficient resonance; the amplitude of the GW spectrum only depends on the strength of the resonance parameter while the frequency range depends on the initial energy scale.Consequently, for the same initial energy, changing the coupling will produce a 'tower' of GW spectrum ET aLIGO CE BBN bound on the integrated GW spectrum FIG.2: Spectrum of gravitational waves generated during preheating projected on the detection sensitivity of various experiments.With take the dimensionless resonance parameter as (g 2 Φ 2 0 )/m 2 Φ = 10 4 and start our simulation when Φ0 = 1 M Pl for all the cases.We change the energy scale, keeping m as a free parameter with this initial setup.We start with mΦ = 5 × 10 −6 M Pl (red in the extreme right) which corresponds to an energy scale of (4.5 × 10 15 GeV) 4 to mΦ = 5 × 10 −20 M Pl (black in the extreme left) which corresponds to a energy scale of (4.5 × 10 8 GeV) 4 .This final value of mΦ = 120 MeV is chosen such that it is well above the bound on the minimum reheating temperature of 10 MeV from BBN predictions.The intermediate curves are obtained choosing the values of mΦ to be 5 × 10 −9 M Pl (blue), 5 × 10 −12 M Pl (green), 5 × 10 −15 M Pl (yellow) and 5 × 10 −18 M Pl (brown).The gravitational waves generated during low-scale inflation fall within the sensitivity range of detectors such as DECIGO, BBO, LIGO, etc.
within the same frequency range [101][102][103].Thus, we can say that those values of the coupling g 2 along the red line larger than the simulations performed here will also ensure parametric resonance.Consequently, the corresponding GW signals should also be detectable with various detectors (see Fig. 2).This implies that the coherently oscillating scalar field during (p)-reheating that behaves as m 2 Φ Φ 2 type can source a promising candidate for the dark matter.Additionally, we would like to comment that if we consider the production of GWs from coherent scalar while relaxing the scalar field as inflaton, we can have detectable GW signals across a much wider range of frequency bands [104].For instance, fitting the peak frequency of the red-shifted GWs signal with the energy scale and hence the effective mass, we find the following relation: Hence, depending on the scalar field mass, we can move to the low-frequency GWs detector bands.All these scenarios could be interesting to explore in light of recently detected GWs background by the North American Nanohertz Observatory for Gravitational Waves (NANOGrav) [105].A comprehensive analysis of the amplitude of GW spectrum involving various interactions between the inflaton and the Higgs field, along with the effects of various initial conditions and a broad range of parameters of interest, will be examined in a future publication (also see Refs.[101-103, 106, 107] for similar works along this direction).

VII. DISCUSSION AND CONCLUSION
In a nutshell, we studied the inflaton behaving as the dark matter of the Universe.We summarize our main findings below: FIG.3: Plot shows the parameter space in the mΦ-g 2 plane, which will satisfy the observed dark matter to photon ratio (redline).
For values of the parameters, values from hatched are allowed but will produce less Y dm,0 in that case; the inflaton will not account for the full DM.The parameters within the yellow-shaded region correspond to q > 1, which naively tells us the region where broad resonance is possible.However, in practice, for broad resonance, one must have q ≫ 1.We took q = 8 × 10 4 .The black dots correspond to the values used in our simulation.The frequency in the upper axis corresponds to the peak frequency of the GWs signal using the fitted relation in (13).
• We showed in a (p)-reheating scenario, the inflaton starts to oscillate, and its zeroth mode energy remnant behaves as nonthermal dark matter as our observable current Universe satisfying the relic density (cf.Eq. 6).We investigated the impact of resonance parameters on dark matter production via lattice simulation in the oscillating regimes for such an oscillating scenario.
• In such a nonthermal origin of dark matter without any interaction with the visible or the SM sector, it is impossible to have laboratory tests for the DM candidate.We showed that the generated gravitational waves could be a viable detection channel during such an oscillatory period with its typical spectrum, as shown in Fig 2.
• In Fig 3, we showed the allowed parameter space of the coupling g 2 and the effective mass of the coherent oscillating inflaton m Φ , which satisfy the observed DM mass to photon ratio.If we have sufficient resonance, the amplitude of the produced GWs is independent of the energy scale of the inflation (cf Fig ( 2)).Thus, we found that for m Φ ≲ 10 5 GeV, the generated stochastic signals of GWs will be of sufficient amplitude to be detectable in various planned and proposed detectors while also explaining the origin of dark matter as the remnant inflaton component.
• As depicted in Fig 2 , for the higher mass range (m Φ ≳ 10 8 GeV), the coupling required to satisfy the observed dark matter relic is excluded by perturbative limit.Moreover, the allowed coupling range will overproduce dark matter and, hence, not allowed when considering inflaton as dark matter.For the lower end, although the couplings required to satisfy the dark matter relic are within the perturbative limit, the resulting resonance parameters became so large that the simulations became numerically unstable.However, in principle, the GWs produced in such parameter regions should also be detectable.Thus, the lower-mass ranges m Φ ≲ 10 8 GeV provide an interesting parameter space for the DM probe using GWs detectors.
• We can probe mass scales (≥ TeV) otherwise not accessible to colliders or any other laboratory physics.
There is huge interest in the community to find the detection principles of high-frequency gravitational gaves [108][109][110].With more advanced detectors, although we may be able to probe high-scale inflationary scenarios, the paradigm of inflation as the dark matter may be in tension within the scope of simple four-legged interaction as considered here.
We have thus explored a novel alternative to the thermal WIMP dark matter scenario in context to the non-thermal production of dark matter in the early Universe.The inflaton itself starts to behave as dark matter in late times, with the energy density stored during the oscillations of the (p)-reheating era.Due to no other or negligible interactions being present between the inflaton sector (dark sector) and the visible sector, we showed our prospects of detecting such non-thermal dark matter candidates within the measurements of SBGW in current and future GW detectors.We envisage our results will be the commencement of huge activity in this direction, whereas an alternative to hunt for various non-thermal dark matter production mechanisms which are otherwise not traceable in laboratory-based or astrophysical searches and, therefore, not in traditional DM detectors due to very feeble or no coupling between the visible and the dark sectors but in GW detectors in the near future, particularly when LIGO has already put bounds [111], and LIGO-VIRGO has further searched for such SGWB [112] from primordial sources as described in [113].
Appendix A: Towards Inflationary model building for reheating with detectable GWs signal The effective mass in most canonical inflationary models is constrained by the amplitude of the scalar spectrum at the CMB scales.For instance, in the simple oneparameter chaotic m 2 Φ 2 potential, the value of the mass m is fixed at 10 −6 M Pl .Looking beyond, the extra parameter models may provide us with some freedom.The currently observationally favored alpha-attractor models, for instance, the super gravity-inspired T-models having the inflationary potential V (Φ) = Λ 4 tanh 2 (ϕ/M ) reduces to V (Φ) ∝ (Λ 4 /M 2 )Φ 2 around the minima.Hence, the effective mass can be written as m 2 eff = Λ 4 /M 2 , where Λ is the amplitude of the potential's plateau and M is a free mass scale describing the width of the potential.In this case, however, it turns out that the observed amplitude of the scalar perturbations imposes approximately Λ 4 ∝ M 2 , so we have m Φ ≈ 10 −6 M Pl for all values of M [114].Consequently, the effective mass m Φ is the same as the vanilla chaotic inflation3 .Thus, we must expand our search for suitable models in other directions.We list below some models where we can take the effective mass during preheating as a free parameter and describe the associated challenges.

Two-field hybrid inflation
To make the effective mass a free parameter, the authors in [66] assumed that the full potential during inflation is given by: During inflation, when Φ > M/h, σ stays in its false σ = 0. Thus, the effective potential during inflation is The auxiliary field σ becomes tachyonic when Φ = M/h and evolves towards its true vacuum at M √ λ.Assuming σ = ⟨σ⟩, the potential reduces to Now, in the limit m 2 ≪ h 2 M 2 /λ, The effective potential is given by where µ 2 ≡ h 2 M 2 /λ, and, consequently, the effective mass for the inflaton can be taken as a free parameter.However, if m 2 ≪ h 2 M 2 /λ and during inflation when Φ > M/h, the potential in (A3) will result in blue tilted scalar spectra (n s = 1).Although such blue-tilted spectra have implications for resolving the Hubble tension [115][116][117], such spectra are ruled out from observations.

A model in the framework of α-attractor
A broad class of inflationary models defined in terms of a parameter α in the context of supergravity are the well-known α-attractors.The aforementioned parameter determines the curvature and cutoff of these models.These models are known to predict some universal results of the inflationary observables, viz., the spectral index (n s ) and the tensor-to-scalar ratio (r).There are several inequivalent formalisms to introduce the cosmological αattractors which eventually leads to the following Einstein frame Lagrangian [118]: which after rescaling of the field ϕ Finally, using the scalar field redefinition ϕ √ 6α = tanh(Φ/ √ 6α), we arrive at the following form of Lagrangian in the Einstein frame with a canonically normalized kinetic term: which will have n step-like transitions.We will only concentrate on the n = 1 case for the rest of our analysis.We take x CMB = tanh ΦCMB √ 6α that determines the normalization of the potential in the flat region to for matching the predictions of the model with the CMB observables.
x Step = tanh ΦStep √ 6α determines the appearance of the step.The inflationary potential is found to be: We have expressed all the dimensionful quantities in the units of M Pl and V 0 = f 2 0 .In practice, we can use continuous step function, for instance, Θ(x − x 0 ) = (1 + tanh((x − x 0 )/∆))/2 for numerical stability when solving the equations for the background inflationary dynamics and perturbations.The effective mass around the minima: Now, being derivatives of the α-attractor models, these models with steps can easily produce the necessary CMB observables.For instance, taking α = 1 and A = 0.05 we can have n s = 0.9630 and r 0.05 = 0.003 while the amplitude of power spectrum fixes the value of V 0 to be 1.344 × 10 −10 .This will make the effective mass around the minima m Φ = 4×10 −6 M Pl .Increasing A will decrease the m Φ but will require tuning the other parameters.In all cases, the only requirement for the preheating scenario to work is that the position of the step is at higher field values compared to the initial amplitude of the inflaton for simulation (Φ 0 = 1M Pl ).In the above case, we have taken Φ CMB = 10 and Φ Step = 4.Such models can also boost the scalar power spectrum, which has interesting applications for Primordial black holes and secondary gravitational waves [119,120].It is important to mention that to make m Φ as low as in the MeV range for the above model, we will need to take an extremely large height of the step.This may result in a kination phase after the end of inflation that can make the relation between the effective mass of inflaton and the observed GWs signal more convoluted.At this point, we must emphasize that the models listed above are only suggestive ways of making the effective mass a free parameter and looking for a more realistic model is still an open issue.We adjourn the details of model building and parameter explorations for future work.

FIG. 1 :
FIG.1:We plot the comoving density of the inflaton (blue) and radiation (red-dashed) during the perturbative reheating process after preheating.The radiation-dominated universe commences with a stable inflaton relic working as the dark matter candidate.

FIG. 4 :
FIG. 4:An illustration of the scalar potential for α = 1 and A = 0.5 (red), 0.05 (green) and 0.1 (blue).The height of the potential at the position of the step determines the effective mass of the potential around the minima.