1 Introduction

A successful scenario for solving shortcomings of the hot Big Bang model and providing a mechanism for creating the primordial density perturbation is based on the accelerated expansion of the early universe driven by a scalar field known as inflaton. This scenario is the cosmic inflation [1,2,3,4]. Notice that the pioneer work [1], where the first full and internally self-consistent inflationary model with the graceful exit from inflation to the final radiation dominated state was developed in the \(R+R^2\) modified gravity, remains viable by now. In principle, the cosmic inflation can be achieved by introducing an additional degree of freedom, e.g., inflaton, in the universe. Hence, the cosmic inflation can also be realized in modified theories of gravity such as the scalar-tensor theories of gravity or f(R) gravity [1] due to the existence of a scalar degree of freedom for gravity in these theories, (see also review articles, e.g., [5, 6]).

Usually, the modification of the theories of gravity can be done by adding extra degrees of freedom to gravitational interaction in the General relativity (GR) [7]. However, there are many attempts to construct a class of minimally modified theories of gravity in which the degrees of freedom for gravity are two tensorial modes similar to GR [8,9,10,11,12]. This class of theories could be constructed by breaking the temporal diffeomorphism invariance of GR. Cosmology in these theories has been investigated in [13,14,15]. A possible form of minimally modified theories of gravity is the cuscuton theory [16, 17] which is constructed by adding a minimally coupled scalar field to GR. To keep the total degrees of freedom unchanged, the scalar field namely the cuscuton field has to be a non-dynamical field in the unitary gauge. However, the existence of the cuscuton field alters the constraint equation of GR and the dispersion relation of the field coupled to the cuscuton.

Influences of cuscuton on Cosmic Microwave Background (CMB) have been studied in [17]. Bouncing cosmology has been investigated in [18,19,20]. In [21], it has been shown that power-law solutions of inflation become viable in cuscuton gravity. A solution describing an accelerating universe with a stable extra dimension in cuscuton gravity has been studied in [22]. Observable predictions in Generalized cuscuton models of inflation have been analyzed in [23]. In this work, we revisit the model proposed in [23] and estimate observational predictions of the cuscuton inflation with two particular forms of the cuscuton potential. We also study preheating effects in the cuscuton inflation.

Regarding existing literature, the preheating mechanism was examined in [24,25,26]. The properties of resonance of models with non-minimally coupled scalar field \({\chi }\) were carried out in [27]. The inclusion of a non-minimal coupling \(\xi R\chi ^{2}\) term with a sizeable range of parameter \(\xi \) produced a sufficient resonance. Higher-curvature inflation models with \((R+\alpha ^{n}R^{n})\) allowing to study a parametric preheating of a scalar field coupled non-minimally to a spacetime curvature were also investigated [28]. In [29], the authors studied preheating effects in the extended Starobinsky model [1] with an additional scalar field which interacts directly with the inflaton field via a four-leg interaction term. Interestingly, regarding multi-field inflation, a preheating mechanism in asymmetric \(\alpha \)-attractors has been discussed in [30], see also [31] for preheating after multifield inflation; while Palatini formalism of gravitational dark matter production during preheating stage was carried out in [32]. The study of preheating in the Palatini formalism with a quadratic inflaton potential with a \(\alpha \,R^2\) term has been discussed in [33].

In Sect. 2, we present our formalism of cuscuton inflation by considering two types of the cuscuton potentials. In Sect. 3, we constrain the model parameters using observational data. We also investigate the preheating effect in the cuscuton inflation in Sect. 4. Here we demonstrate if parametric resonances in our models are sufficiently broad possible for the exponential growth of the number of particles and discuss the case in which the expansion of space is considered. We conclude our findings in Sect. 5.

2 Cuscuton inflation

We study the cuscuton inflation described by the action

$$\begin{aligned} S{} & {} = \int d^4x\,\sqrt{-g} \nonumber \\{} & {} \quad \times \Bigg ( \frac{m_\textrm{pl}^2}{2} R + \mu ^2 \sqrt{- \partial _\mu \varphi \partial ^\mu \varphi } - V(\varphi ) - \frac{1}{2} \partial _\mu \phi \partial ^\mu \phi - U(\phi ) \Bigg ),\nonumber \\ \end{aligned}$$
(1)

where \(m_\textrm{pl}= 1/ \sqrt{8\pi G}\) is the reduced Planck mass, \(\mu \) is a constant with dimension of mass, V and U are the potentials of the cuscuton field \(\varphi \) and inflaton field \(\phi \).

The cuscuton can be a non-dynamical field and consequently a theory of gravity has two degrees of freedom at non-linear level if the cuscuton is homogenous [34]. Nevertheless, if the cuscuton can have inhomogeneous part, there will be instantaneous modes which enhance the degrees of freedom of gravity [23]. Hence, in the following consideration, the cuscuton is supposed to be a homogeneous field.

The spacetime of the background universe is described by the spatially flat Friedmann-Lemaître-Robertson-Walker (FLRW) metric,

$$\begin{aligned} ds^2 = - dt^2 + a^2(t) \delta _{ij} dx^i dx^j\,, \end{aligned}$$
(2)

where a(t) is the cosmic scale factor. The homogeneity and isotropy of the background universe require that the field \(\phi \) depend only on time, so that the equations of motion are [21]

$$\begin{aligned}&3 m_\textrm{pl}^2 H^2 = V + U + \frac{1}{2} \dot{\phi }^2\,, \end{aligned}$$
(3)
$$\begin{aligned}&\textrm{sign}(\dot{\varphi })\, 3 \mu ^2 H + V_\varphi = 0\,, \end{aligned}$$
(4)
$$\begin{aligned}&\ddot{\phi }+ 3 H \dot{\phi }+ U_\phi =0\,, \end{aligned}$$
(5)

where \(H=\dot{a}/a\) is the Hubble parameter. During inflation, the inflaton slowly evolves so that Eqs. (5) and (3) become

$$\begin{aligned} \frac{d \phi }{d N} \simeq \frac{U_\phi }{3H^2}\,, \quad \text{ and }\quad 3 m_\textrm{pl}^2 H^2 \simeq V + U\,, \end{aligned}$$
(6)

where \(N \equiv \ln (a_e/a)\) is the number of e-folding measured from the end of inflation at scale factor \(a = a_e\). The first equation in Eq. (6) can be integrated as

$$\begin{aligned} N_* = \int _{\phi _e}^{\phi _*} \frac{3 H^2(\phi )}{U_{\phi }} d \phi \,, \end{aligned}$$
(7)

where subscripts \({}_e\) and \({}_*\) denote evaluation at the end of inflation and at the time when the perturbation modes relevant to CMB anisotropy cross outside the Hubble horizon. The Hubble parameter in Eq. (7) can be expressed in terms of the inflaton field if the form of the cuscuton potential V is specified. According to the analysis in [16], the cuscuton field affects dynamics of the background universe by modifying the Friedmann equation. It follows from Eq. (3) that this modification depends on the form of the potential V of the cuscuton. We are interested in the quadratic and exponential forms of the cuscuton potential, because the cuscuton has a scaling solution when its potential takes the quadratic form, while for the exponential potential, the modified Friedmann equation takes the form as in the Dvali-Gabadadze-Poratti (DGP) brane-world model [17]. In the following consideration, we choose the inflaton potential in the form

$$\begin{aligned} U(\phi ) = \frac{1}{2} m_\phi ^2 \phi ^2\,, \end{aligned}$$
(8)

where the inflaton mass \(m_\phi \) is constant.

The scalar spectral index \(n_s\) and the tensor-to-scalar ratio r in the cuscuton inflation are [23]

$$\begin{aligned} n_s -1 = -2 \epsilon - \beta \,, \quad \text{ and }\quad r = 16 \alpha \,, \end{aligned}$$
(9)

where

$$\begin{aligned} \epsilon = - \frac{\dot{H}}{H^2}\,, \quad \alpha = \frac{\dot{\phi }^2}{2 m_\textrm{pl}^2 H^2}\,, \quad \beta = \frac{\dot{\alpha }}{\alpha H}\,. \end{aligned}$$
(10)

Note that the prediction from standard inflationary model is \(r = 16\epsilon \) rather than that in Eq. (9).

3 Observational constraints

3.1 Quadratic potential

We first choose the quadratic potential for the cuscuton in the form

$$\begin{aligned} V(\varphi ) = \frac{1}{2}m^{2}\varphi ^2(t)\,, \end{aligned}$$
(11)

where the mass m is constant. Substituting Eq. (11) into Eq. (4), we get

$$\begin{aligned} 9 \mu ^4 H^2 = m^{4}\varphi ^2\,. \end{aligned}$$
(12)

Inserting Eq. (12) into Eq. (3), the friedmann equation becomes

$$\begin{aligned} \left( 3 m_\textrm{pl}^2 - \frac{9}{2} \Lambda \right) H^2 = U + \frac{1}{2} \dot{\phi }^2\,, \end{aligned}$$
(13)

where \(\Lambda \equiv \mu ^4 / m^2\). To estimate observable quantities, we compute a slow-roll parameter \(\epsilon \) by differentiating Eq. (13) with respect to time and applying the slow-roll approximation as

$$\begin{aligned} \epsilon \simeq \left( m_\textrm{pl}^2 - \frac{3}{2} \Lambda \right) \frac{U_\phi ^2}{2 U^2}\,, \end{aligned}$$
(14)

where the slow-roll approximation \(U \gg \dot{\phi }^2\) is used. The number of e-folding can be computed from Eq. (6) as

$$\begin{aligned} N_*= & {} \int _{\phi _e}^{\phi _*} \frac{1}{m_\textrm{pl}^2 - 3 \Lambda / 2} \frac{U}{U_\phi } d \phi \nonumber \\= & {} \frac{1}{2\left( m_\textrm{pl}^2 - 3 \Lambda / 2\right) } \left[ \left. \frac{U^2}{U_\phi ^2}\right| _{\phi _e}^{\phi _*} + \int _{\phi _e}^{\phi _*} \frac{U^3}{U_\phi ^3}\frac{U_{\phi \phi }}{U} d\phi \right] \nonumber \\\sim & {} \frac{1}{2\left( m_\textrm{pl}^2 - 3 \Lambda / 2\right) } \frac{U^2}{U_\phi ^2} - \tilde{N}, \end{aligned}$$
(15)

where \(\tilde{N}\) is a constant and only the dominant term is keeped in the last line. Here, we use the approximation \(U/ (m_\textrm{pl}U_{\phi }) \sim \mathcal{O}(1/ \sqrt{\epsilon })\) and \(m_\textrm{pl}^2 U_{\phi \phi } / U \sim \mathcal{O}(\epsilon )\) as well as \(\epsilon \ll 1\). Inserting Eq. (15) into Eq. (14), we get

$$\begin{aligned} \epsilon _* \sim \frac{1}{4 (N_* + \tilde{N})}\,. \end{aligned}$$
(16)

We see that for slow-roll inflation the slow-roll parameter \(\epsilon \) is not influenced by cuscuton at the leading order if the cuscuton potential takes the quadratic form. From Eq. (10), the parameter \(\alpha \) can be estimated as

$$\begin{aligned} \alpha _*{} & {} = \left. \frac{\dot{\phi }^2}{2 m_\textrm{pl}^2 H^2}\right| _* = \frac{1}{2} \left( \frac{d \phi _*}{d N_*}\right) ^2 \simeq \frac{1}{2}\nonumber \\{} & {} \quad \times \left( m_\textrm{pl}^2 - \frac{3}{2} \Lambda \right) ^2 \frac{U_\phi ^2}{U^2}\, \sim \left( m_\textrm{pl}^2 - \frac{3}{2} \Lambda \right) \frac{1}{4 \left( N_* + \tilde{N}\right) }, \nonumber \\ \end{aligned}$$
(17)

and hence the parameter \(\beta \) is

$$\begin{aligned} \beta _* = - \frac{1}{\alpha _*}\frac{d \alpha _*}{d N_*} \sim \frac{1}{N_* + \tilde{N}}\,. \end{aligned}$$
(18)

It follows from the above estimations and Eq. (9) that \(n_s\) is not affected by the cuscuton at the leading order, while r can be suppressed.

We now consider the specific form of the inflaton potential. Inserting Eqs. (13) and (8) into Eq. (7), the number of e-folding can be computed as

$$\begin{aligned} N_* = \int _{\phi _e}^{\phi _*} \frac{1}{2\left( m_\textrm{pl}^2 - 3 \Lambda / 2\right) } \phi d \phi \,, \end{aligned}$$
(19)

Applying the slow-roll approximation to Eq. (13) and differentiating the result with respect to time, we get

$$\begin{aligned} \epsilon = 2 \left( m_\textrm{pl}^2 - \frac{3}{2} \Lambda \right) \frac{1}{\phi ^2}\,. \end{aligned}$$
(20)

At the end of inflation, \(\epsilon = 1\) so that the above equation yields

$$\begin{aligned} \phi _e^2 = 2 \left( m_\textrm{pl}^2 - \frac{3}{2} \Lambda \right) \,. \end{aligned}$$
(21)

This equation shows that \(\phi _e\) is suppressed in cuscuton inflation compared with the standard one. Substituting this result into Eq. (19), Eq. (19) gives relation between the inflaton field and the number of e-folding as

$$\begin{aligned} \phi _*^2 = 4\left( m_\textrm{pl}^2 - \frac{3}{2} \Lambda \right) \left( N_* +\frac{1}{2}\right) \,. \end{aligned}$$
(22)

Hence, Eq. (20) gives

$$\begin{aligned} \epsilon _* = \frac{1}{2 N_* + 1}\,. \end{aligned}$$
(23)

The parameter \(\alpha \) is computed by diffferentiating Eq. (22) with respect to \(N_*\) as

$$\begin{aligned} \alpha _* = \frac{1}{2 m_\textrm{pl}^2}\left( \frac{d\phi _* }{d N_*}\right) ^2 = \frac{1 - 3 \Lambda / (2 m_\textrm{pl}^2)}{2 N_* + 1}\,. \end{aligned}$$
(24)

From the above expression for \(\alpha _*\), we obtain

$$\begin{aligned} \beta _* = \frac{2}{2 N_* + 1}\,. \end{aligned}$$
(25)

Hence, the scalar spectral index and the tensor-to-scalar ratio at the horizon crossing are

$$\begin{aligned} n_s -1 = - \frac{4}{2 N_* + 1}\,, \quad \text{ and }\quad r = \frac{16 - 24 \Lambda / m_\textrm{pl}^2}{2 N_* + 1}\,, \end{aligned}$$
(26)

The observational bound on these parameters from Planck [35] are

$$\begin{aligned} n_s \sim 0.96\,, \qquad r < 0.1\,. \end{aligned}$$
(27)

Setting \(n_s\) in Eq. (26) according to the above observational value, and plugging the obtained \(N_*\) in the expression of r, we get \(r < 0.1\) if \(\Lambda / m_\textrm{pl}^2 > 0.25\). To avoid negative r, we require \(\Lambda / m_\textrm{pl}^2 < 2 / 3\).

After the end of inflation, the inflaton oscillates around the minimum of its potential. The oscillating solution can be estimated by using variable \(\phi =a^{-3/2}Y\), so that Eq. (5) becomes

$$\begin{aligned} \ddot{Y} +\left( m_\phi ^{2}-\frac{9}{4}H^{2}-\frac{3}{2}\dot{H}\right) Y = 0\,. \end{aligned}$$
(28)

According to the approximations \(m_\phi ^{2}\gg H^{2}\) and \(m_\phi ^{2}\gg |\dot{H}|\) after inflation, Eq. (28) is satisfied by the solution in the form

$$\begin{aligned} Y = A \sin (m_\phi t)\,, \end{aligned}$$
(29)

where A is a constant. Hence, the solution of \(\phi \) is

$$\begin{aligned} \phi = \frac{A}{a^{3/2}} \sin (m_\phi t)\,. \end{aligned}$$
(30)

It is worth mentioning here that a massive minimally coupled scalar field behaves like dust-like matter after inflation (\(mt\gg 1\)), so that \(a(t)\propto t^{2/3}\) after averaging over time interval much larger than \(m^{-1}\), was derived analytically in [36]. Averaging over several oscillations of the inflaton, see also Ref. [37], the energy and pressure densities associated to the solution (30) are given by

$$\begin{aligned} P_{\phi } \approx \left<\frac{1}{2}{\dot{\phi }}^{2} - \frac{1}{2}m^{2}_{\phi }\phi ^{2}\right> \approx 0, \end{aligned}$$
(31)

where we suppose \(m_\phi > H\) after inflation. This implies that the averaged energy density of the inflaton behaves like matter, and therefore Eq. (13) gives \(a \propto t^{2/3}\). As a result, Eq. (30) becomes

$$\begin{aligned} \phi (t) \approx \frac{\phi _{e}}{m_{\phi }t}\,\sin (m_{\phi }t)\equiv \Phi (t) \sin (m_{\phi }t), \end{aligned}$$
(32)

where the amplitude A of the oscillation is computed by matching \(\phi (t \rightarrow 0) = \phi _e\) and \(\phi _e\) in this case is given by Eq. (21). It follows from Eqs. (21) and (26) that the amplitude of inflaton oscillation decreases if r decreases.

3.2 The exponential potential

We now consider the exponential potential in the form

$$\begin{aligned} V(\varphi ) = V_0 e^{-\kappa \varphi }\,, \end{aligned}$$
(33)

where \(V_0\) and \(\kappa \) are constant with dimension of mass\({}^4\) and mas\({}^{-1}\). For the exponential potential, Eq. (4) yields

$$\begin{aligned} \textrm{sign}(\dot{\varphi }) 3 \mu ^2 H = \kappa V\,. \end{aligned}$$
(34)

According to Eqs. (33) and (34), if we deman that the Hubble parameter is positive, the field \(\varphi \) has to increase with time for a positive \(\kappa \) and has to decrease with time for a negative \(\kappa \). Hence, we can drop \( \textrm{sign}(\dot{\varphi })\) and suppose that \(\kappa > 0\) in the following calculations. The relation in Eq. (34) can be used to eliminate V from Eq. (3), so that the friedmann equation becomes

$$\begin{aligned} 3 \left( m_\textrm{pl}^2 H^2 - \frac{\mu ^2}{\kappa } H\right) = U + \frac{1}{2} \dot{\phi }^2\,. \end{aligned}$$
(35)

Hence, the Hubble parameter is given by

$$\begin{aligned} H = \Xi \pm \sqrt{\Xi ^2 + \frac{1}{3m_\textrm{pl}^2} \left( U + \frac{1}{2} \dot{\phi }^2\right) }\,, \end{aligned}$$
(36)

where \(\Xi \equiv \mu ^2 / (2 \kappa m_\textrm{pl}^2)\). In the following calculations we choose the solution with plus sign to avoid the negative H. From Eq. (36), we have

$$\begin{aligned} 3H^2 = 6\Xi ^2 + \frac{\rho _\phi }{m_\textrm{pl}^2} + 2 \Xi \sqrt{9 \Xi ^2 + \frac{3}{m_\textrm{pl}^2} \rho _\phi }\,, \end{aligned}$$
(37)

where \(\rho _\phi \equiv \dot{\phi }^2 / 2 + U\) is the energy density of the inflaton. To investigate situation in which the inflationary universe can be realized, we compute the slow-roll parameter \(\epsilon \) by differentiating Eq. (36) with respect to time:

$$\begin{aligned} \epsilon = \frac{\sqrt{3}}{2 m_\textrm{pl}H} \frac{\dot{\phi }^2}{\sqrt{3 m_\textrm{pl}^2 \Xi ^2 + \rho _\phi }}\,. \end{aligned}$$
(38)

The slow-roll parameter \(\epsilon \) can be less than unity corresponding to inflationary phase and can be larger than one corresponding to epoch after inflation, i.e., the model has graceful exit, if the energy density of inflaton \(\rho _\phi \) is larger than \(m_\textrm{pl}^2 \Xi ^2\). In this situation, the inflationary phase can be realized if the inflaton field slowly evolves. Supposing that the energy density of inflaton is significantly larger than \(m_\textrm{pl}^2 \Xi ^2\), Eq. (37) can be approximated as

$$\begin{aligned} 3H^2 \simeq \frac{2\sqrt{3}}{m_\textrm{pl}} \sqrt{\rho _\phi } \Xi + \frac{\rho _\phi }{m_\textrm{pl}^2}\,. \end{aligned}$$
(39)

Applying the slow-roll approximation to Eq. (39), we can write

$$\begin{aligned} 3H^2 \simeq \sqrt{6} \bar{m}_\phi \phi \Xi + \frac{\bar{m}_\phi ^2}{2} \phi ^2\,, \end{aligned}$$
(40)

where \(\bar{m}_\phi = m_\phi / m_\textrm{pl}\) and we have used the inflaton potential from Eq. (8). Differentiating Eq. (40) with respect to time, the slow-roll parameter is obtained under our assumption as

$$\begin{aligned} \epsilon = \frac{2 m_\textrm{pl}^2}{\phi } \frac{\sqrt{6} \bar{m}_\phi \Xi + \bar{m}_\phi ^2 \phi }{\left( 2 \sqrt{6} \Xi + \bar{m}_\phi \phi \right) ^2}\,. \end{aligned}$$
(41)

Setting \(\epsilon = 1\), the inflaton field at the end of inflation can be computed from the relation

$$\begin{aligned} \phi _e \left( 2 \sqrt{6} \Xi + \bar{m}_\phi \phi _e\right) ^2 = 2 m_\textrm{pl}^2 \left( \sqrt{6} \bar{m}_\phi \Xi + \bar{m}_\phi ^2 \phi _e\right) \,. \end{aligned}$$
(42)

If we suppose that \(\phi \sim \mathcal{O}(m_\textrm{pl})\), the approximation \(\rho _\phi \gg m_\textrm{pl}^2 \Xi ^2\) is equivalent to \(m_\phi \gg \Xi \). Further assuming that \(\bar{m}_\phi < 1\), the approximation becomes \(m_\textrm{pl}\gg \Xi \). The solutions for Eq. (42) under this approximation are

$$\begin{aligned}{} & {} \phi _e = \sqrt{2} m_\textrm{pl}- \frac{3 \sqrt{3/2}}{\bar{m}_\phi } \Xi , \quad \phi _e = - \sqrt{2} m_\textrm{pl}- \frac{3 \sqrt{3/2}}{\bar{m}_\phi } \Xi , \nonumber \\{} & {} \quad \phi _e = - \frac{\sqrt{6}}{\bar{m}_\phi }\Xi . \end{aligned}$$
(43)

We are interested in the first solution. The above approximated \(\phi _e\) differs by less than a few percent from that numerically computedfrom Eq. (38) under the slow-roll approximation if \(\Xi / m_\phi \le 0.05\). Supstituting Eq. (40) into Eq. (7), the number of e-folding for this case can be computed as

$$\begin{aligned} N_* = \int _{\phi _e}^{\phi _*} \frac{1}{m_\phi m_\textrm{pl}} \left( \sqrt{6} \Xi + \frac{\bar{m}_\phi }{2} \phi \right) d \phi \,. \end{aligned}$$
(44)

The above equation yields

$$\begin{aligned} \phi _* = \frac{2}{\bar{m}_\phi }\left( - \sqrt{6} \Xi \pm \sqrt{6 \Xi ^2 + m_\phi ^2 \left( N_* + \tilde{N}\right) }\right) \,, \end{aligned}$$
(45)

where

$$\begin{aligned} \tilde{N}\equiv \frac{4 \bar{m}_\phi ^2 m_\textrm{pl}^2 + 4 \sqrt{3} \bar{m}_\phi m_\textrm{pl}\Xi - 45 \Xi ^2}{8 m_\phi ^2}\,. \end{aligned}$$
(46)

We will consider the solution with the plus sign in the following analysis. Substituting Eq. (45) into Eq. (41), we get

$$\begin{aligned}{} & {} \epsilon _* \nonumber \\{} & {} {=}\frac{m_\phi ^2}{4\left( {-} \sqrt{6} \Xi {+} \sqrt{6 \Xi ^2 {+} m_\phi ^2 \left( N_* {+} \tilde{N}\right) }\right) \sqrt{6 \Xi ^2 {+} m_\phi ^2 \left( N_* {+} \tilde{N}\right) }} \nonumber \\{} & {} \quad + \frac{ m_\phi ^2 }{4 \left( 6 \Xi ^2 + m_\phi ^2 \left( N_* + \tilde{N}\right) \right) }. \end{aligned}$$
(47)

The parameter \(\alpha \) can be computed by differentiating Eq. (45) with respect to \(N_*\) as

$$\begin{aligned} \alpha _* = \frac{m_\phi ^2}{12 \Xi ^2 + 2 m_\phi ^2 \left( N_* + \tilde{N}\right) }\,. \end{aligned}$$
(48)

This expression of \(\alpha _*\) yields

$$\begin{aligned} \beta _* = \frac{m_\phi ^2}{6 \Xi ^2 + m_\phi ^2 \left( N_* + \tilde{N}\right) }\,. \end{aligned}$$
(49)

Hence, the predictions for this model are

$$\begin{aligned} n_s -1{} & {} = - \frac{m_\phi ^2}{2\left( - \sqrt{6} \Xi + \sqrt{6 \Xi ^2 + m_\phi ^2 \left( N_* + \tilde{N}\right) }\right) \sqrt{6 \Xi ^2 + m_\phi ^2 \left( N_* + \tilde{N}\right) }}\nonumber \\{} & {} \quad - \frac{ 3 m_\phi ^2 }{2 \left( 6 \Xi ^2 + m_\phi ^2 \left( N_* + \tilde{N}\right) \right) } \simeq - \frac{4}{2 N_*+ 1} + 2 \sqrt{3}\frac{\Xi }{m_\phi } \frac{2 - \sqrt{2 N_* + 1}}{\left( 2 N_* + 1\right) ^2}, \end{aligned}$$
(50)
$$\begin{aligned} r{} & {} = 16 \frac{m_\phi ^2}{12 \Xi ^2 + 2 m_\phi ^2 \left( N_* + \tilde{N}\right) } \simeq \frac{16}{2 N_* + 1} - \frac{16 \sqrt{3} \Xi }{m_\phi \left( 2 N_* + 1\right) ^2}. \end{aligned}$$
(51)

We estimate how much the observational quantities in this model of cuscuton inflation differ from those in standard inflation by computing the factions:

$$\begin{aligned} \frac{n_s - n_s^{(s)}}{n_s^{(s)} - 1}\sim & {} -\frac{\sqrt{3}}{2\sqrt{2}}\frac{\Xi }{m_\phi \sqrt{N_*}} < - 0.01, \end{aligned}$$
(52)
$$\begin{aligned} \frac{r - r^{(s)}}{r^{(s)}}\sim & {} - \frac{\sqrt{3} \Xi }{2 m_\phi N_*} < - 0.01, \end{aligned}$$
(53)

where superscript \({}^{(s)}\) denotes quantities from standard inflation, i.e., \(\Xi = 0\). The upper bound \(\Xi / m_\phi < 0.05\) is used in the estimation of Eqs. (52) and (53) to ensure that the above approximations have error up to a few percent.

After inflation, the inflaton field starts to oscillate around local minimum of its potential. To estimate the evolution of the scale factor, we have to average Eq. (39) over several oscillations of the inflaton

$$\begin{aligned} 3 \langle H^2 \rangle \simeq \frac{2\sqrt{3}}{m_\textrm{pl}} \langle \sqrt{\rho _\phi }\rangle \Xi + \frac{\langle \rho _\phi \rangle }{m_\textrm{pl}^2}\,. \end{aligned}$$
(54)

According to Eq. (5), the oscillation of inflaton depends on the form of its potential while the Hubble parameter plays a role of the dampping term. Hence, the averaged evolution of inflaton for this case is the same as that in the previous section such that averaged energy density of inflaton behaves like matter. According to Eq. (30), we suppose that the damped oscillation of the inflaton can be described by

$$\begin{aligned} \phi (t) = \phi _0(t) \sin (m_{\phi }t)\,, \end{aligned}$$
(55)

so that the energy density of inflaton is

$$\begin{aligned} \rho _{\phi } {=} \frac{1}{2} \dot{\phi }_{0}^{2} \sin ^2(m_{\phi }t) {+} \dot{\phi }_{0} m_{\phi } \phi _{0} \sin (m_{\phi }t)\cos (m_{\phi }t) {+} \frac{1}{2} m_{\phi }^{2} \phi _{0}^{2}. \end{aligned}$$
(56)

Since \(\phi _{0}(t)\) is a function of a scale factor, one can suppose that \(\dot{\phi }_{0}/ \phi _{0} \sim \mathcal{O}(H)\). Using \(m_{\phi } > H\) during preheating, we can approximate

$$\begin{aligned} \sqrt{\rho }_{\phi } \approx \sqrt{\frac{1}{2}} \dot{\phi }(t)_{0}\sin (m_{\phi }t)\cos (m_{\phi }t) + \sqrt{\frac{1}{2}} m_{\phi } \phi _{0}(t). \end{aligned}$$
(57)

Averaging Eq. (57) over several oscillations of the inflaton, we get \(\langle \sqrt{\rho }_\phi \rangle \approx \sqrt{1/2} m_\phi \phi _0(t)\) and consequently Eq. (54) yields

$$\begin{aligned} 3 H^2 \approx \sqrt{6} \bar{m}_\phi \Xi \phi _0(t) + \frac{1}{2} \bar{m}_\phi ^2 \phi _0^2(t)\,. \end{aligned}$$
(58)

Setting \(\phi _0(t) = A / a^{(3/2)}\), the above equation can be integrated as

$$\begin{aligned} 4 \frac{\sqrt{\sqrt{6} \bar{m}_\phi \Xi A a^{3/2} + \bar{m}_\phi ^2 A^2 / 2}}{3\sqrt{2} \bar{m}_\phi \Xi A} \approx t + C\,. \end{aligned}$$
(59)

The integration constant C is chosen such that the above equation reduces to \(a \propto t^{2/3}\) when \(\Xi \rightarrow 0\), so that we get \(C = 2 / (3 \Xi )\) and therefore

$$\begin{aligned} a^{3/2}\approx & {} \frac{3 \sqrt{6} \bar{m}_\phi \Xi A}{16} t^2 + \frac{\sqrt{3} \bar{m}_\phi A}{2 \sqrt{2}} t \nonumber \\\approx & {} \frac{\sqrt{3} \bar{m}_\phi A}{2 \sqrt{2}} t\left( 1 + \frac{3 \Xi }{4} t\right) . \end{aligned}$$
(60)

Inserting this result into Eq. (30), we get

$$\begin{aligned} \phi \approx B \left( 1 + \frac{3 \Xi }{4} t\right) ^{-1} \frac{2 \sqrt{2}}{\sqrt{3} \bar{m}_\phi t}\sin (m_\phi t)\,, \end{aligned}$$
(61)

where B is a constant. Matching \(\phi (m_\phi t \rightarrow 0) \rightarrow \phi _e\), we get

$$\begin{aligned} \phi \approx \phi _e \left( 1 + \frac{3 \Xi }{4} t\right) ^{-1} \frac{\sin (m_\phi t)}{m_\phi t} \equiv \bar{\Phi }(t) \sin (m_\phi t)\,, \end{aligned}$$
(62)

where \(\phi _e\) is given by Eq. (43). After the end of inflation, the averaged energy density of the inflaton evolves as \(\rho _\phi \propto a^{-3}\). Hence, the energy density of the inflaton can drop below \(m_\textrm{pl}^2 \Xi ^2\). Consequently, the universe could start an inflationary phase again as followed from Eq. (38). This implies that \(\Xi \) has to be suppressed as the universe evolves which could be achieved if \(\mu \) or \(\kappa \) depend on time, otherwise \(\Xi \) has to be small such that it has no effect during inflationary phase.

4 Preheating

In this section, we study parametric resonances of models when an inflaton field \(\phi \) coupled to another scalar field \(\chi \) with the interaction term \(g^{2}\phi ^{2}\chi ^{2}\), so that the action for our model is

$$\begin{aligned} S{} & {} = \int d^4x\,\sqrt{-g} \Bigg ( \frac{m_\textrm{pl}^2}{2} R + \mu ^2 \sqrt{- \partial _\mu \varphi \partial ^\mu \varphi } - V(\varphi ) \nonumber \\{} & {} \quad - \frac{1}{2} \partial _\mu \phi \partial ^\mu \phi - U(\phi ) \nonumber \\{} & {} \quad -\frac{1}{2} \partial _\mu \chi \partial ^\mu \chi - W(\chi ) -\frac{1}{2} g^{2}\phi ^{2}\chi ^{2} \Bigg ), \end{aligned}$$
(63)

where \(W(\chi )\) is the potential of the field \(\chi \) and g is the coupling constant. In this case scenario, we will choose the potential \(W(\chi )\) of the form:

$$\begin{aligned} W(\chi ) = \frac{1}{2}m^{2}_{\chi }\chi ^{2} , \end{aligned}$$
(64)

where the mass \(m_{\chi }\) is constant. Models with different form of \(W(\chi )\) may also receive physical interest, e.g., Refs. [38, 39] for a non-minimally coupled scenario of the form \(\sim \xi R\chi ^{2}\). The time evolution of the quantum fluctuations of the field \(\chi \) is ruled by the classical equation of motion, a.k.a., the Klein-Gordan equation in an expanding flat FLRW universe. From Eq. (63), it is rather straightforward to derive the equation of motion for the field \(\chi \) to obtain

$$\begin{aligned} \ddot{\chi }+3H\dot{\chi } - \frac{1}{a^{2}}\nabla ^{2}\chi +\Big [m^{2}_{\chi } + g^{2}\phi ^{2}\Big ] \chi = 0, \end{aligned}$$
(65)

where an effective \(\chi \) mass \(m^{2}_\mathrm{eff.}=m^{2}_{\chi }+g^{2}\phi ^{2}\). We then expand the scalar fields \(\chi \) in terms of the Heisenberg representation to yield

$$\begin{aligned} \chi (t,\textbf{x}) \sim \int \left( a_{k}\chi _{k}(t)e^{-i\textbf{k}\cdot \textbf{x}} + a^{\dagger }_{k}\chi ^{*}_{k}(t)e^{i\textbf{k}\cdot \textbf{x}}\right) d^{3}{} \textbf{k} , \end{aligned}$$
(66)

where \(a_{k}\) and \(a^{\dagger }_{k}\) are annihilation and creation operators. We can show that \(\chi _{k}\) obeys the following equation of motion:

$$\begin{aligned} \ddot{\chi }_{k}+3H\dot{\chi }_{k} +\Big [ \frac{k^{2}}{a^{2}} + m^{2}_{\chi } + g^{2}\phi ^{2}\Big ]\chi _{k} = 0. \end{aligned}$$
(67)

Performing Fourier transformation to this equation and rescaling the field using \(Y_{k} = a^{3/2}\chi _{k}\), we have

$$\begin{aligned} \ddot{Y}_{k} + \omega ^{2}_{k}Y_{k} = 0, \end{aligned}$$
(68)

where a time dependent frequency of \(Y_{k}\) is given by

$$\begin{aligned} \omega ^{2}_{k} = \frac{k^{2}}{a^{2}} + m^{2}_{\chi }+ g^{2}\phi ^{2}. \end{aligned}$$
(69)

As is expected, Eq. (68) describes an oscillator with a periodically changing frequency \(\omega ^{2}_{k} = \frac{k^{2}}{a^{2}} + m^{2}_{\chi }+ g^{2}\phi ^{2}\). The physical momentum \(\textbf{p}\) coincides with \(\textbf{k}\) for Minkowski space such that \(k=\sqrt{\textbf{k}^{2}}\). The periodicity of Eq. (68) may drive the parametric resonance for modes with certain values of k. We suppose that the effective mass of the field \(\chi \) vanishes \(m_{\chi }=0\) and neglect for a moment the expansion of the universe, taking \(a=1\) in Eq. (67).

Let us consider for the first scenario the quadratic cuscuton potential examined in Sect. 3.1. The simplest way to describe this important effect is to make a change of variables \(m_{\phi } t \equiv z\). This simplifies Eq. (68) to the well-known Mathieu equation governed by

$$\begin{aligned} \frac{d^{2}Y_{k}}{dz^{2}} + \left( A_{k} - 2q\cos (2z)\right) Y_{k} = 0. \end{aligned}$$
(70)

The frequency terms are

$$\begin{aligned} A_{k} = \frac{k^{2}}{m_{\phi }^{2}}+2q,\,\, q = \frac{g^{2}\Phi ^{2}(t)}{4m_{\phi }^{2}} . \end{aligned}$$
(71)

An important feature of the solution of Mathieu equation is the existence of an exponential instability \(Y_{k} \propto \exp (\mu _{k}z)\). This instability corresponds to an exponential growth of occupation number of quantum fluctuation \(n_{k}(t) \propto \exp (2\mu _{k}z)\) [38, 39]. The modes \(Y_{k}\) with momentum corresponding to the center of the resonance at \(k\sim m_{\phi }\) grow as \(e^{qz/2}\) which in this work equals \(e^{qz/2}\sim e^{g^{2}\Phi ^{2}t/ 8m_{\phi }}\). Then the number of \(\chi \)-particles grows as \(e^{2\mu _{k}z}\sim e^{qz}\sim e^{g^{2}\Phi ^{2}t/4m_{\phi }}\) [38, 39]. From Eq. (32), we have \(q= g^{2}\phi _e^{2}/4 t^2\,m^{4}_{\phi }\), so that q can be very large if \(m_\phi \ll m_\textrm{pl}\). In this regime, the resonance occurs for a broad range of values of k. Therefore, the parameter \(\mu _{k}\) can be also large in which the resonance occurs for modes with \(k^{2}/m^{2}_{\phi }=A-2q\). Specifically, the broad resonance occurs above the line \(A=2q\) on the stability/instability chart for the Mathieu equation. We recommend the readers to follow the work done by [38, 39]. for detailed analysis on the stability/instability chart for the Mathieu equation [40].

In the second case in which the cuscuton potential takes the exponential form discussed in Sect. 3.2, we can also have an equation which describes an oscillator with a periodically changing frequency \(\omega ^{2}(t) = k^{2} + g^{2}\bar{\Phi }^{2} \sin ^{2} (m_{\phi }t)\). One can have a Mathieu equation with \(A_{k} = \frac{k^{2}}{m_{\phi }^{2}}+2q\), where \(q = \frac{g^{2}\bar{\Phi }^{2}}{4m_{\phi }^{2}}\). Similar to the case of the quadratic potential, a broad parametric resonance for the field \(\chi \) in Minkowski space for \(q \gg 1\) in this scenario can be obtained if \(m_\phi \ll m_\textrm{pl}\). However, q decreases faster than the case of the quadratic potential due to a factor \(f \equiv (1 + 3 \Xi t/ 4)^{-1}\). For a rough estimation, \(H \sim 1/t\) so that \(3 \Xi t/ 4 \sim 3 \Xi / (4\,H) \sim (3/4) (\Xi /m_\phi ) (m_\phi / H)\) which can be in order of unity. This suggests that a factor f can significantly alter the decreasing rate of q. To investigate how does the factor f influence the broad parametric resonance, we compute q for two cases of potentials. For the quadratic potential, we have

$$\begin{aligned} q = \frac{g^{2}m_\textrm{pl}^{2}}{2m_{\phi }^4 t^{2}}\left( 1 - \frac{3}{2} \frac{\Lambda }{m_\textrm{pl}^{2}} \right) ^{2}\,, \end{aligned}$$
(72)

while for the exponential potential, we have

$$\begin{aligned} q = \frac{g^{2}m_\textrm{pl}^{2}}{2 m_{\phi }^{4} t^{2}\Big (1-\frac{3\Xi t}{4}\Big )^{2}} \Bigg (1+\frac{3\sqrt{3}\Xi }{2m_{\phi }}\Bigg )^{2}\,. \end{aligned}$$
(73)

If we suppose that the broad resonance stops at \(q \gtrsim 1\), the broad resonance stops at

$$\begin{aligned} t_{s\,q}{} & {} \sim \frac{g m_\textrm{pl}}{\sqrt{2} m^{2}_{\phi }}, \quad t_{s\,e} \sim \frac{2}{3 \Xi }\left( - 1 + \sqrt{1 + 3 \Xi t_{s\,q}}\right) \nonumber \\{} & {} =\frac{2\,t_{s\,q}}{3 \Xi t_{s\,q}}\left( - 1 + \sqrt{1 + 3 \Xi t_{s\,q}}\right) , \end{aligned}$$
(74)

for the quadratic (\(t_{s\,q}\)) and exponential (\(t_{s\,e}\)) potentials, respectively. For the exponential potential, when \(3\Xi t_{s\,q}< 1\), we can expand \(\sqrt{1 + 3 \Xi t_{s\,q}}\sim 1+3 \Xi t_{s\,q}/2\) to simply find that \(t_{s\,e} \sim t_{s\,q}\). When \(3\Xi t_{s\,q}\sim 1\), we get \(t_{s\,e} \sim 0.8 t_{s\,q}\). However, in case of \(3\Xi t_{s\,q}> 1\), we instead find that \(t_{s\,e}<t_{s\,q}\). The later case implies that the broad resonance in the exponential potential stops earlier than that of the quadratic one.

Nevertheless, a sufficient broad parametric resonance could be achieved if q is initially large. Consider the case where the expansion of space is taken into account. Let’s take only the first model and we see that the equation of motion for the scalar matter fields, Eq. (68), can still be written in the form of a simple harmonic oscillator with a time varying frequency. Using the the same technique of performing Fourier transformation and taking \(Y_{k}(t) = a(t)^{3/2}\chi _{k}(t)\), we have

$$\begin{aligned} \ddot{Y}_{k} + \omega ^{2}(k,t)\,Y_{k} = 0, \end{aligned}$$
(75)

where a time dependent frequency of \(Y_{k}(t)\) in this case is given by

$$\begin{aligned} \omega ^{2}(k,t) = \frac{k^{2}}{a^{2}} + m^{2}_{\chi }+ g^{2}\phi ^{2}-\frac{9}{4}H^{2}-\frac{3}{2}{\dot{H}}. \end{aligned}$$
(76)

It is clear that the above equation is not the Hill’s equation anymore. To quantify the evolution of particular comoving modes, the effects of the space expansion on particle production can be traced by considering the flow lines to the Floquet chart [41]. Note that in the matter-dominated background with \(U(\phi )= m_{\phi }^{2} \phi ^{2}/2\), and \(\Phi \propto a^{-3/2},\,k\propto a^{-1}\) and \(3\,H^{2}=-2{\dot{H}}\), the last two terms cancel. A given comoving mode k flows along \(\Phi \sim k^{3/2}_{phys}\equiv (k/a)^{3/2}\) curve in \(k_{phys}-\Phi \) plane, see for example [41].

The parameter \(q = g^{2}\Phi ^{2}/4m^{2}_{\phi }\) depends on time via \(\Phi \propto a^{-3/2}\). Additionally, we can check in which a resonance band in our process develops. Following [39], the number of the band in the theory of the Mathieu equation is given by \(n=\sqrt{A}\). In our case, reheating occurs for \(A\sim 2q\), i.e. \(n\sim \sqrt{2q}\sim \sqrt{g^{2}\Phi ^{2}/2\,m^{2}_{\phi }}\). Suppose we have an inflationary theory with \(m_{\phi }\sim 10^{-6} m_\textrm{pl}\), and let us take as an example \(g\sim 10^{-2}\). Then after the first oscillation of the field, we have \(\Phi (t)\sim 0.2 m_\textrm{pl}\sqrt{\left( 1 - \frac{3}{2} \frac{\Lambda }{m_\textrm{pl}^2}\right) }\sim 0.2 m_\textrm{pl}\) for \(\Lambda \ll m_\textrm{pl}\), which corresponds to \(q\sim 10^{6}\). This gives the band number \(n \sim 1414\).

5 Conclusions

In this work, we investigate observational predictions and preheating in cuscuton inflation. For the case where the potential of the cuscuton takes the quadratic form, the scalar spectral index is the same as for the standard inflationary model at the leading order if the slow-roll approximation is applied. Under the slow-roll approximation, the tensor-to-scalar ratio can be reduced such that its value agrees with the Planck data for any form of inflaton potential. If the cuscuton potential has exponential form, the energy density of the inflaton is required to be larger than the energy scale \(m_\textrm{pl}^2 \Xi ^2\) to have a graceful exit. According to this condition, the scalar spectral index and tensor-to-scalar ratio from this model deviate by a few percent compared with those from standard inflation. The existence of the constant energy scale \(\Xi \) leads to the problem when the energy density of the inflaton drops below this energy scale after inflation because the universe could accelerate again. This problem could be alleviated if \(\Xi \) can vary with time or its value is small such that it has no effect during inflation but could drive cosmic acceleration at late time.

We have also examined the particle production due to parametric resonances in both models. We demonstrate that in Minkowski space the stage of parametric resonances can be described by the Mathieu equation. For the case of quadratic potential, the amplitude of the driving force in the Mathieu equation has a similar form as that in the standard inflationary framework. Nevertheless, in the case of exponential potential, the amplitude of the driving force decreases faster than in the standard case. We find that the parametric resonances in our models can be sufficiently broad possible for the exponential growth of the number of particles. We also discuss the case in which the expansion of space is considered.