# The Linear Regime of Tachyonic Preheating

Niko Koivunen  Eemeli Tomberg  and Hardi Veermäe

###### Abstract

Tachyonic preheating is realized when the inflaton repeatedly returns to a
convex region of the potential during the post-inflationary oscillating phase.
This will induce a strong tachyonic instability and lead to a rapid
fragmentation of the coherent field that can complete within a fraction of an
$e$-fold. In this paper, we study the linear regime of this process in a
model-independent way. To this purpose, we construct simplified models that
provide an analytic Floquet theoretic description of mode growth. This
approach captures the essential features of well-motivated tachyonic
preheating scenarios, including scenarios in which the inflaton is part of a
larger scalar multiplet. We show that tachyonic preheating is efficient if the
field excursions are sub-Planckian, can produce gravitational waves in the
frequency range of current and future gravitational wave interferometers, and
can be consistent with any experimentally allowed tensor-to-scalar ratio.

##  1 Introduction

The inflationary paradigm resolves several outstanding problems in Big Bang
cosmology [1, 2, 3, 4, 5, 6, 7, 8] and successfully explains the primordial
density perturbations inferred from the cosmic microwave background (CMB)
measurements by the Planck satellite [8].

During inflation, cosmic expansion is driven by the inflaton’s potential
energy. This energy must eventually be transferred to Standard Model (SM)
particles. This process is called reheating [9, 10]. The initial non-
perturbative phase of reheating is dubbed preheating111The nomenclature
differs in literature. Sometimes the distinction is made between the non-
perturbative stage of particle production (preheating) and perturbative
inflaton decay that produces the thermal bath of SM particles (reheating). and
can proceed through parametric resonance [11, 12, 9, 13, 10, 14] or through a
tachyonic instability [15, 16, 17]. In the first case, the post-inflationary
oscillating inflaton induces time-dependence in the masses of the fields it is
coupled to. This periodically changing mass leads to explosive particle
production through a resonant instability. On the other hand, tachyonic
preheating occurs when the squared mass of the fields becomes negative,
triggering particle production due to a tachyonic instability.

Tachyonic preheating has been studied in the context of hybrid [15, 16, 18,
19, 20, 21, 22], hilltop [23, 24, 25, 26, 27, 28], small field [29], and
plateau inflation [30, 31, 32, 17]. In single field inflation with a potential
that is convex at inflationary scales, tachyonic preheating generally occurs
when the Hubble scale at the end of inflation is much smaller than the
inflaton’s effective mass. To have such a sub-mass scale inflation, the
potential should become flatter as the field moves away from the minimum. Such
setups are naturally realized in plateau inflation. When the condition for the
Hubble scale is satisfied, the field repeatedly returns to the tachyonic
plateau region as it oscillates. These scenarios also tend to be in good
agreement with the increasingly tightening experimental bounds on inflationary
parameters [33]. Flattened plateau potentials are quite common and appear, for
instance, in alpha-attractor models [34, 35, 36, 37, 38, 39], and in Palatini
formulations with non-minimal gravitational couplings [40, 41, 42, 43, 44].
Also models in [1, 45, 46, 47] give rise to flat potentials, but these are
incapable of producing tachyonic preheating.

The current study provides a simple analytic model-independent description of
the initial linear phase of inflaton fragmentation due to a dominant tachyonic
instability. It complements our earlier numerical study [17]. Since the
timescales associated with fragmentation in tachyonic preheating are much
shorter than the Hubble time, we can neglect the expansion of space. The
coherent background then oscillates with an almost constant frequency and
amplitude, and the instability of each mode can be understood using Floquet
theory. The modes then follow the flat space equation

$\delta\ddot{\phi}_{k}(t)+\quantity(m_{\rm
eff}^{2}(t)+k^{2})\delta\phi_{k}(t)=0\,,$ (1.1)

with a time-dependent effective mass $m_{\rm eff}^{2}$ that becomes negative
during a part of each period. Since $m_{\rm eff}^{2}$ changes sign, the models
considered here contain elements of both parametric resonance and tachyonic
instability. Moreover, these situations do not obey the Mathieu or the Lame
equation, commonly encountered in discussions of parametric resonance [12, 9,
13, 10, 14]. We note that similar scenarios involving periodic sign-flips of
the effective mass, dubbed the flapping resonance, are encountered in the
context of axion models [48, 49].

Even in the linear phase, describing fragmentation can be a computationally
demanding task, since in the absence of a general analytical solution to (1.1)
the evolution of each mode must be tracked numerically. For this reason, we
propose simple analytically solvable models constructed by approximating
$m_{\rm eff}^{2}$ with simpler time-dependence, where $m_{\rm eff}^{2}$ is
constant and negative on the plateau and positive near the minimum of the
potential. We choose the positive part as a delta function or a top-hat like
box. These simplified models can be understood as temporal analogues to the
Kronig-Penney model [50] that describes electron conduction in a periodic one-
dimensional spatial lattice. The obtained analytic picture illuminates how
quantitative features of a given inflation model translate into the spectrum
of leading instability bands. In particular, this approach provides model-
independent analytic estimates for the characteristics of the instability,
such as the fastest growing mode and its growth rate.

Our approach does not describe processes during the final non-linear stage of
preheating, such as the creation of secondary peaks due to rescattering,
formation of oscillons or the eventual thermalization of the inflaton
particles. These processes have so far been studied using lattice methods [51,
52, 39, 53, 54]. We also do not consider the decay of the fragmented inflaton
into the SM thermal bath. Nevertheless, the ideas discussed here may be
extended to include some of the aforementioned effects. As an example, we
consider preheating when the inflaton is embedded into an $O(N)$ multiplet.

This paper is structured as follows. In section 2, we give the general set-up
required for tachyonic preheating with an oscillating inflaton and outline the
Floquet theoretical basis for perturbation growth. Section 3 introduces the
simplified models for which the perturbation growth is analytically solvable.
The structure of instability bands is discussed in detail in section 4 with
the aid of the simplified models as well as a numerically worked out realistic
scenario. A multi-field scenario is studied in section 5 and phenomenological
implications of our results are discussed in section 6. We conclude in section
7. Some technical details are given in the appendix. Throughout this paper we
use natural units $\hbar=c=1$ and set the reduced Planck mass to unity,
$M_{\rm P}=1$.

##  2 General considerations

We will study the dynamics of the inflaton $\phi$ arising from the Einstein
frame action

$S=\int{\rm
d}^{4}x\sqrt{-g}\left[\frac{1}{2}R-\frac{1}{2}(\partial\phi)^{2}-U(\phi)\right]\,,$
(2.1)

where $g$ denotes the determinant of the metric tensor, $R$ is the Ricci
scalar, and $U$ is the scalar potential. In a
Friedmann–Lemaître–Robertson–Walker spacetime, the equation of motion of the
coherent background field $\bar{\phi}$ and the Friedmann equation are

$\ddot{\bar{\phi}}+3H\dot{\bar{\phi}}+U^{\prime}(\bar{\phi})=0\,,\qquad
3H^{2}=\bar{\rho}\,,$ (2.2)

respectively. Above, $\bar{\rho}\equiv\dot{\bar{\phi}}^{2}/2+U(\bar{\phi})$ is
the energy density of the scalar, $H\equiv\dot{a}/a$ is the Hubble parameter,
$a$ is the scale factor, and a dot denotes a derivative with respect to the
cosmic time. These equations describe both inflation and the post-inflationary
oscillations of $\bar{\phi}$ around the minimum of $U$. However, to transition
into the hot Big Bang era, the energy density of the background field then has
to be transferred into radiation and, ultimately, into the SM degrees of
freedom.

The early stages of this process are characterized by the growth of the scalar
perturbations $\delta\phi\equiv\phi-\bar{\phi}$ due to the time-dependence of
the background. This non-perturbative particle production after inflation is
called preheating, and it ultimately leads into the fragmentation of the
background field. We focus on the early linear stages of this process during
which each mode of the perturbed field evolves according to

$\delta\ddot{\phi}_{k}+3H\delta\dot{\phi}_{k}+\omega_{k}^{2}\delta\phi_{k}=0\,,\qquad\omega_{k}^{2}\equiv
k^{2}/a^{2}+U^{\prime\prime}(\bar{\phi})\,.$ (2.3)

As these perturbations grow, interactions between the modes become
increasingly important, and the remaining evolution has to be studied using
non-linear methods, such as lattice simulations or thermal field theory [51,
52, 39, 53, 54].

Figure 1: A schematic of the potential _(left panel)_ and the effective
squared mass _(right panel)_ exhibiting a post-inflationary tachyonic
instability. In the blue region, the tachyonic instability is active, while
field excursions through the red region induce rapid changes in the mass and
contribute to the parametric resonance. The gray region in the left panel
correspond to the inflationary regime and there the potential does not have to
be specified for the purposes of the current study.

Preheating is more often studied in a context where the inflaton $\phi$
oscillates in a quadratic potential, $U(\phi)\sim\phi^{2}$, and decays to
another scalar field $\chi$ through a coupling of the form
$g\phi^{2}\chi^{2}$, or the inflaton potential is quartic,
$U(\phi)\sim\phi^{4}$, and the coherent inflaton fragments. The perturbations
will then obey a Mathieu or Lame equation, respectively, and the instability
responsible for particle production is the parametric resonance [12, 9, 13,
10, 14]. In both cases, perturbations predominantly grow when $\phi$ crosses
zero in the effective broad resonance regime. Instead, we are interested in a
different situation, in which $U^{\prime\prime}<0$ is possible away from the
minimum of the potential and leads to a tachyonic instability [15, 16]. Notice
that, around the minimum, we must still have $U^{\prime\prime}>0$ by
construction, so the tachyonic instability cannot be active at all times if
the background is to oscillate.

A schematic of a typical potential that gives rise to tachyonic preheating is
shown in the left panel of Fig. 1. The right panel of Fig. 1 depicts the
corresponding temporal evolution of the effective mass
$U^{\prime\prime}(\bar{\phi}(t))$. One can observe long stretches in the
tachyonic $U^{\prime\prime}<0$ region followed by brief transitions through
the non-tachyonic $U^{\prime\prime}>0$ minimum. The height $\Delta U_{i}$ and
width $\Delta\phi_{i}$ of the corresponding regions give a rough estimate for
the effective mass: $U^{\prime\prime}\approx 2\Delta
U_{1}/\Delta\phi_{1}^{2}\equiv m_{1}^{2}$ and $U^{\prime\prime}\approx-2\Delta
U_{2}/\Delta\phi_{2}^{2}\equiv-m_{2}^{2}$ in the non-tachyonic and the
tachyonic region, respectively. The time $T_{1}$ spent transitioning through
the non-tachyonic region is of the order $m_{1}$ and must satisfy $(1+\Delta
U_{2}/\Delta U_{1})^{-1/2}\leq m_{1}T_{1}/2\leq(\Delta U_{1}/\Delta
U_{2})^{1/2}$, computed from the maximal and minimal velocities in the limit
of negligible Hubble friction. Throughout the rest of the paper we assume that
the potential is symmetric, so that $U^{\prime\prime}$ is periodic with period
$T=T_{1}+T_{2}$, while the background field $\bar{\phi}$ has a period $2T$.

We do not consider the details of the inflationary epoch, and thus, the shape
of the potential beyond the oscillating region is not specified in Fig. 1.
Nevertheless, the Hubble rate during the oscillating period is roughly $\Delta
U_{1}+\Delta U_{2}\approx 3H^{2}$. In these scenarios, slow-roll ends due to
the violation of the second slow-roll condition $U^{\prime\prime}/U\lesssim 1$
implying that the potential should satisfy $\Delta\phi_{2}^{-2}\Delta
U_{2}/(\Delta U_{1}+\Delta U_{2})\gg 1$. This condition will generally
guarantee that Hubble friction is weak in the tachyonic region. For fast
preheating, we want this to be true also near the minimum, so we require that
the time spent in the non-tachyonic region does not exceed the Hubble time,
$T_{1}\ll H^{-1}$. Using estimates from the discussion above, we find that
this is satisfied when $\Delta\phi_{1}(1+\Delta U_{1}/\Delta U_{2})^{1/2}\ll
1$. This condition is less relevant because, in most cases, the time spent in
the tachyonic region is significantly longer, i.e., $T_{1}\ll T_{2}$, since
the field velocity is much lower there. In all, we can formulate the
sufficient condition

$(\Delta\phi_{1}+\Delta\phi_{2})(1+\Delta U_{1}/\Delta U_{2})^{1/2}\ll 1$
(2.4)

for effective tachyonic preheating. In particular, for $\Delta
U_{1}\lesssim\Delta U_{2}$, it is sufficient that the field value at the end
of inflation, given roughly by $\Delta\phi_{1}+\Delta\phi_{2}$, is sub-
Planckian.

Preheating can be extremely rapid if the field repeatedly returns to the
tachyonic region, causing the complete fragmentation of the coherent
background within a fraction of an $e$-fold. As the processes under study take
place in sub-Hubble timescales, we can neglect cosmic expansion altogether,
and the perturbation equation (2.3) simplifies to the Hill equation [55]

$\delta\ddot{\phi}_{k}+\omega_{k}^{2}\delta\phi_{k}=0\,,$ (2.5)

while the background (2.2) obeys $\ddot{\bar{\phi}}+U^{\prime}(\bar{\phi})=0$.
We remark that $\omega_{k}^{2}(t)$ depends on the energy density of the
background $\bar{\rho}$, which, even in Minkowski space, will decrease due to
the feedback from the perturbation growth. Nevertheless, if $\bar{\rho}$ is
not damped considerably during a single oscillation, its evolution can be
studied in the adiabatic approximation by defining an effective equation of
state by taking suitable averages over a single period [56, 32, 17]. Below, we
will work in the limit where both the Hubble friction and backreaction can be
neglected during an oscillation.

###  2.1 Floquet theory

For a periodic $\omega_{k}^{2}$, the solutions of Eq. (2.5) are governed by
Floquet theory, which invites us to look for quasiperiodic solutions of the
form

$\delta\phi_{k}(t+T)=e^{\lambda_{k}T}\,\delta\phi_{k}(t)\,,$ (2.6)

where $\lambda_{k}$ are the Floquet exponents for the mode $k$. We define the
growth rate for the mode $k$ as

$\mu_{k}\equiv\max\real\lambda_{k}.$ (2.7)

When $\mu_{k}>0$, the leading component of the general solution will grow
exponentially, which manifests physically as the production of
$\phi$-particles. Solving (2.5) amounts to finding $\lambda_{k}$ and, in
particular, $\mu_{k}$, for all values of $k$.

Given two independent solutions $u_{1}$ and $u_{2}$ of (2.5) we can compute
the growth exponents using the monodromy matrix

$G=w^{-1}(0)w(T)\,,\qquad\mbox{where}\qquad
w(t)\equiv\begin{pmatrix}u_{1}(t)&u_{2}(t)\\\ \dot{u}_{1}(t)&\dot{u}_{2}(t)\\\
\end{pmatrix}$ (2.8)

is the Wronskian matrix. The $G$ is a constant matrix with eigenvalues
$e^{\lambda_{k}T}$, with the two exponents $\lambda_{k}$ being additive
inverses due to the conservation of the Wronskian. Thus, they can be solved as

$\lambda_{k}\,T=\pm{\rm acosh}\left(\frac{1}{2}{\rm tr}\,G\right).$ (2.9)

To find $\lambda_{k}$ in the numerical examples studied in section 4, we solve
$u_{1}$ and $u_{2}$ from $0$ to $T$ with initial conditions
$u_{1}(0)=\dot{u}_{2}(0)=1$, $u_{2}(0)=\dot{u}_{1}(0)=0$, so ${\rm
tr}\,G=u_{1}(T)+\dot{u}_{2}(T)$.

An important property of growth rates for inflaton perturbations is that the
$k=0$ mode must be stable, that is, $\mu_{k=0}=0$ [17]. In single field
inflation, the wavenumber at which the first instability band begins depends
on how the period reacts to changes in the energy density of the background
field: when $\partial_{\bar{\rho}}T>0$, then the first instability band begins
at $k=0$, while for $\partial_{\bar{\rho}}T<0$, there is a stability band
beginning at $k=0$ [17]. In detail, for small $k$, the Floquet exponent can be
expanded as

$\lambda_{k}=-\frac{i\pi}{T}+(k/a)\frac{\sqrt{W\partial_{\bar{\rho}}T}}{T}+\mathcal{O}(k/a)^{3}\,,$
(2.10)

and is thus completely determined by the abbreviated action (see also appendix
A)

$\displaystyle W(\rho)$ $\displaystyle=2\int^{\phi_{\text{amp}}}_{0}{\rm
d}\phi\sqrt{2(\bar{\rho}-U(\phi))}\,,$ (2.11)

where $\bar{\rho}=U(\phi_{\text{amp}})$ is the background energy density and
$\phi_{\text{amp}}$ is the field oscillation amplitude. In terms of the
abbreviated action, the duration of a half-oscillation of $\bar{\phi}$ is
$T=\partial_{\bar{\rho}}W$. The presently considered models possess an
instability band starting at $k=0$, that is, $\mu_{k}>0$ immediately as $k>0$
entering a wide band comprising of all the modes for which $\omega_{k}^{2}$ is
negative on the plateau. This situation resembles the more commonly studied
broad parametric resonance regime [9].

##  3 Simplified models

The common characteristics of the potentials supporting tachyonic preheating
scenarios suggest that the essential features of this process can be studied
by generic simplified models. In the following sections, we construct such
models for the growth of linear perturbations in Eq. (2.5) by studying
different approximations of $\omega_{k}^{2}$. As depicted in Fig. 1, in the
models of interest the background is characterized by long stretches in the
tachyonic regime followed by rapid crossings of the origin where the sign of
$\omega_{k}^{2}$ is briefly flipped. This behaviour can be qualitatively
captured by the ansatz

$\omega_{k}^{2}=-\Gamma_{k}^{2}+\sum_{j\in\mathbb{Z}}f(t-jT)\,,$ (3.1)

where $f(t)$ is a positive $k$-independent function peaked at $t=0$ describing
the temporal evolution of $U^{\prime\prime}$ around origin crossings and222To
include expansion, one must change $k\to k/a$ to account for redshifting of
the modes.

$\Gamma_{k}^{2}=\Gamma_{0}^{2}-k^{2}\,,$ (3.2)

where $\Gamma_{0}$ is a complex constant.333Where it matters, we choose
$\real\Gamma_{k}>0$. However, the system is invariant under
$\Gamma_{k}\leftrightarrow-\Gamma_{k}$. The background mode $k=0$ should be
stable by Eq. (2.10). This provides an additional constraint between
$\Gamma_{0}$, $T$, and $f(t)$. If the ansatz (3.1) describes evolution arising
from a specific potential, this determines $\Gamma_{0}$ and $T$, which then
depend on $\bar{\rho}$ (or $\phi_{\text{amp}}$).

Although we focus on models exhibiting tachyonic instabilities, by extending
our results to an imaginary $\Gamma_{0}$, our analytic approximations may be
adapted to preheating via parametric resonance only. This approach may be
useful when the linear regime is not governed by the Mathieu or Lame equation.

In the following, we will consider two separate ansätze for $f(t)$: the delta
model and its generalization, the box model.

###  3.1 Delta model

In the simplest scenario, we take the non-tachyonic phase to be infinitely
short, so the peak in $\omega_{k}^{2}$ is described by
$f(t)=\Lambda\delta(t)$, that is,

$\omega_{k}^{2}=-\Gamma_{k}^{2}+\Lambda\sum_{j\in\mathbb{Z}}\delta(t-jT)\,,$
(3.3)

where $\Lambda$ is positive and $j$ enumerates the periods. We will refer to
this ansatz as the _delta model_. Between the peaks,
$\omega_{k}^{2}=-\Gamma_{k}^{2}$ is constant and, by the mode equation (2.5),
the perturbations evolve as

$\delta\phi_{k}=A_{j}e^{\Gamma_{k}(t-jT)}+B_{j}e^{-\Gamma_{k}(t-jT)}\,.$ (3.4)

Continuity of $\delta\phi_{k}$ together with the mode equation implies that
the first derivative must jump as
$\Delta\delta\dot{\phi}_{k}=-\Lambda\delta\phi_{k}$ after encountering each
peak. Then the coefficients in Eq. (3.4) are joined together as

$\begin{pmatrix}A_{j+1}\\\ B_{j+1}\end{pmatrix}=G\begin{pmatrix}A_{j}\\\
B_{j}\end{pmatrix}\,\equiv\begin{pmatrix}\left(1-\frac{\Lambda}{2\Gamma_{k}}\right)e^{\Gamma_{k}T}&-\frac{\Lambda}{2\Gamma_{k}}e^{\Gamma_{k}T}\\\
\frac{\Lambda}{2\Gamma_{k}}e^{-\Gamma_{k}T}&\left(1+\frac{\Lambda}{2\Gamma_{k}}\right)e^{-\Gamma_{k}T}\end{pmatrix}\begin{pmatrix}A_{j}\\\
B_{j}\end{pmatrix}.$ (3.5)

Although $G$ is not constructed via the Wronskian matrix as in Eq. (2.8), its
eigenvalues correspond to solutions that change by a constant factor after
each period, so it is similar to the monodromy matrix. From Eq. (2.9) we
obtain the spectrum of growth rates

$\mu_{k}=\frac{1}{T}{\rm Re}\,{\rm acosh}\left(\frac{1}{2}{\rm
tr}\,G\right)=\frac{1}{T}{\rm Re}\,{\rm
acosh}\left(\cosh(\Gamma_{k}T)-\frac{\Lambda}{2\Gamma_{k}}\sinh(\Gamma_{k}T)\right)\,.$
(3.6)

Regardless of other details such as the inflaton’s potential, perturbation
growth must satisfy the condition (2.10) at $k=0$. This fixes

$\Lambda=2\Gamma_{0}\coth\left(\frac{\Gamma_{0}T}{2}\right)\overset{\Gamma_{0}T\gg
1}{\approx}2\Gamma_{0}\,,$ (3.7)

so model depends on two dimensional parameters, $T$ and $\Gamma_{0}$. The
number of parameters can be further reduced by appropriate rescalings to a
single dimensionless parameter, e.g. , $\Gamma_{0}T$ [17]. Other dimensionless
quantities such as $\lambda_{k}T$ can be expressed in terms of it.

Figure 2: Delta model growth rate of modes for $\Gamma_{0}T=2$ (blue) and
$\Gamma_{0}T=10$ (green) and $\Gamma_{0}T\to\infty$ (orange). The gray dashed
line shows the approximation (3.8) of the first peak for $\Gamma_{0}T=10$. To
the right of the dashed vertical line, the instability stops being tachyonic.
The dashed green line shows the box model with $T_{1}=T/4$, $\Gamma_{0}T=10$
and ${\tilde{\omega}}_{0}$ fixed by (3.19).

The growth rates of different modes for different values of $\Gamma_{0}T$ are
shown in Fig. 2. We observe a band structure in which unstable bands
($\mu_{k}>0$) alternate with stable ones ($\mu_{k}=0$). The tachyonic
instability appears in the first and dominant band that begins at $k=0$. In
the other bands, $k>\Gamma_{0}$ and thus $\omega_{k}^{2}>0$ at all times, so
the instability is not tachyonic but due to parametric resonance. Since
$\Gamma_{k}$ is imaginary there, the hyperbolic functions inside the brackets
of Eq. (3.6) become trigonometric, causing $\cosh(\Gamma_{k}T)$ to oscillate
and creating unstable resonance bands near the peaks where the subdominant
$\sinh$ term lifts the absolute value of the sum above one. Eq. (3.6) implies
that each band is terminated at $\Gamma_{k}T=j\pi i$, with $j$ a positive
integer. Therefore, the first band lies in the range
$0<k<\sqrt{\Gamma_{0}^{2}+\pi^{2}/T^{2}}$. Although there are infinitely many
unstable bands at higher $k$, their growth rate is lower than that of the
leading peaks, and they can often be neglected in practice. We discuss the UV
behaviour in detail in section 4.2.

The higher $k$ instability bands also become weaker with respect to the
tachyonic band when $\Gamma_{0}$ increases. In the limit
$\Gamma_{0}T\to\infty$, shown by the orange curve in Fig. 2, they disappear
altogether, and the growth rate approaches the fully tachyonic
$\mu_{k}\to\Gamma_{k}$, as expected. This limit corresponds to fast
fragmentation due to a constant tachyonic mass term appearing like, e.g. , in
hilltop inflation [15]. In the opposite limit $\Gamma_{0}T\ll 1$, the first
band gets broader but has a weaker instability.

In the mostly tachyonic case $\Gamma_{0}T\gg 1$, which describes well the
models we wish to study, the shape of the first peak is nicely approximated by

$\lambda_{k}\approx\Gamma_{k}+T^{-1}\ln\left[1-\frac{\Gamma_{0}}{\Gamma_{k}}\right],$
(3.8)

as can be seen from Fig. (2). The fastest growing mode lies at

$k_{\text{peak}}\approx\sqrt{2\Gamma_{0}/T}\left[1+\mathcal{O}\left(\Gamma_{0}T\right)^{-1}\right]$
(3.9)

and grows with the rate

$\mu_{\text{peak}}\approx\Gamma_{0}-T^{-1}(1+\ln(\Gamma_{0}T))+\mathcal{O}\left(\Gamma_{0}^{-1}T^{-2}\right).$
(3.10)

These approximations agree with our previous observations. In particular, we
can use them to compare a pure tachyonic instability, that is, a constant but
negative mass term $-\Gamma_{0}^{2}$, with an oscillating mass with a strong
tachyonic region. The purely tachyonic limit $\Gamma_{0}T\to\infty$ can be
obtained by sending the period $T$ to infinity while keeping $\Gamma_{0}$
constant, giving $\mu_{\text{peak}}\to\Gamma_{0}$, $k_{\text{peak}}\to 0$. In
comparison, once we start decreasing $T$, the first peak’s position shifts
towards larger $k$, and its height decreases. That is, the instability bands
become wider and flatter, as shown in Fig. 2.

####  3.1.1 Generating potential

Above, we postulated an idealized time-dependence (3.3) for $\omega_{k}^{2}$
that allowed us to study mode growth analytically. We can go one step further
and ask whether we can construct a potential that induces the time dependence
(3.3) of $\omega_{k}^{2}$. The answer is positive: (3.3) is generated with the
potential

$U=-\frac{\Gamma_{0}^{2}}{2}\left((|\phi|-\phi_{2})^{2}-\phi_{2}^{2}\right),$
(3.11)

if one also imposes the condition (3.7), i.e., that the $k=0$ mode is stable.
To see this, first note that the effective mass corresponding to the potential
(3.11) is
$U^{\prime\prime}=-\Gamma_{0}^{2}+2\Gamma_{0}^{2}\phi_{2}\delta(\phi)$. As the
background crosses the origin periodically, we have

$\delta(\bar{\phi}(t))=|\dot{\bar{\phi}}|^{-1}\sum_{j\in\mathbb{Z}}\delta(t-jT)\,,$
(3.12)

and comparing $U^{\prime\prime}$ to the ansatz (3.3), we have a match if
$\Lambda=2\Gamma_{0}^{2}\phi_{2}/|\dot{\bar{\phi}}(0)|$.

To make the match explicit, we want to express all other quantities in terms
of the parameters $\Gamma_{0}$ and $T$. As we work in the limit of negligible
Hubble friction, the background’s energy density
$\bar{\rho}=U(\phi_{\text{amp}})$ is conserved and thus
$|\dot{\bar{\phi}}(\bar{\phi})|=\sqrt{2(\bar{\rho}-U(\bar{\phi}))}$.
Integrating this gives the half-period of $\bar{\phi}$,

$T=\frac{2}{\Gamma_{0}}{\rm
acoth}\left(\frac{\Gamma_{0}\phi_{2}}{\sqrt{2\bar{\rho}}}\right)\,.$ (3.13)

The velocity at the minimum is $|\dot{\bar{\phi}}(0)|=\sqrt{2\bar{\rho}}$. By
combining everything, we find

$\Lambda=2\Gamma_{0}\coth\left(\frac{T\Gamma_{0}}{2}\right),\,$ (3.14)

which matches the condition (3.7) for the stability of the $k=0$ mode. We find
that, for the delta model, the stability of the $k=0$ mode is equivalent to
the existence of a potential from which the ansatz (3.3) for $\omega_{k}^{2}$
can be generated.

By relating the oscillation time $T$ and the amplitude, Eq. (3.13) gives the
final ingredient for determining background evolution. Nevertheless, by
starting from $\omega_{k}^{2}$, it is not possible to fix the absolute scale
of the field, but only the ratio $\phi_{\text{amp}}/\phi_{2}$ can be fixed.
Varying $\phi_{2}$ affects the height of the potential (3.11), as can be seen
in the bottom right panel of Fig. 3.

###  3.2 Box model

Generalizing the previous case, we can account for the duration of the non-
tachyonic phase $T_{1}$. This leads us to model the peak in $\omega_{k}^{2}$
with the ansatz

$f_{\rm
box}(t)=({\tilde{\omega}}_{0}^{2}-\Gamma_{0}^{2})\theta\left(T_{1}/2-|t|\right)-{\tilde{\Lambda}}\,\delta\left(|t|-T_{1}/2\right)\,,$
(3.15)

where $\theta$ denotes the unit step function. In the tachyonic phase with
length $T_{2}$, we now have
$\omega_{k}^{2}=k^{2}-\Gamma_{0}^{2}=-\Gamma_{k}^{2}$ as before, and in the
non-tachyonic phase with length $T_{1}$,
$\omega_{k}^{2}=k^{2}+{\tilde{\omega}}_{0}^{2}\equiv{\tilde{\omega}}_{k}^{2}$.
The $\delta$-function is now located at the transition between the phases. It
mimics the dip in $\omega_{k}^{2}$ before the onset of the tahyonic phase
seen, for instance, in Figs. 1 and 3. The stability of the $k=0$ mode imposes
a condition that fixes ${\tilde{\Lambda}}$ in terms of the other 4 free
parameters of the model: $T_{1}$, $T_{2}\equiv T-T_{1}$,
${\tilde{\omega}}_{0}$, and $\Gamma_{0}$.

The mode equation (2.5) can be solved using standard methods for a piecewise
constant $\omega_{k}^{2}$. The modes evolve as

$\delta\phi_{k}=\left\\{\begin{array}[]{ll}A^{(1)}_{j}e^{i{\tilde{\omega}}_{k}(t-jT)}+B^{(1)}_{j}e^{-i{\tilde{\omega}}_{k}(t-jT)},&\qquad|t-jT|<T_{1}/2\\\
A^{(2)}_{j}e^{\Gamma_{k}(t-jT)}+B^{(2)}_{j}e^{-\Gamma_{k}(t-jT)},&\qquad
T_{1}/2\leq|t-jT|<T/2\end{array}\right.\,.$ (3.16)

As in Eq. (3.5), the monodromy matrix is given by the transformation
connecting the $A^{(i)}_{j}$ coefficients. Appropriate matching, demanding the
continuity of $\delta\phi_{k}$ and the discontinuity
$\Delta\delta\dot{\phi}_{k}={\tilde{\Lambda}}\delta\phi_{k}$ during
transitions, then gives

$\displaystyle\frac{1}{2}{\rm tr}\,G$
$\displaystyle=\left[\cos(T_{1}{\tilde{\omega}}_{k})+\frac{{\tilde{\Lambda}}}{{\tilde{\omega}}_{k}}\sin(T_{1}{\tilde{\omega}}_{k})\right]\cosh(\Gamma_{k}T_{2})$
(3.17)
$\displaystyle+\left[\frac{\Gamma_{k}^{2}-{\tilde{\omega}}_{k}^{2}+{\tilde{\Lambda}}^{2}}{2\Gamma_{k}{\tilde{\omega}}_{k}}\sin(T_{1}{\tilde{\omega}}_{k})+\frac{{\tilde{\Lambda}}}{\Gamma_{k}}\cos(T_{1}{\tilde{\omega}}_{k})\right]\sinh(\Gamma_{k}T_{2})$

and, by Eqs. (2.7), (2.9), it gives the growth rates as
$\mu_{k}=T^{-1}\real\,{\rm acosh}({\rm tr}\,G/2)$. The stability of the $k=0$
mode, i.e., ${\rm tr}\,G(k=0)=-2$, gives

${\tilde{\Lambda}}={\tilde{\omega}}_{0}\tan\left(\frac{T_{1}{\tilde{\omega}}_{0}}{2}\right)-\Gamma_{0}\coth\left(\frac{T_{2}\Gamma_{0}}{2}\right)\,,$
(3.18)

eliminating one of the 5 parameters of the initial ansatz (3.15).

The band structure (3.17) of the box model is relatively similar to the one of
the delta model (3.6). The most noticeable difference is the $k$-space
modulation with period $k\approx\pi/T_{1}$ when $k\gtrsim{\tilde{\omega}}_{0}$
that arises due the introduction of a new scale $T_{1}$. A comparison between
the delta and the box models is shown in Fig. 2 by the solid and the dashed
green lines, respectively. Although these lines almost overlap for the first
three peaks, the fourth instability band is seen to be much weaker due to the
modulating terms. Since the leading peak in $\mu_{k}$ tends to be similar in
both models, the added complexity of the box model is most useful for
obtaining a more accurate description of the subleading peaks.

####  3.2.1 The $T_{1}\to 0$ limit

As a consistency check, we consider the limiting case of an instantaneous non-
tachyonic phase, $T_{1}\to 0$. In this process, we fix $\int f_{\rm
box}(t)\,{\rm
d}t=({\tilde{\omega}}_{0}^{2}-\Gamma_{0}^{2})T_{1}-2{\tilde{\Lambda}}$ and
allow for the possibility that ${\tilde{\omega}}_{0}\to\infty$. We find that
this limit reproduces both the Floquet spectrum (3.6) and the stability
condition (3.7) of the delta-model with

$({\tilde{\omega}}_{0}^{2}-\Gamma_{0}^{2})T_{1}-2{\tilde{\Lambda}}=\Lambda$
(3.19)

or equivalently, with $\int f_{\rm box}(t)\,{\rm d}t=\int f_{\rm
delta}(t)\,{\rm d}t$.

Moreover, this matching permits us to compute the small $T_{1}$ corrections to
the delta model: the trace of the monodromy matrix is

$\frac{1}{2}\Tr
G=\left(1-\frac{T_{1}\Lambda}{2}\right)\cosh(\Gamma_{k}T)-\frac{\Lambda}{2\Gamma_{k}}\left(1-\frac{T_{1}\Lambda}{4}\right)\sinh(\Gamma_{k}T)+\mathcal{O}(T_{1}^{2})\,,$
(3.20)

while zero-mode stability gives

$\Lambda=2\Gamma_{0}\coth\left(\frac{\Gamma_{0}T}{2}\right)-T_{1}\Gamma_{0}^{2}+\mathcal{O}(T_{1}^{2}).$
(3.21)

Notice that these approximations hold as long as
$T_{1}\tilde{\omega}_{k}\approx T_{1}k\ll 1$. Thus, it can fail for modes that
oscillate with frequencies similar or larger to $1/T_{1}$, as these modes
begin to probe the structure of the peak. In fact, this is a general limit of
validity for these approximate models. As we will show in section 4.2, their
growth rates $\mu_{k}$ decrease much slower at large $k$ than in scenarios
based on smooth physical potentials.

####  3.2.2 Generating potential

The potential that generates the box model (3.15) can be constructed by
glueing together quadratic functions,

$U(\phi)=\left\\{\begin{array}[]{lr}\frac{1}{2}{\tilde{\omega}}_{0}^{2}\phi^{2},&\qquad|\phi|\leq\phi_{1}\\\
-\frac{1}{2}\Gamma_{0}^{2}(|\phi|-\phi_{2})^{2}+C,&\qquad|\phi|>\phi_{1}\\\
\end{array}\right.\,,$ (3.22)

where the constant
$C=\frac{1}{2}m_{1}^{2}\phi_{1}^{2}+\frac{1}{2}m_{2}^{2}(\phi_{1}-\phi_{2})^{2}$
is fixed by demanding continuity at $|\phi|=\phi_{1}$. For the second
derivative of the potential

$U^{\prime\prime}(\phi)=\left(\Gamma_{0}^{2}(\phi_{1}-\phi_{2})-{\tilde{\omega}}_{0}^{2}\phi_{1}\right)\delta(|\phi|-\phi_{1})+\left\\{\begin{array}[]{l}{\tilde{\omega}}_{0}^{2},\qquad|\phi|\leq\phi_{1}\\\
-\Gamma_{0}^{2},\qquad|\phi|>\phi_{1}\\\ \end{array}\right.$ (3.23)

to give rise to (3.15), we must demand that

${\tilde{\Lambda}}=\frac{{\tilde{\omega}}_{0}^{2}\phi_{1}-\Gamma_{0}^{2}(\phi_{1}-\phi_{2})}{|\dot{\bar{\phi}}|(\phi_{1})}={\tilde{\omega}}_{0}\tan\left(\frac{T_{1}{\tilde{\omega}}_{0}}{2}\right)-\Gamma_{0}\coth\left(\frac{T_{2}\Gamma_{0}}{2}\right)\,,$
(3.24)

equivalent to (3.18), where, as with the delta model, we used
$|\dot{\bar{\phi}}(\bar{\phi})|=\sqrt{2(\bar{\rho}-U(\bar{\phi}))}$ and
integrated $|\dot{\bar{\phi}}(\bar{\phi})|$ to get the durations of each
epoch,

$T_{1}=\frac{2}{{\tilde{\omega}}_{0}}\atan\left(\frac{{\tilde{\omega}}_{0}\,\phi_{1}}{|\dot{\bar{\phi}}|(\phi_{1})}\right),\qquad
T_{2}=\frac{2}{\Gamma_{0}}{\rm
acoth}\left(\frac{\Gamma_{0}(\phi_{2}-\phi_{1})}{|\dot{\bar{\phi}}|(\phi_{1})}\right)\,.$
(3.25)

The first derivative of the potential $U^{\prime}$ is continuous if and only
if ${\tilde{\Lambda}}=0$. Again we find that imposing the condition (3.18)
allows us to consistently construct a potential from which the ansatz (3.15)
follows. Analogously to the delta model, (3.25) fixes the ratios
$\phi_{\text{amp}}/\phi_{2}$ and $\phi_{\text{amp}}/\phi_{1}$, completing the
matching between $\omega_{k}^{2}$ and the background evolution. The overall
scale of the field, and thus the height of the potential, is again left free.

##  4 Modelling the band structure

We aim to use the simplified models derived in the previous section to provide
a detailed description of tachyonic preheating that is applicable to realistic
scenarios on both the qualitative and the quantitative level.

As an example, we will use the simplified setups to model the instability
bands of preheating with the potential

$U(\phi)=U_{0}\tanh^{2}(\phi/\phi_{0})\,,$ (4.1)

to which we refer to as the $\tanh^{2}$ potential. This potential represents
the more general class of potentials with an exponentially flat plateau, that
is, $U\sim U_{0}(1-A\exp(-2\phi/\phi_{0}))$, when $\phi\gtrsim\phi_{0}$. Here
$A$ and $\phi_{0}$ are constants. All such potentials lead to similar
preheating phases [17]. Plateau potentials can arise, e.g. , from scalar field
models of Palatini gravity with non-minimal gravitational couplings [42, 43,
44], or string theory inspired alpha-attractor models [34, 35, 36, 37, 38].
The perturbations have a tachyonic instability when $\phi\gtrsim 0.67\phi_{0}$
with a corresponding mass scale

$m_{\rm th}^{2}\equiv U_{0}/\phi_{0}^{2}\,.$ (4.2)

Applying the condition (2.4) for an effective tachyonic instability, we find
that $\Delta\phi_{1}\sim\Delta\phi_{2}\sim\phi_{0}$, $\Delta U_{1}\sim\Delta
U_{2}\sim U_{0}$. As a result, we find that efficient tachyonic preheating is
possible when

$\phi_{0}\ll 1\,,$ (4.3)

consistent with the numerical analysis of [17] where the more precise limit
$\phi_{0}\lesssim 0.01$ was given. In this case, the potential (4.1) predicts
the scalar power spectrum $A_{s}\approx N^{2}m_{\rm th}^{2}/(3\pi^{2})$, the
spectral index $n_{s}\approx 1-2/N$, and the tensor-to-scalar ratio $r\approx
2\phi_{0}^{2}/N^{2}$, when there are $N$ e-folds of inflaton left after the
CMB pivot scale exits the cosmic horizon. With $N\approx 50$ and $m_{\rm
th}\approx 5\times 10^{-6}$, the model fits the latest observational
constraints [8, 57]. In such plateau models, the tachyonicity requirement
$\phi_{0}\lesssim 0.01$ together with the fixed $m_{\rm th}$ sets the limits
$r\lesssim 10^{-7}$ and $U_{0}^{1/4}\lesssim 10^{14}$ GeV, giving an upper
limit for the reheating temperature.

###  4.1 Leading peaks

Figure 3: Comparison of the delta (blue) and box model (purple) fitted to the
$\tanh^{2}$ potential (4.1) (orange) with amplitude
$\phi_{\text{amp}}=3\phi_{0}$. _Left panel:_ The corresponding growth rates
$\mu_{k}$. _Upper right panel:_ The time-dependence of $\omega^{2}_{0}$. The
arrows represent delta functions. _Lower right panel:_ The generating
potential. The solid lines indicate the region where the field oscillates.
Parameters for the delta model fit: $\Gamma_{0}T=6.72$ (giving
$\phi_{\text{amp}}=0.93\,\phi_{2}$), and for the box model fit:
$\Gamma_{0}T=4.88$, $\tilde{\omega}_{0}T=11.61$, $T_{1}=0.22\,T$,
$T_{2}=0.78\,T$ (giving $\phi_{\text{amp}}=0.88\,\phi_{2}$ and
$\phi_{2}=1.70\,\phi_{1}$). The $\tanh^{2}$ potential uses
$\phi_{0}=1.29\times 10^{-3}$, $U_{0}=4.15\times 10^{-17}$ (corresponding to
an initial oscillation amplitude $\phi_{\text{amp}}=3\phi_{0}$ [17], $m_{\rm
th}=5\times 10^{-6}$, and $T=4.47\times 10^{6}$). Fixing $T$ determines the
generating potentials up to a total scaling of $\phi$. In the lower right
panel, three alternatives are depicted: one where the height of the potential
is equal that of the $\tanh^{2}$ case ($\phi_{2}=6.06\times 10^{-3}$ for the
delta and $\phi_{2}=5.72\times 10^{-3}$ for the box), and two where the
$\phi_{2}$ is either $20\%$ higher or lower.

The speed of fragmentation dynamics, as well as the initial spectrum of
perturbations, is determined mainly by the highest peaks in $\mu_{k}$, which
lie at low $k$.444In the tachyonic cases considered here, the first $\mu_{k}$
peak (with the lowest $k$) is always the highest one. However, this is not
generally true, e.g. in the low amplitude non-tachyonic limits of the
$\tanh^{n}$ models studied in [17] and for the orthogonal fields considered in
section 5 one of the subsequent peaks rises above the first. In this section,
we will compare the leading instability bands of the $\tanh^{2}$ model (4.1)
to the simplified scenarios. To this purpose, we fix the parameters of the
simplified models as follows:

  1. 1.

We demand that all models in the comparison have the same period $T$. For the
$\tanh^{2}$ model, the period and the amplitude are related by $m_{\rm
th}T=\pi{\rm cosh}(\phi_{\rm amp}/\phi_{0})$. Matching the period $T$ leads to
similar $k$-space periodicity of the instability bands. In fact, this
agreement between the periods in $\mu_{k}$ becomes exact in the limit
$k\to\infty$.

  2. 2.

We fix the initial slope $\partial_{k}\mu_{k=0}$ of the growth exponents. For
the $\tanh^{2}$ model, it can be obtained analytically from Eq. (2.10) as

$\partial_{k}\mu_{k=0}=\sqrt{\sqrt{2}m_{\rm th}T/\pi-1}\,$ (4.4)

even though the rest of the $\mu_{k}$ must be computed numerically. For the
delta and the box models, the initial slope is easily obtained from Eqs. (3.6)
and (3.17), respectively.

This procedure completely fixes the parameters of the delta model while it
leaves two of the parameters free for the box model. The remaining parameters
of the box model were fitted by eye to match the leading peaks of the box and
$\tanh^{2}$ models as closely as possible.

An example is shown in Fig. 3, where we compare the $\mu_{k}$ spectrum,
$\omega_{k}^{2}$, and the generating potentials of the three cases. The
$\tanh^{2}$ spectrum is computed numerically. The delta model reproduces the
leading first peak reasonably well, but the following secondary peaks are
lower and narrower than the $\tanh^{2}$ counterparts. The additional freedom
of the box model allows for a practically perfect match for the leading peak
and yields a better fit for the secondaries, but we see additional modulation
in their heights, arising from the two competing time scales $T_{1}$ and
$T_{2}$. Such modulation is absent in the smooth $\tanh^{2}$ case. Also, the
models differ considerably in their large-$k$ behaviour, as is discussed in
the following subsection. In the best fit box model, the non-tachyonic regime
is relatively wide compared to the $\tanh^{2}$ model, $T_{1}=0.28\,T_{2}$.
This is necessary to differentiate the box from the delta and to make the
secondary peaks higher and wider.

Figure 4: The height and location of the leading first $\mu_{k}$ peak for
$\tanh^{2}$ model and for the corresponding delta fit, as explained in the
text. The asymptotic behaviour (dotted black lines) for the $\tanh^{2}$
potential is $\mu_{\text{peak}}T\approx
0.92+1.00\,\phi_{\text{amp}}/\phi_{0}$, $k_{\text{peak}}T/\pi\approx 1.13$,
and for the delta model, $\mu_{\text{peak}}T\approx
0.89+1.10\,\phi_{\text{amp}}/\phi_{0}$, $k_{\text{peak}}T/\pi\approx
1.02+0.083\,\phi_{\text{amp}}/\phi_{0}$. For comparison, the value of $T$
varies from $2.22\,m_{\rm th}^{-1}$ to $1218.05\,m_{\rm th}^{-1}$ as
$\phi_{\text{amp}}/\phi_{0}$ changes from $0$ to $7$.

In the $\tanh^{2}$ model, preheating dynamics is regulated almost entirely by
the first peak when $10^{-4}\lesssim\phi_{0}\lesssim 10^{-2}$ [17]. In Fig. 4,
we have compared the predictions for this peak for the $\tanh^{2}$ model and
the delta fit over a wide range of $\phi_{\mathrm{amp}}$ values. The
analytical delta model agrees remarkably well with the numerical results of
the $\tanh^{2}$ model. This implies that the first peak is largely determined
by the behaviour of $\omega_{k}^{2}$ on the tachyonic plateau, which the delta
model captures well. The delta model is completely analytic and relies mainly
on Eq. (4.4) for fixing $\Gamma_{0}$ from the slope of $\mu_{k}$ at $k=0$.
This method provides a robust way for estimating the quantities relevant for
tachyonic preheating with minimal computational effort. The box model does not
provide such simple fits because of its additional free parameters.

Figure 5: Floquet charts for the $\tanh^{2}$ potential _(left panel)_ , the
delta model _(middle panel)_ and the box model _(right panel)_. The red dashed
curve indicates tachyonicity on the plateau: above it,
$k^{2}+U^{\prime\prime}(\phi_{\text{amp}})<0$. For the $\tanh^{2}$ potential,
the orange curve indicates overall tachyonicity: above it $\omega_{k}^{2}<0$
at the edge of the plateau where $U^{\prime\prime}$ is most negative. Notice
the broadening of the secondary instability bands once a tachyonic
contribution is present. The corresponding potentials are shown in Fig. 3 and
the black dotted line indicates the amplitudes used there.

Fig. 5 shows the Floquet charts for the potentials of Fig. 3 for all three
models.555Notice a difference in methodology between Figs. 4 and 5: in 4, we
fit the delta model to the $\tanh^{2}$ model separately for each
$\phi_{\text{amp}}$ and plot the best-fit results, while in 5, the fit is only
done at one amplitude, the delta and box potentials are inferred from this,
and then $\phi_{\text{amp}}$ is varied within the fixed potential. The delta
and box Floquet charts do not aim to match the $\tanh^{2}$ chart beyond the
marked horizontal fit line. We see that, in all models, only the first peak is
fully tachyonic. The secondary peaks always have $\omega_{k}^{2}>0$ on the
plateau at $\bar{\phi}=\phi_{\text{amp}}$. The secondary $\tanh^{2}$ bands get
wider as $\phi_{\text{amp}}$ grows, consistently with the results of [17],
while the delta bands are always narrow, and the box bands exhibit modulation
that varies with $\phi_{\text{amp}}$. The strength of the secondary bands is
tied to the properties of the dip in $\omega_{k}^{2}$ between the plateau and
the non-tachyonic region. Due to this dip, the modes in these bands can still
be partly tachyonic in the $\tanh^{2}$ model, leading to stronger peaks. The
number of partly tachyonic peaks in the limit $\phi_{\text{amp}}\gg\phi_{0}$
is

$N_{\rm th}\approx\frac{1}{2\sqrt{3}}e^{\phi_{\text{amp}}/\phi_{0}}\,,$ (4.5)

where we used the estimate $n\approx kT/\pi$ of the number of peaks below $k$
and that $k^{2}<-U^{\prime\prime}_{\mathrm{min}}=2/3$ at these peaks. For
instance, when $\phi_{\text{amp}}/\phi_{0}\lesssim 2$, only the first peak is
tachyonic, but for $\phi_{\text{amp}}=5\phi_{0}$, over $100$ peaks are partly
tachyonic. This is reflected in the Floquet chart in Fig. 5. As was noted in
[17], in the limit of high $\phi_{\text{amp}}$, the secondary peaks play a
major role in preheating and cannot be neglected. In comparison, there is no
dip in the delta model, so the peaks are narrower and weaker. In the box
model, the negative delta functions mimic the dip and produce stronger
instability bands. Nevertheless, the match to the $\tanh^{2}$ model is not
perfect. The detailed shape of the dip is essential for accurate modelling of
the secondary peaks.

When $\phi_{\text{amp}}$ is increased, the growth rate of modes within a
single period, i.e., $\mu_{k}T$ gets larger. For the delta and box models,
$\phi_{\text{amp}}$ is limited from above by $\phi_{2}$, where the potential
has a local maximum. In particular, the half-period $T$ diverges as
$\phi_{\text{amp}}\to\phi_{2}$. In the opposite limit $\phi_{\text{amp}}\to
0$, the instability bands become increasingly narrower and weaker for the
$\tanh^{2}$ potential. Beyond the leading band, the delta model is
qualitatively similar. In the delta model, in the $T\ll 1/\Gamma_{0}$ limit we
find a first band at $k\in(0,\pi/T)$ with height $\mu_{k}T=0.79$. Curiously,
the instability of the box model first increases as
$\phi_{\text{amp}}\to\phi_{1}=0.59\,\phi_{2}$ because $\tilde{\Lambda}$ from
(3.24) increases as $|\dot{\bar{\phi}}|(\phi_{1})$ decreases, but then shuts
off completely for $\phi_{\text{amp}}<\phi_{1}$, since the potential is
quadratic.

###  4.2 UV instability bands

Having discussed the first instability bands, let us turn our attention to the
instability bands at large $k$. First, since the large $k$ modes behave almost
adiabatically with $\dot{\omega}_{k}/\omega_{k}^{2}\ll 1$, we have

$\frac{1}{2}{\rm
tr}\,G\sim\cos(\alpha_{k}(T)),\qquad\alpha_{k}(t)\equiv\int_{0}^{t}{\rm
d}t^{\prime}\,\omega_{k}(t^{\prime})\approx kt\,,$ (4.6)

that is, it oscillates between $-1$ and $1$. This result is derived by
considering the independent solutions $u_{\pm}(t)\propto e^{\pm
i\alpha(t)}/\sqrt{\omega_{k}}$ to the mode equation and then applying Eq.
(2.8) to compute the monodromy matrix. Since, by Eq. (2.9), instability bands
$\mu_{k}>0$ are possible only when $|{\rm tr}\,G/2|>1$, the adiabatic limit is
stable, as expected. In other words, $\mu_{k}\to 0$ as $k\to\infty$.

According to the general theory of the Hill equation [55], there is an
infinite number of instability bands and $|{\rm tr}\,G|-2$ must decrease as
$k^{-2}$ or faster when $k\to\infty$, implying that $\mu_{k}$ should diminish
at least as fast as $k^{-1}$.666This follows from $\mu_{k}T\sim{\rm
acosh}|{\rm tr}\,G/2|\sim\sqrt{|{\rm tr}\,G|-2}$ as $|{\rm tr}\,G|\sim 2$.
However, the $\mu_{k}\sim k^{-1}$ scaling in the UV would be too slow to be
physical as it would lead to a catastrophic growth of the energy density of UV
perturbations. To demonstrate this, consider initial growth of perturbations
sourced by vacuum fluctuations, in which case the energy density grows as [17]

$\delta\dot{\rho}\approx j\equiv\int\frac{{\rm
d}^{3}k}{(2\pi)^{3}}2\mu_{k}\delta\rho^{\rm vac}_{k}\,,$ (4.7)

where $\delta\rho^{\rm vac}_{k}=k/2$ is the vacuum energy density. This
integral must be finite in well-behaved models. Otherwise, the linear
description cannot be valid up to arbitrarily high $k$, and there must exist a
cut-off scale $k_{\mathrm{cut}}$ regulating the growth of
perturbations.777This implies a breakdown of the perturbative vacuum at high
$k$, though this could be cured by non-linear effects. Approximating the
$n$-th peak in the spectrum of growth rates as

$\mu_{k}\approx\mu_{\text{peak},n}\left(1-4\frac{(k-k_{n})^{2}}{\Delta_{\text{peak},n}^{2}}\right),$
(4.8)

where $\mu_{\text{peak},n}$ $k_{n}$, and $\Delta_{\text{peak},n}$ give the
height, position and width of the peak, respectively, we can approximate

$j\sim\frac{1}{3\pi^{2}}\sum_{n}\mu_{\text{peak},n}\Delta_{\text{peak},n}k_{n}^{3}\,.$
(4.9)

As implied from extrema of ${\rm tr}\,G$ in Eq. (4.6), it approximately holds
that $k_{n}\propto n$. Thus, for the sum to converge,
$\mu_{\text{peak},n}\Delta_{\text{peak},n}$ must die away faster than
$k^{-4}$. In the following, we will show that this requirement is not
satisfied by the simplified models. Therefore, they must be equipped with a UV
cut-off on $\mu_{k}$. In practice, this is a benign issue, as the scales at
which the problematic UV behaviour shows up are well separated from the scale
of the tachyonic instability. Thus, since the fragmentation dynamics is
governed by the first instability bands, imposing a UV cut-off has a
negligible impact on the applicability of the simplified models beyond curing
the far-UV divergence.

Furthermore, although Hubble friction can be mostly neglected for the leading
bands, expansion may soften the UV instability by causing the growing high-$k$
modes to redshift out of the increasingly narrow instability bands. For any
amount of expansion, modes with a sufficiently high $k$ will not stay on a
given resonance band for even within a single oscillation of the background
field – periodicity of the mode equation cannot be assumed for such modes and
the Floquet analysis breaks down. Particles may still be produced at high $k$
in the presence of a strong feature such as the delta peaks, so UV divergences
may still appear, but they should be analysed using different techniques.
Since the $k$ value where such effects kick in varies from case to case, we
ignore this issue for the rest of the section. We proceed to study the UV
behaviour of our models in the flat space limit, believing it to be of
theoretical interest, if not directly relevant for the preheating dynamics.

In the UV, it is the small deviations from perfect adiabaticity that cause
instability. As these deviations approach zero when $k\to\infty$, the maxima
of $\mu_{k}$ must lie close to the extrema of $\cos(\alpha_{k}(T))$ from
(4.6), that is, at

$\alpha_{k_{n}}(T)=\pi n+\epsilon_{n}\,,\quad n\in\mathbb{Z}_{+},\qquad n\gg
1\,,$ (4.10)

where $\epsilon_{n}$ denotes small deviations from the adiabatic expectation.
Indeed, the peaks get increasingly narrower so the difference is important
especially in numerical searches, for which $\alpha_{k}(T)=\pi n$ also gives a
considerably better initial guess than the simpler limiting case $kT=\pi n$.
We look for the extrema of ${\rm tr}\,G/2$ of the form (4.10) by considering
the leading terms in the limits $n\gg 1$ and $\epsilon_{n}\ll 1$ and then
maximizing with respect to $\epsilon_{n}$. A general expansion around the
peaks then reads

$\quantity|\frac{1}{2}{\rm tr}\,G|\approx
1-\frac{1}{2}\epsilon_{n}^{2}+F_{0,n}+\epsilon_{n}F_{1,n}\,,$ (4.11)

where the first term arises from expanding the leading adiabatic contribution
(4.6) and the last two terms represent small non-adiabatic corrections to
zeroth and first orders in $\epsilon_{n}$. At the leading order, we can ignore
corrections of the form $\epsilon_{n}^{2}F_{2,n}$, since we expect $F_{2,n}\ll
1/2$. The maximum of (4.11) then lies at $\epsilon_{n}=F_{1,n}$, with $|{\rm
tr}\,G/2|=1+F_{0,n}+F_{1,n}^{2}/2$. Expanding ${\rm
acosh}(1+x)\approx\sqrt{2x}$ to obtain $\mu_{k}$ from (2.9) and mapping this
into the parameters of (4.8) we find that, at the leading order,

$\mu_{\text{peak},n}=\frac{1}{T}\sqrt{2F_{0,n}+F_{1,n}^{2}}\,,\qquad\Delta_{\text{peak},n}=2\mu_{\text{peak},n}\,,\qquad
k_{n}=\pi n/T+F_{1,n}\,.$ (4.12)

Thus the widths of the instability bands are diminished at the same pace as
their heights. In particular,
$j\propto\sum_{n}\mu_{\text{peak},n}^{2}k_{n}^{3}$, so it is sufficient to
require that $\mu_{\text{peak},n}$ decreases faster than $k^{-2}$ in the UV
for $j$ to converge.

Let us now give concrete examples of the UV instability bands in the
simplified models and the $\tanh^{2}$ model.

Figure 6: Positions of the peaks in $\mu_{k}$ in the large-$k$ tail. The
example models are the same as those in Fig. 3, except for the
$\tilde{\Lambda}=0$ box, which is obtained from the box model by adjusting
$T_{1}/T$ so that $\tilde{\Lambda}=0$, achieved with $T_{1}=0.070\,T$. The
grey lines give fits to the data: for the $\tanh^{2}$ model, the numerical fit
$\mu_{\text{peak}}T=3.63\,e^{-1.47k/m_{\rm th}}$; for the delta model, the
$k^{-1}$ scaling of (4.14); and for the $\tilde{\Lambda}=0$ box model, the
modulated $k^{-2}$ scaling of (4.17). The dashed black line gives the envelope
of the modulated peaks for the box from the magnitude of (4.17). Notice the
excellent matches between the peaks and the analytic estimates.

#### UV instability of the delta model

For the delta model, $\alpha(T)=-i\Gamma_{k}T$. Expanding the trace of the
monodromy matrix (3.6) around $\alpha(T)=\pi n$ then gives

$\quantity|\frac{1}{2}{\rm tr}\,G|\approx
1-\frac{1}{2}\epsilon_{n}^{2}-\frac{\Lambda}{2k}\epsilon_{n}\,,$ (4.13)

where $k$ corresponds to the value at $\epsilon_{n}=0$. Comparing with Eq.
(4.11), we find $F_{0,n}=0$ and $F_{1,n}=\Lambda/(2k)$, so that by Eq. (4.12)
we have

$\mu_{\text{peak}}=\frac{\Lambda}{2kT}+\mathcal{O}(k^{-2})\,.$ (4.14)

This behaviour is shown in Fig. 6. The peak decreases too slowly, so the $j$
integral (4.7) diverges. The source of this slow damping is the delta function
in the $\omega_{k}^{2}$ ansatz (3.3), which can efficiently excite arbitrarily
high $k$-bands through parametric resonance. In particular, the adiabaticity
condition $\dot{\omega}_{k}/\omega_{k}^{2}\ll 1$ is never completely satisfied
due to this delta function.

#### Softening the UV instability with the box model

For the box model, $\alpha(T)=\alpha_{1}+\alpha_{2}$, where $\alpha_{1}\equiv
T_{1}{\tilde{\omega}}_{k}\sim kT_{1}$, $\alpha_{2}\equiv T_{2}\Gamma_{k}\sim
kT_{2}$. At $k\gg 1$, the trace of the monodromy matrix (3.17) an be arranged
as

$\frac{1}{2}{\rm
tr}\,G=\cos(\alpha(T))+\frac{{\tilde{\Lambda}}}{k}\sin(\alpha(T))-\frac{(\Gamma_{0}^{2}+{\tilde{\omega}}^{2})^{2}}{8k^{4}}\sin(\alpha_{1})\sin(\alpha_{2})+\mathcal{O}(\Lambda^{2},k^{-6})\,,\\\
$ (4.15)

where we included the leading term that survives it the ${\tilde{\Lambda}}\to
0$ limit. Replacing $\alpha(T)=\pi n+\epsilon_{n}$ and expanding in
$\epsilon_{n}$ gives

$\quantity|\frac{1}{2}{\rm tr}\,G|\approx
1-\frac{1}{2}\epsilon_{n}^{2}+\frac{{\tilde{\Lambda}}}{k}\epsilon_{n}+\frac{(\Gamma_{0}^{2}+{\tilde{\omega}}^{2})^{2}}{8k^{4}}\sin^{2}(T_{1}k)\,,$
(4.16)

where the last two terms give us $F_{1,n}$ and $F_{0,n}$, respectively. If
${\tilde{\Lambda}}\neq 0$, then
$\mu_{\text{peak}}\sim|{\tilde{\Lambda}}|/(kT)$ resembling the UV-behaviour of
the delta model. However, when ${\tilde{\Lambda}}=0$, we find

$\mu_{\text{peak},{\tilde{\Lambda}}=0}=\frac{\Gamma_{0}^{2}+{\tilde{\omega}}^{2}}{2k^{2}T}\quantity|\sin(T_{1}k)|+\mathcal{O}(k^{-3})\,,$
(4.17)

so the $j$ integral (4.7) diverges only logarithmically. Notice that the
adiabaticity condition $\dot{\omega}_{k}/\omega_{k}^{2}\ll 1$ is not satisfied
at the transition between the tachyonic and non-tachyonic phases, i.e., at the
boundary of the box, even if ${\tilde{\Lambda}}=0$. Although this violation is
weaker, it allows the parametric resonance to stay active at very high $k$ and
provides a physical explanation for the divergence of $j$.

Another interesting feature is the modulation of
$\mu_{\text{peak},{\tilde{\Lambda}}=0}$ due to the secondary scale $T_{1}$,
also shown in Fig. 6. Moreover, (4.17) reproduces the delta model (4.14) in
the $T_{1}\to 0$ limit when we make the identification
$T_{1}{\tilde{\omega}}^{2}=\Lambda$. Notice that this matching differs
slightly from Eq. (3.19) since we set ${\tilde{\Lambda}}=0$. In particular,
the $k\to 0$ and $T_{1}\to 0$ limits of $\mu_{k}$ do not commute.

#### The $\tanh^{2}$ potential and UV stability of smooth potentials

For the $\tanh^{2}$ model (4.1), we must resolve the UV instability bands
numerically. A numerical scan of amplitudes in the range
$\phi_{\text{amp}}\in(\phi_{0},7\phi_{0})$ reveals an exponential
suppression888We expect this suppression to become effective in the fully non-
tachyonic region which, for high $\phi_{\text{amp}}$, means a high band number
$n$, as discussed around Eq. (4.5).

$\mu_{\text{peak}}\approx\Delta_{\text{peak}}/2\sim e^{-1.5k/m_{\rm th}}\,,$
(4.18)

with some percent level variation in the exponential slope lying within the
numerical uncertainty of the computation. Consistent with Eq. (4.12), the
widths obey roughly $\Delta_{\text{peak}}\approx 2\mu_{\text{peak}}$. The case
$\phi_{\text{amp}}=3\phi_{0}$ is shown in Fig. 6. Due to the exponential
damping of $\mu_{k}$, the $j$ integral (4.7) converges, and the model is free
of instabilities in the UV. The damping scale $m_{\rm th}$ can be understood
intuitively – it is the highest scale contributing to the time-dependence of
$U^{\prime\prime}$. Effectively, $m_{\rm th}$ acts as the cut-off scale of the
model—modes with $k\gg m_{\rm th}$ are increasingly adiabatic and stable.

Above we have observed that the degree of discontinuity in $\omega_{k}^{2}$
corresponds how fast $\mu_{k}$ decreases with $k$: if $\omega_{k}^{2}$
contains a $\delta$-function or a simple discontinuity, then the scalings
$\mu_{k}\sim k^{-1}$ and $\mu_{k}\sim k^{-2}$ are obtained, respectively. A
smooth $\omega_{k}^{2}$, on the other hand, leads to the exponential
suppression of $\mu_{k}$. One can further show that a continuous
$\omega_{k}^{2}$ leads to a decrease at least as fast as $\mu_{k}\sim k^{-3}$
[55], suggesting that in a model where the $n^{\mathrm{th}}$ derivative of
$\omega_{k}^{2}$ is the first discontinuous one, the heights and widths of the
peaks decrease as $k^{-n-2}$.999Although the UV-instability can be removed by
using a piecewise continuous ansatz for $\omega_{k}^{2}$, the resulting UV-
instability bands would still significantly differ from smooth physical
models. We further expect that $\mu_{k}$ is exponentially damped whenever
$\omega_{k}^{2}$ is smooth, though we are not aware of a general proof.

##  5 Multi-field scenarios

Figure 7: Comparison of mode growth of the inflation (orange) and orthogonal
(green) perturbations for $O(N)$ inflation assuming the $\tanh^{2}$ potential.
_Left panel:_ An example of the $\mu_{k}$ spectra for
$\phi_{\text{amp}}=3\phi_{0}$. _Upper right panel:_ The function
$\omega_{0}^{2}$ for $\phi_{\text{amp}}=3\phi_{0}$. Note that the orthogonal
component is never tachyonic. _Lower right panel:_ The Floquet chart for the
orthogonal modes. The dashed horizontal line indicates
$\phi_{\text{amp}}=3\phi_{0}$ depicted in the other panels.

Preheating can be more complicated when the inflaton is coupled to other
fields since these too can be excited due to the oscillating background. It is
also possible that the inflaton itself is not a simple real scalar. This
section focuses on the latter possibility and considers an inflaton symmetric
under a global $O(N)$ group as a simple multi-field example of tachyonic
preheating.

As above, we assume fast preheating and neglect the effect of expansion. The
$O(N)$ symmetry dictates that the potential $V(\phi)$ is a function of the
field modulus $\phi\equiv\sqrt{\sum_{i}\phi^{2}_{i}}$ only and the background
field equations read

$\ddot{\bar{\phi}}_{i}=-\hat{\phi}_{i}V^{\prime},$ (5.1)

where $\hat{\phi}_{i}\equiv\bar{\phi}_{i}/\bar{\phi}$, the prime denotes a
derivative with respect to $\phi$, and a bar refers to the background
field.101010In contrast to the rest of the article, in this section
$\bar{\phi}$ is the modulus of the background field and thus always positive.
We assume that the potential has a minimum at $\phi=0$ and that it is a
monotonously growing function of $\phi$, so that $V^{\prime}>0$. Expansion
during inflation damps all angular motion in the field space, so the
background field rolls down in a radial direction. We choose the basis so that
the background is non-vanishing only along the $\phi_{1}$ component,
$\bar{\phi}_{i}\propto\delta_{1i}=\hat{\phi}_{i}$, which we will refer to as
the inflaton direction. The perturbations
$\delta\phi\equiv\phi_{i}-\bar{\phi}_{i}$ evolve as

$\delta\ddot{\phi}_{i,k}+(m^{2}_{ij}+k^{2}\delta_{ij})\delta\phi_{j,k}=0\,,\qquad\mbox{where}\quad
m^{2}_{ij}\equiv\partial_{\phi_{j},\phi_{i}}V\,.$ (5.2)

For $O(N)$ symmetric fields, the mass matrix

$m^{2}_{ij}=\hat{\phi}_{i}\hat{\phi}_{i}V^{\prime\prime}+\left(\delta_{ij}-\hat{\phi}_{i}\hat{\phi}_{i}\right)V^{\prime}/\bar{\phi}={\rm
diag}(V^{\prime\prime},V^{\prime}/\bar{\phi},V^{\prime}/\bar{\phi},\ldots)$
(5.3)

is diagonal and the directions decouple from each other. Moreover, only the
inflaton direction possesses a tachyonic instability through
$V^{\prime\prime}<0$—for monotonous potentials, $V^{\prime}>0$, so the mass
squared of orthogonal perturbations is positive. These directions are still
subject to parametric instabilities. The corresponding mode equations are

$\delta\ddot{\phi}_{\perp,k}+(V^{\prime}/\bar{\phi}+k^{2})\delta\phi_{\perp,k}=0\,,$
(5.4)

leading to some growth rates $\mu_{\perp,k}$. The symbol $\perp$ refers to an
orthogonal component with $i>1$. These components are identical to each other.
These modes have the following general features:

  * •

The $k=0$ mode is not growing. This can be seen by noting that
$V^{\prime}/\bar{\phi}=-\ddot{\bar{\phi}}_{1}/\bar{\phi}_{1}$ so that
$\delta\phi_{\perp,0}(t)\propto\bar{\phi}_{1}(t)$. Since the background is
periodic in $\bar{\phi}(t)$, then, analogously to the inflaton mode, its
Floquet exponents must vanish [17, 55].

  * •

The first instability band of orthogonal modes must begin at some finite, non-
zero $k$ (for a derivation, see appendix A). This is to be contrasted with the
perturbations in the tachyonic direction for which the first instability band
begins already at $k=0$.111111This holds when the oscillation period decreases
with decreasing energy density (see Eq. (2.10)). As the band structure of both
the inflaton and the orthogonal modes have a similar period in $k$-space, the
first orthogonal modes must then be narrower. An example of these features is
shown in Fig. 7.

  * •

In the limit $k\to\infty$, the perturbations are adiabatic and thus the growth
rates vanish. How fast $\mu_{\perp,k}$ approaches 0 depends on
$V^{\prime}/\bar{\phi}$. If $V^{\prime}/\bar{\phi}$ is smooth, as for the
$\tanh^{2}$ model, it’s expected that $\mu_{\perp,k}$ is damped exponentially
for large $k$ as was the case with the inflaton mode. By Eq. (4.6), the
corresponding monodromy matrix oscillates as $\cos(kT)$ in the $k\to\infty$
limit.

One could construct simplified models for the orthogonal modes, for instance,
by relying on the generating potentials (3.11) and (3.22). However, the
orthogonal mode equation (5.2) does not have a simple analytic solution for
these potentials. Thus we found it more illuminating to rely on general
analytic derivations and numerical methods. As an example, we worked out the
case of the $\tanh^{2}$ potential (4.1) shown in Fig. 7.

In all, the orthogonal modes are subject to a narrow parametric resonance that
is generally weaker than the broad tachyonic instability of inflaton modes.
Thus, at least at the leading order, it is sufficient to consider only the
inflaton modes when studying inflaton fragmentation in the linear regime.

##  6 Discussion

In a large class of inflationary models with an exponentially flat plateau,
including the $\tanh^{2}$ model, tachyonic preheating is efficient only when
the tensor-to-scalar ratio is extremely small, $r<10^{-7}$ [17]. Preheating is
also a potential source of gravitational waves (GWs) [58, 21, 22, 59, 60, 61,
62, 63, 64, 53, 65, 66], and for tachyonic preheating in these potentials, the
GW signal is expected to peak at frequency $f_{\mathrm{GW}}=10^{9}$ Hz or
above [17]. Such predictions are useful for differentiating models of
inflation from each other, but in this case, they also lie well beyond our
current observational capabilities. It is then interesting to ask whether
these predictions are generic to all inflationary models that transition into
an epoch of rapid tachyonic preheating. We show that this is not the case, but
instead, it is possible to construct inflationary models with fast tachyonic
preheating that can produce virtually any $r$ and $f_{\mathrm{GW}}$.

##### Inflationary observables.

Let us first review how the limits on the $\tanh^{2}$ preheating arise. As
discussed in section 4, the CMB scalar perturbation strength $A_{s}$ fixes the
mass parameter $m_{\rm th}^{2}=U_{0}/\phi_{0}^{2}\approx 5\times 10^{-6}$,
which sets also the scale of tachyonicity during preheating,
$|U^{\prime\prime}|\sim m_{\rm th}^{2}$. Thus, in this case, the shape of the
potential relates the scales relevant to CMB and preheating. Fast tachyonic
preheating means that the background oscillation time $T$ during preheating
must satisfy $T\ll H^{-1}$. Since $H\sim\sqrt{U_{0}}$ from the first Friedmann
equation and $T\sim m_{\rm th}^{-1}\phi_{0}^{-1/2}$ [17], we get the condition
$\phi_{0}\ll 1$. A more precise analysis shows that $\phi_{0}\lesssim 0.01$ is
necessary and sufficient for exponentially flat potentials [17]. As $m_{\rm
th}$ is fixed, we obtain that $U_{0}\lesssim(5\times 10^{14}\,\text{GeV})^{4}$
and $r\lesssim 8\times 10^{-8}$. To put it briefly, an effective tachyonic
preheating requires that the scale of inflation must be lower than the mass.
On the other hand, if the mass is fixed as for the $\tanh^{2}$ potential, it
may infer a small scale of inflation and thus a small $r$.

Let us consider a more general model, where the CMB and preheating scales are
not intimately related. For instance, as a proof of concept, we may construct
a potential by glueing good inflationary plateaus to the delta model
potential. The shape of the plateau and the minimum can then be tuned
independently. The condition (2.4) for tachyonic preheating, $T\ll H^{-1}$,
reads $\phi_{2}\ll 1$, which simply means that _efficient preheating requires
sub-Planckian field excursions_. We should further impose that the mass scales
are sub-Planckian, $\Gamma_{0}\ll 1$, which gives the hierarchy

$H\ll\phi_{2}\ll 1\,.$ (6.1)

At this very general level of discussion, this is the strongest theoretical
bound and leaves plenty of room for the current experimental bounds $H<2\times
10^{-5}$ [8] to be satisfied121212The bound on $H$ arises from $r<0.036$ [57]
in single field slow-roll inflation.. The smaller the field scale $\phi_{2}$,
the more difficult it will be for $r$ (or $H$) to saturate the observational
bounds. Although this argument does not give an explicit physically viable
inflationary model, it provides a proof of concept showing that tachyonic
preheating can be consistent with any experimentally allowed set of
inflationary parameters.

##### Induced gravitational waves.

Let us now consider GWs. Their frequency can be estimated as
$f_{\mathrm{GW}}\approx(zT)^{-1}$, where the redshift is $z\approx T_{\rm
reh}/T_{\mathrm{CMB}}$. Here $T_{\rm reh}$ is the reheating temperature, we
assumed instant inflaton fragmentation into radiation, the number of effective
degrees of freedom is $g_{*}\approx 100$, and $T_{\mathrm{CMB}}\approx 2.7$ K
is the present CMB temperature. First, in the $\tanh^{2}$ and similar models,
the period scales as $T\sim m_{\rm th}^{-1}\phi_{0}^{-1/2}$ (assuming
$\phi_{0}\gtrsim 10^{-4}$). Then, with $U_{\text{preh}}\approx U_{0}$ we
obtain that the frequency $f_{\mathrm{GW}}\approx T_{\mathrm{CMB}}\,m_{\rm
th}^{1/2}\approx 8\times 10^{8}$ Hz, which is out of the range of GW
interferometers. Since this only depends on $m_{\rm th}$, the CMB constraint
on $A_{s}$ fixes $f_{\mathrm{GW}}$.

In contrast, the tachyonicity scale is not constrained in the delta model, and
we can decrease $f_{\mathrm{GW}}\sim T^{-1}\sim\Gamma_{0}$ by making
$\Gamma_{0}$ small. Since we require $T\ll H^{-1}$, lowering $H$ (and $r$)
will lower the bound on $f_{\mathrm{GW}}$ and vice versa. In detail, we find
the lower bound

$f_{\mathrm{GW}}\gg 5\times 10^{-7}{\rm Hz}\times\left[\frac{T_{\rm reh}}{{\rm
GeV}}\right]\,.$ (6.2)

For example, for $T_{\rm reh}=1\,\text{TeV}$, we find $f_{\mathrm{GW}}\gg 0.5$
mHz, while at the current observational upper limit $T_{\rm reh}\approx
5\times 10^{15}{\rm GeV}$, we find $f_{\mathrm{GW}}\gg 2\times 10^{9}$ Hz.
Therefore, tachyonic preheating can generate GWs within the frequency range of
LIGO-Virgo-Kagra [67, 68, 69] and future GW interferometers such as ET [70]
and LISA [71] if the scale of inflation is sufficiently low.

##### Broad vs narrow resonances.

Finally, let us compare our results to prior literature. It is worth noting
that the band structure in tachyonic preheating is quite different from the
more commonly studied parametric resonance. Early on [9, 10], parametric
resonance models were divided into two categories: those with narrow and broad
resonance. In narrow resonance, the mode functions evolve almost
adiabatically, and the resonance consists of multiple weak, narrow bands
located roughly at integer values of $kT/\pi$ that do not extend to $k=0$.
This resembles the behaviour of the orthogonal fields in our multi-field setup
in section 5. In the broad resonance regime, the resonance is dominated by one
dominant band starting from $k=0$ and extending to high values of $kT$. Broad
resonance is usually observed in multi-field setups for fields coupled to the
inflaton. Our tachyonic resonance shares features from both: the dominant band
always starts at $k=0$ and reaches up to $k=2\pi/T$. The subsequent secondary
bands appear periodically in $kT/\pi$, and the lowest ones, still partly
tachyonic, are broad and strong and can significantly affect the preheating
dynamics. Bands in the UV tail are weak and narrow.

In [72, 48, 49], the authors discussed a ‘flapping resonance’ similar to ours,
where the inflaton repeatedly passes through tachyonic regions. Using an
effective parameter ${\tilde{q}}$, they classified the resulting spectra to
‘broad’ (${\tilde{q}}\ll 1$), ‘intermediate’ (${\tilde{q}}\simeq 1$), and
‘narrow’ (${\tilde{q}}\ll 1$) [49], mimicking the classification discussed in
the previous paragraph. In particular, they found that tachyonicity produces
an intermediate value of ${\tilde{q}}$. However, the value of ${\tilde{q}}$
approaches zero for our $\tanh^{2}$ model in the limit of high
$\phi_{\text{amp}}$, for instance, ${\tilde{q}}=0.034$ for
$\phi_{\text{amp}}=7\phi_{0}$, even though these models do not have narrow
resonances. At the same time, for the delta model, formally
${\tilde{q}}=\infty$, even though the resonance is not really of the broad
form. This shows that the classification scheme of [49] is not universal. The
difference between our results arises from the length of the tachyonic period:
in the models considered in [49], the tachyonic region was relatively short,
whereas, in our models, the modes in the leading instability band are
tachyonic for most of the time.

##  7 Conclusions

In this paper, we studied the linear regime of tachyonic preheating. The
considered scenarios take place if the inflaton repeatedly returns to a
tachyonic region of its potential, e.g. , a plateau during post-inflationary
oscillations. The presence of a tachyonic instability leads to a rapid
fragmentation of the coherent background field within a fraction of an
$e$-fold. Complementing our recent numerical study [17], we constructed
simplified models by postulating the time-dependence of the effective mass of
inflation perturbations. With these simplifications, the Floquet exponents
determining the growth of each mode can be found analytically. As an applied
example, we use these models to study tachyonic preheating for a physically
well-motivated $\tanh^{2}$ potential.

As for general inflationary potentials leading to tachyonic preheating, the
analytic models assume that the effective mass alternates periodically between
tachyonic and non-tachyonic phases. We considered the following cases:

  * •

In the minimal simplified model— _the delta model_ —we neglected the duration
of the non-tachyonic phase and described the effective mass squared by a delta
function and a negative constant. Mode growth is characterized by two
parameters: the period of oscillations and the constant negative mass squared.
We showed that this ansatz is equivalent to assuming the idealized two-
parameter potential $-m^{2}(|\phi|-\phi_{0})^{2}/2$.

  * •

As a generalization of the delta model, we considered _the box model_ which
also accounts for the duration of the tachyonic period. In this set-up, mode
growth is determined by four parameters: the constant effective masses and the
durations of the tachyonic and non-tachyonic epochs. Motivated by numerical
observations of the oscillating behaviour of realistic potentials, we further
allow for a drop in the effective mass squared during transitions between the
tachyonic and non-tachyonic phases, which we model by a delta function. As in
the delta model, the postulated time dependence of the effective mass can be
realized by an idealized piecewise quadratic continuous potential. The box
model reduces to the delta model in the limit when the duration of the non-
tachyonic period approaches zero.

In both models, the ansätze for the time-dependence of the effective mass
squared contains an extra parameter related to the others by the requirement
that the zero momentum mode must not grow.

Comparing the idealized analytic spectra of growth rates $\mu_{k}$ to
numerical $\tanh^{2}$ results, we find that both models capture well the
leading tachyonic peak in $\mu_{k}$. This justifies the usefulness of this
simplification, as, in tachyonic preheating, the overall growth of the energy
density, as well as the spectrum of the fragmented field component, is
determined chiefly by the first instability band. The delta model, in
particular, is useful for estimating the rate of tachyonic preheating since it
provides an excellent analytic model of the first peak. The subleading peaks
are better approximated by the box model. All in all, the simplified models
provide a good analytical understanding of the band structure in tachyonic
preheating. The most significant discrepancy with realistic models was found
in the UV, where, for smooth potentials, the instability bands get
exponentially narrower and weaker, while, for the simplified models, we
observe only a power-law reduction with growing wavenumber. This discrepancy
can be traced back to the discontinuities introduced by hand. This observation
suggests that one can construct simplified models with improved UV behaviour
by using ansätze with continuous higher time derivatives.

We also considered preheating in multi-field scenarios in the example of an
$O(N)$ inflaton. In such scenarios, one direction of the field space
corresponds to the inflaton degree of freedom while the rest behave as
spectator fields. They are excited during preheating due to their coupling to
the inflaton. We showed that only the inflaton is subject to a tachyonic
instability as long as the potential is monotonous in $|\phi|$. The
perturbations of the spectator components grow due to a less effective
parametric resonance. Therefore, in the linear regime, the inflaton fragments
into itself and the expected equipartitioning of the energy density among all
components of the $O(N)$ multiplet completes during the subsequent non-linear
evolution.

Based on our results, we make the following general observations on tachyonic
preheating:

  * •

In typical tachyonic inflationary potentials, tachyonic preheating is
effective if the field excursions are sub-Planckian.

  * •

Tachyonic preheating is expected to complete within less than an $e$-fold in
all cases. Comparing the characteristic timescales $1/H=\sqrt{3/\rho}$,
$1/\mu_{\text{peak}}$ and $T$, we find that $1/H\gg 1/\mu_{\text{peak}}>T$ in
the tachyonic regime—damping inflaton oscillation amplitude is dominated by
fragmentation/instability and the inflaton will oscillate at least a few times
before fragmentation.

  * •

Effective tachyonic preheating can be possible for any set of experimentally
allowed inflationary parameters. In particular, it does not imply an extremely
small tensor-to-scalar ratio $r$.

  * •

The frequency of GWs produced during tachyonic preheating can range from nHz
to GHz and may thus be accessible to future GW experiments.

These results can be used as guidelines for future model building with
tachyonic preheating.

## Acknowledgments

This work was supported by the Estonian Research Council grants PRG803,
PRG1055, MOBTP135, MOBJD381 and MOBTT5 and by the EU through the European
Regional Development Fund CoE program TK133 “The Dark Side of the Universe.”

##  Appendix A Infrared stability of O(N) field perturbations

In this appendix, we will consider the infrared behaviour of the Floquet
exponents in the scenario in which the inflaton belongs to a multiplet with an
$O(N)$ global symmetry. At $k=0$, the equations for the inflaton mode (2.5)
and the orthogonal modes (5.4) can be recast as
$\delta\ddot{\phi}_{1,0}-(\dddot{\bar{\phi}}_{1}/\dot{\bar{\phi}}_{1})\delta\ddot{\phi}_{1,0}=0$
and
$\delta\ddot{\phi}_{\perp,0}-(\ddot{\bar{\phi}}_{1}/\bar{\phi}_{1})\delta\ddot{\phi}_{\perp,0}=0$,
respectively. Their two independent solutions read

$u_{1}=\left\\{\begin{array}[]{ll}\dot{\bar{\phi}}_{1}&\qquad\mbox{(inflaton)}\\\
\bar{\phi}_{1}&\qquad\mbox{(orthogonal)}\end{array}\right.,\qquad
u_{2}=u_{1}\int\frac{{\rm d}t}{u_{1}^{2}}\,.$ (A.1)

These solutions have a unit Wronskian.

In both cases, we can choose the point $t=0$ so that $\dot{u}_{1}(0)=0$
without of loss of generality: for the inflaton modes, we choose $t=0$ when
the inflaton crosses the minimum of the potential, so that
$\dot{u}_{1}(0)=\ddot{\bar{\phi}}_{1}(0)=0$, and for the orthogonal modes,
$t=0$ corresponds to a turning point since
$\dot{u}_{1}(0)=\dot{\bar{\phi}}_{1}(0)=0$. In both cases, the $u_{2}$
integrand in (A.1) is singular at $T/2$, since $u_{1}(T/2)=0$. For the
inflaton modes, this corresponds to a turning point with
$u_{1}(T/2)=\dot{\bar{\phi}}_{1}(T/2)=0$, while for orthogonal modes we have a
zero-crossing, i.e., $u_{1}(T/2)=\bar{\phi}_{1}(T/2)=0$. Notice that here we
use the assumption that the potential is symmetric—otherwise the field at the
half period may not lie at the origin.

To address this singularity we follow the derivation in Ref. [17] and
construct the $u_{2}$ solution as [17]

$u_{2}(t)=\left\\{\begin{array}[]{lc}\tilde{u}_{2}(t)&,\quad 0\leq t\leq
T/2\\\ \tilde{u}_{2}(T-t)+\beta u_{1}(t)&,\quad T/2\leq t\leq
T\end{array}\right.\,,$ (A.2)

with $\tilde{u}_{2}(t)$ defined by (A.1) when $0\leq t\leq T/2$. In this
construction, $u_{2}$ is continuous while the continuity of $\dot{u}_{2}$
gives

$\displaystyle\beta$
$\displaystyle=\frac{2\dot{\tilde{u}}_{2}(T/2)}{\dot{u}_{1}(t)}=\lim_{t\to
T/2}\frac{2}{\dot{u}_{1}(T/2)}\left[\dot{u}_{1}(t)\int^{t}_{0}\frac{{\rm
d}t^{\prime}}{u_{1}(t^{\prime})^{2}}+\frac{1}{u_{1}(t)}\right]$ (A.3)

so that the monodromy matrix around $k=0$ reads

$\displaystyle\frac{1}{2}\Tr G$
$\displaystyle=-1+\frac{k^{2}}{2}\beta\int^{T}_{0}{\rm
d}t^{\prime}\,u_{1}(t^{\prime})^{2}+\mathcal{O}(k^{4}).$ (A.4)

The derivation is similar for asymmetric potentials, although slightly more
tedious. In the latter case, one also has $\Tr G/2=1$ instead of $\Tr G/2=-1$
at $k=0$.

Plugging in $u_{1}$ from Eq. (A.1) and eliminating the derivatives using the
equations of motion and energy conservation for the background field, we find
that

$\frac{1}{2}\Tr
G=-1-\frac{k^{2}}{2}W\partial_{\bar{\rho}}T+\mathcal{O}(k^{4})$ (A.5)

for the inflaton modes. The abbreviated action $W$ was defined in Eq. (2.11).
For the orthogonal modes

$\displaystyle\frac{1}{2}\Tr G_{\perp}$
$\displaystyle=-1+k^{2}\frac{\beta_{\perp}}{\sqrt{2\bar{\rho}}}\int^{\phi_{\rm
amp}}_{0}\frac{{\rm
d}\phi\,\phi^{2}}{\sqrt{1-V(\phi)/\bar{\rho}}}+\mathcal{O}(k^{4})\,,$ (A.6)

where $\bar{\rho}=V(\phi_{\rm amp})$ and $\beta_{\perp}$ does not have a short
analytical expression. However, since $\beta_{\perp}\geq 0$, we have that
$|\Tr G_{\perp}|\leq 1$ around $k=0$ and thus there must exist a $k_{1}$ such
that $k\in[0,k_{1}]$ is stable. In other words, the $\mu_{k}$ spectrum for
orthogonal perturbations begins with a stability band. The instability bands
of orthogonal modes thus resemble a narrow parametric resonance.

## References

  * [1] A. A. Starobinsky, A New Type of Isotropic Cosmological Models Without Singularity, Phys. Lett. B 91 (1980) 99–102.
  * [2] A. A. Starobinsky, The Perturbation Spectrum Evolving from a Nonsingular Initially De-Sitter Cosmology and the Microwave Background Anisotropy, Sov. Astron. Lett. 9 (1983) 302.
  * [3] A. H. Guth, The Inflationary Universe: A Possible Solution to the Horizon and Flatness Problems, Phys. Rev. D 23 (1981) 347–356.
  * [4] A. D. Linde, A New Inflationary Universe Scenario: A Possible Solution of the Horizon, Flatness, Homogeneity, Isotropy and Primordial Monopole Problems, Phys. Lett. B 108 (1982) 389–393.
  * [5] A. Albrecht and P. J. Steinhardt, Cosmology for Grand Unified Theories with Radiatively Induced Symmetry Breaking, Phys. Rev. Lett. 48 (1982) 1220–1223.
  * [6] A. D. Linde, Chaotic Inflation, Phys. Lett. B 129 (1983) 177–181.
  * [7] D. H. Lyth and A. Riotto, Particle physics models of inflation and the cosmological density perturbation, Phys. Rept. 314 (1999) 1–146, [hep-ph/9807278].
  * [8] Planck Collaboration, Y. Akrami et al., Planck 2018 results. X. Constraints on inflation, Astron. Astrophys. 641 (2020) A10, [arXiv:1807.06211].
  * [9] L. Kofman, A. D. Linde, and A. A. Starobinsky, Reheating after inflation, Phys. Rev. Lett. 73 (1994) 3195–3198, [hep-th/9405187].
  * [10] L. Kofman, A. D. Linde, and A. A. Starobinsky, Towards the theory of reheating after inflation, Phys. Rev. D 56 (1997) 3258–3295, [hep-ph/9704452].
  * [11] A. D. Dolgov and D. P. Kirilova, ON PARTICLE CREATION BY A TIME DEPENDENT SCALAR FIELD, Sov. J. Nucl. Phys. 51 (1990) 172–177.
  * [12] J. H. Traschen and R. H. Brandenberger, Particle Production During Out-of-equilibrium Phase Transitions, Phys. Rev. D 42 (1990) 2491–2504.
  * [13] Y. Shtanov, J. H. Traschen, and R. H. Brandenberger, Universe reheating after inflation, Phys. Rev. D 51 (1995) 5438–5455, [hep-ph/9407247].
  * [14] P. B. Greene, L. Kofman, A. D. Linde, and A. A. Starobinsky, Structure of resonance in preheating after inflation, Phys. Rev. D 56 (1997) 6175–6192, [hep-ph/9705347].
  * [15] G. N. Felder, J. Garcia-Bellido, P. B. Greene, L. Kofman, A. D. Linde, and I. Tkachev, Dynamics of symmetry breaking and tachyonic preheating, Phys. Rev. Lett. 87 (2001) 011601, [hep-ph/0012142].
  * [16] G. N. Felder, L. Kofman, and A. D. Linde, Tachyonic instability and dynamics of spontaneous symmetry breaking, Phys. Rev. D 64 (2001) 123517, [hep-th/0106179].
  * [17] E. Tomberg and H. Veermäe, Tachyonic Preheating in Plateau Inflation, arXiv:2108.10767.
  * [18] J. Garcia-Bellido and E. Ruiz Morales, Particle production from symmetry breaking after inflation, Phys. Lett. B 536 (2002) 193–202, [hep-ph/0109230].
  * [19] E. J. Copeland, S. Pascoli, and A. Rajantie, Dynamics of tachyonic preheating after hybrid inflation, Phys. Rev. D 65 (2002) 103517, [hep-ph/0202031].
  * [20] N. Barnaby and J. M. Cline, Nongaussian and nonscale-invariant perturbations from tachyonic preheating in hybrid inflation, Phys. Rev. D 73 (2006) 106012, [astro-ph/0601481].
  * [21] J. Garcia-Bellido, D. G. Figueroa, and A. Sastre, A Gravitational Wave Background from Reheating after Hybrid Inflation, Phys. Rev. D 77 (2008) 043517, [arXiv:0707.0839].
  * [22] J.-F. Dufaux, G. Felder, L. Kofman, and O. Navros, Gravity Waves from Tachyonic Preheating after Hybrid Inflation, JCAP 03 (2009) 001, [arXiv:0812.2917].
  * [23] M. Desroche, G. N. Felder, J. M. Kratochvil, and A. D. Linde, Preheating in new inflation, Phys. Rev. D 71 (2005) 103516, [hep-th/0501080].
  * [24] J. F. Dufaux, G. N. Felder, L. Kofman, M. Peloso, and D. Podolsky, Preheating with trilinear interactions: Tachyonic resonance, JCAP 07 (2006) 006, [hep-ph/0602144].
  * [25] S. Antusch, D. Nolde, and S. Orani, Hill crossing during preheating after hilltop inflation, JCAP 06 (2015) 009, [arXiv:1503.06075].
  * [26] S. Antusch, F. Cefala, D. Nolde, and S. Orani, Parametric resonance after hilltop inflation caused by an inhomogeneous inflaton field, JCAP 02 (2016) 044, [arXiv:1510.04856].
  * [27] Y. Ema, K. Mukaida, and K. Nakayama, Electroweak Vacuum Metastability and Low-scale Inflation, JCAP 12 (2017) 030, [arXiv:1706.08920].
  * [28] S. Antusch, F. Cefala, and S. Orani, What can we learn from the stochastic gravitational wave background produced by oscillons?, JCAP 03 (2018) 032, [arXiv:1712.03231].
  * [29] P. Brax, J.-F. Dufaux, and S. Mariadassou, Preheating after Small-Field Inflation, Phys. Rev. D 83 (2011) 103510, [arXiv:1012.4656].
  * [30] J. Rubio and E. S. Tomberg, Preheating in Palatini Higgs inflation, JCAP 04 (2019) 021, [arXiv:1902.10148].
  * [31] A. Karam, M. Raidal, and E. Tomberg, Gravitational dark matter production in Palatini preheating, JCAP 03 (2021) 064, [arXiv:2007.03484].
  * [32] A. Karam, E. Tomberg, and H. Veermäe, Tachyonic preheating in Palatini R 2 inflation, JCAP 06 (2021) 023, [arXiv:2102.02712].
  * [33] J. Martin, C. Ringeval, R. Trotta, and V. Vennin, The Best Inflationary Models After Planck, JCAP 03 (2014) 039, [arXiv:1312.3529].
  * [34] S. Ferrara, R. Kallosh, A. Linde, and M. Porrati, Minimal Supergravity Models of Inflation, Phys. Rev. D 88 (2013), no. 8 085038, [arXiv:1307.7696].
  * [35] R. Kallosh, A. Linde, and D. Roest, Superconformal Inflationary $\alpha$-Attractors, JHEP 11 (2013) 198, [arXiv:1311.0472].
  * [36] J. J. M. Carrasco, R. Kallosh, and A. Linde, Cosmological Attractors and Initial Conditions for Inflation, Phys. Rev. D 92 (2015), no. 6 063519, [arXiv:1506.00936].
  * [37] M. Galante, R. Kallosh, A. Linde, and D. Roest, Unity of Cosmological Inflation Attractors, Phys. Rev. Lett. 114 (2015), no. 14 141302, [arXiv:1412.3797].
  * [38] R. Kallosh and A. Linde, Universality Class in Conformal Inflation, JCAP 07 (2013) 002, [arXiv:1306.5220].
  * [39] T. Krajewski, K. Turzyński, and M. Wieczorek, On preheating in $\alpha$-attractor models of inflation, Eur. Phys. J. C 79 (2019), no. 8 654, [arXiv:1801.01786].
  * [40] F. Bauer and D. A. Demir, Inflation with Non-Minimal Coupling: Metric versus Palatini Formulations, Phys. Lett. B 665 (2008) 222–226, [arXiv:0803.2664].
  * [41] J. Rubio, Higgs inflation, Front. Astron. Space Sci. 5 (2019) 50, [arXiv:1807.02376].
  * [42] V.-M. Enckell, K. Enqvist, S. Rasanen, and L.-P. Wahlman, Inflation with $R^{2}$ term in the Palatini formalism, JCAP 02 (2019) 022, [arXiv:1810.05536].
  * [43] I. Antoniadis, A. Karam, A. Lykkas, and K. Tamvakis, Palatini inflation in models with an $R^{2}$ term, JCAP 11 (2018) 028, [arXiv:1810.10418].
  * [44] T. Tenkanen, Tracing the high energy theory of gravity: an introduction to Palatini inflation, Gen. Rel. Grav. 52 (2020), no. 4 33, [arXiv:2001.10135].
  * [45] D. S. Salopek, J. R. Bond, and J. M. Bardeen, Designing Density Fluctuation Spectra in Inflation, Phys. Rev. D 40 (1989) 1753.
  * [46] D. I. Kaiser, Primordial spectral indices from generalized Einstein theories, Phys. Rev. D 52 (1995) 4295–4306, [astro-ph/9408044].
  * [47] F. L. Bezrukov and M. Shaposhnikov, The Standard Model Higgs boson as the inflaton, Phys. Lett. B 659 (2008) 703–706, [arXiv:0710.3755].
  * [48] N. Kitajima, J. Soda, and Y. Urakawa, Gravitational wave forest from string axiverse, JCAP 10 (2018) 008, [arXiv:1807.07037].
  * [49] H. Fukunaga, N. Kitajima, and Y. Urakawa, Efficient self-resonance instability from axions, JCAP 06 (2019) 055, [arXiv:1903.02119].
  * [50] R. Kronig and W. G. Penney, Quantum mechanics of electrons in crystal lattices, Proc. R. Soc. Lond. A 130 (1931) 499–513.
  * [51] K. D. Lozanov and M. A. Amin, Equation of State and Duration to Radiation Domination after Inflation, Phys. Rev. Lett. 119 (2017), no. 6 061301, [arXiv:1608.01213].
  * [52] K. D. Lozanov and M. A. Amin, Self-resonance after inflation: oscillons, transients and radiation domination, Phys. Rev. D 97 (2018), no. 2 023533, [arXiv:1710.06851].
  * [53] K. D. Lozanov and M. A. Amin, Gravitational perturbations from oscillons and transients after inflation, Phys. Rev. D 99 (2019), no. 12 123504, [arXiv:1902.06736].
  * [54] A. Bhoonah, J. Bramante, S. Nerval, and N. Song, Gravitational Waves From Dark Sectors, Oscillating Inflatons, and Mass Boosted Dark Matter, JCAP 04 (2021) 043, [arXiv:2008.12306].
  * [55] R. Reissig, Magnus, W. und S. Winkler, Hill’s Equation (Interscience Tracts in Pure and Applied Mathematics, Number 20). VIII + 127 S. New York/London/Sydney 1966. Interscience Publishers. Preis geb. 68 s.net, Zeitschrift Angewandte Mathematik und Mechanik 48 (Jan., 1968) 138–139.
  * [56] M. S. Turner, Coherent Scalar Field Oscillations in an Expanding Universe, Phys. Rev. D 28 (1983) 1243.
  * [57] BICEP, Keck Collaboration, P. A. R. Ade et al., Improved Constraints on Primordial Gravitational Waves using Planck, WMAP, and BICEP/Keck Observations through the 2018 Observing Season, Phys. Rev. Lett. 127 (2021), no. 15 151301, [arXiv:2110.00483].
  * [58] J. Garcia-Bellido and D. G. Figueroa, A stochastic background of gravitational waves from hybrid preheating, Phys. Rev. Lett. 98 (2007) 061302, [astro-ph/0701014].
  * [59] J.-F. Dufaux, D. G. Figueroa, and J. Garcia-Bellido, Gravitational Waves from Abelian Gauge Fields and Cosmic Strings at Preheating, Phys. Rev. D 82 (2010) 083518, [arXiv:1006.0217].
  * [60] S. Antusch, F. Cefala, and S. Orani, Gravitational waves from oscillons after inflation, Phys. Rev. Lett. 118 (2017), no. 1 011303, [arXiv:1607.01314]. [Erratum: Phys.Rev.Lett. 120, 219901 (2018)].
  * [61] J. Liu, Z.-K. Guo, R.-G. Cai, and G. Shiu, Gravitational Waves from Oscillons with Cuspy Potentials, Phys. Rev. Lett. 120 (2018), no. 3 031301, [arXiv:1707.09841].
  * [62] P. Adshead, J. T. Giblin, and Z. J. Weiner, Gravitational waves from gauge preheating, Phys. Rev. D 98 (2018), no. 4 043525, [arXiv:1805.04550].
  * [63] M. A. Amin, J. Braden, E. J. Copeland, J. T. Giblin, C. Solorio, Z. J. Weiner, and S.-Y. Zhou, Gravitational waves from asymmetric oscillon dynamics?, Phys. Rev. D 98 (2018) 024040, [arXiv:1803.08047].
  * [64] P. Adshead, J. T. Giblin, M. Pieroni, and Z. J. Weiner, Constraining axion inflation with gravitational waves from preheating, Phys. Rev. D 101 (2020), no. 8 083534, [arXiv:1909.12842].
  * [65] T. Hiramatsu, E. I. Sfakianakis, and M. Yamaguchi, Gravitational wave spectra from oscillon formation after inflation, JHEP 03 (2021) 021, [arXiv:2011.12201].
  * [66] Y. Cui and E. I. Sfakianakis, Detectable Gravitational Wave Signals from Inflationary Preheating, arXiv:2112.00762.
  * [67] LIGO Scientific Collaboration, J. Aasi et al., Advanced LIGO, Class. Quant. Grav. 32 (2015) 074001, [arXiv:1411.4547].
  * [68] VIRGO Collaboration, F. Acernese et al., Advanced Virgo: a second-generation interferometric gravitational wave detector, Class. Quant. Grav. 32 (2015), no. 2 024001, [arXiv:1408.3978].
  * [69] KAGRA Collaboration, K. Somiya, Detector configuration of KAGRA: The Japanese cryogenic gravitational-wave detector, Class. Quant. Grav. 29 (2012) 124007, [arXiv:1111.7185].
  * [70] M. Punturo et al., The Einstein Telescope: A third-generation gravitational wave observatory, Class. Quant. Grav. 27 (2010) 194002\.
  * [71] LISA Collaboration, P. Amaro-Seoane et al., Laser Interferometer Space Antenna, arXiv:1702.00786.
  * [72] J. Soda and Y. Urakawa, Cosmological imprints of string axions in plateau, Eur. Phys. J. C 78 (2018), no. 9 779, [arXiv:1710.00305].