1 Introduction

The Starobinsky model of inflation, proposed by Starobinsky [1], is a geometric model that incorporates linear and quadratic terms of the scalar curvature in its action, distinguishing it from other inflationary models (for reviews on inflation, see e.g., [2,3,4,5]). By expressing the action in the Einstein frame, a scalar field potential emerges, aligning the model with typical inflation models. However, despite being proposed more than 40 years ago, most of the literature has focused on determining ranges for cosmological quantities like the scalar spectral index \(n_s\), tensor-to-scalar ratio r, and the number of e-folds during inflation, denoted by \(N_k\). In this study, we adopt an approach that allows us to derive solutions for the observables, as well as the number of e-folds during inflation, reheating and radiation, with minimal assumptions. We impose reheating conditions on inflation and obtain equations that enable us to solve for the desired quantities (for reviews on reheating, see e.g., [6,7,8]). To understand the impact of reheating on inflation, we constrain the equation of state parameter (EoS) at the end of reheating, denoted as \(\omega _{re}\). By establishing a connection between inflation and reheating, we derive an equation that determines the scalar spectral index, \(n_s\). Using consistency relations within the model, we determine the remaining observables. Furthermore, we calculate the number of e-folds during inflation, reheating, and radiation, denoted as \(N_k\), \(N_{re}\), and \(N_{rd}\), respectively.

The Starobinsky model modifies Einstein’s theory of general relativity by introducing additional terms in the action. The action consists of the Einstein-Hilbert term, which is proportional to the Ricci scalar R, and an extra term proportional to the square of the Ricci scalar, \(R^2\). The original Starobinsky model’s action is given by

$$\begin{aligned} S =\int d^4 x \sqrt{-g}\left[ \frac{M_{Pl}^2}{2}\left( R + \frac{1}{6M^2}R^2\right) + L_m \right] , \end{aligned}$$
(1.1)

where g is the determinant of the metric tensor, and \(M_{Pl}=2.44\times 10^{18} \,\textrm{GeV}\) is the reduced Planck mass. The parameter M is related to the energy scale of inflation, and \(L_m\) represents the Lagrangian density for matter fields. The inclusion of the \(R^2\) term leads to modified equations of motion and a modified theory of gravity. During inflation, the \(R^2\) term dominates over the Einstein-Hilbert term, resulting in exponential expansion. The \(R^2\) term introduces a repulsive interaction that counteracts the attractive gravitational behavior, driving the accelerated expansion of the universe. This unique feature of the Starobinsky model provides a mechanism for inflation based solely on modifications to the gravitational sector, without the need for additional scalar fields. The dominance of the \(R^2\) term during inflation has implications for the dynamics of gravity and the resulting accelerated expansion. It allows for a prolonged period of inflation, resolving significant cosmological puzzles such as the horizon problem and the flatness problem. This extended exponential expansion is responsible for the observed large-scale homogeneity and isotropy of the universe.

The model, originally defined in the Jordan frame, is equivalent to a single-field model with an asymptotically flat potential when transformed through a conformal transformation to the Einstein frame. Additionally, by considering the Standard Model fields as minimally coupled to gravity in the Jordan frame, the transformation to the Einstein frame induces a coupling between these fields and the inflaton, which provides a natural mechanism for graceful exit and reheating [9,10,11].

The model in the Einstein frame is obtained through a conformal transformation applied to the metric. This transformation is given by

$$\begin{aligned} g_{\mu \nu }\rightarrow e^{\sqrt{\frac{2}{3}}\frac{\phi }{M_{Pl}}}g_{\mu \nu }. \end{aligned}$$
(1.2)

Applying this transformation to the Starobinsky model, we obtain the action

$$\begin{aligned} S=\int d^{4}x \sqrt{-g}\left( \frac{M_{Pl}^{2}}{2}R-\frac{1}{2}\partial _{\mu }\phi \partial ^{\mu }\phi -V(\phi )\right) . \end{aligned}$$
(1.3)

Within this formulation, the potential of the scalar field takes the following form

$$\begin{aligned} V(\phi )=V_{0}\left( 1-e^{-\sqrt{\frac{2}{3}}\frac{\phi }{M_{Pl}}}\right) ^{2}, \end{aligned}$$
(1.4)

where \(V_{0}\) is defined as \(\frac{3}{4}M_{Pl}^{2}M^2\). By having the model expressed in the Einstein frame and identifying a scalar potential \(V(\phi )\), we can employ the usual expressions for single-field inflation under the slow-roll approximation.

The Starobinsky model is remarkably consistent with current measurements of the power spectrum of primordial curvature fluctuations and the constraints on primordial gravitational waves. These measurements are derived from collaborations such as Planck and BICEP/Keck [12,13,14]. We consider the bounds provided by the Table 3 of [12] for the cosmological model \(\Lambda \)CDM\(+r+dn_s/d\ln k\) with the data set Planck TT,TE,EE+lowE+lensing +BK15+BAO. These bounds offer constraints on the parameters and observables within the specific cosmological model and data combination

$$\begin{aligned} n_s =0.9658\pm 0.0040 \quad (68\%\,\, C.L.), \end{aligned}$$
(1.5)
$$\begin{aligned} r < 0.068 \quad (95\%\,\, C.L.). \end{aligned}$$
(1.6)

The Starobinsky model predicts specific patterns in the anisotropies of the cosmic microwave background (CMB), including a characteristic damping scale in the power spectrum. Future observations from experiments such as the Simons Observatory [15] and the CMB-S4 [16] collaboration are anticipated to yield more accurate measurements of the CMB. The Simons Observatory comprises four telescopes positioned at an elevation of 5200 ms in the Atacama Desert of Chile, while the CMB-S4 collaboration aims to enhance our understanding of key cosmological parameters by conducting precise measurements of the CMB’s temperature and polarization. Specifically, the improved measurements of the CMB are expected to provide additional constraints on the Starobinsky model and other inflationary models, refining our understanding of the early universe and its evolution. These experiments aim to study the CMB with high sensitivity and precision using advanced detectors, larger arrays, multiple frequencies, improved angular resolution, and careful site selection. These experiments will provide valuable data to refine our understanding of the early universe, test inflationary models, and investigate fundamental physics.

The paper is organized as follows: Sect. 2 introduces the concept of reheating and establishes a general equation that connects inflation with reheating. This equation is then used to solve for \(n_s\) in the context of the Starobinsky model, as presented in Sect. 3, where we also provide definitions for relevant quantities. Additionally, we derive an expression for the number of e-folds during inflation in terms of the spectral index \(n_s\). Furthermore, we employ the consistency relations of the model to determine other observables. The analysis also encompasses the calculation of the number of e-folds during inflation, reheating, and the radiation era. Finally, Sect. 4 provides the conclusion of our paper.

2 Reheating constraints

Models of inflation can be related to cosmological observables, which, to first order in the slow-roll (SR) approximation, are expressed as (see, for example, [3, 17])

$$\begin{aligned} n_{t}= & {} -2\epsilon = -\frac{r}{8}, \end{aligned}$$
(2.1)
$$\begin{aligned} n_{s}= & {} 1+2\eta -6\epsilon , \end{aligned}$$
(2.2)
$$\begin{aligned} n_{sk}\equiv \frac{d n_s}{d \ln k}= & {} 16\epsilon \eta -24\epsilon ^{2}-2\xi _2, \end{aligned}$$
(2.3)
$$\begin{aligned} n_{tk}\equiv \frac{d n_t}{d \ln k}= & {} 4\epsilon \left( \eta -2\epsilon \right) , \end{aligned}$$
(2.4)
$$\begin{aligned} A_s= & {} \frac{1}{24\pi ^{2}} \frac{V}{\epsilon \, M_{Pl}^{4}}. \end{aligned}$$
(2.5)

Here, \(M_{Pl}=2.43568\times 10^{18} \textrm{GeV}\) is the reduced Planck mass, r denotes the tensor-to-scalar ratio, \(n_s\) represents the scalar spectral index, \(n_{sk}\) its running (which is commonly denoted as \(\alpha \)), \(n_t\) represents the tensor spectral index, and \(n_{tk}\) represents its running, in a self-explanatory notation. The amplitude of density perturbations at a particular wave number k is denoted by \(A_s\). All quantities are evaluated at the moment of horizon crossing at wavenumber \(k=0.05\)/Mpc. The SR parameters involved in the above expressions are

$$\begin{aligned} \epsilon \equiv \frac{M_{Pl}^{2}}{2}\left( \frac{V^{\prime }}{V }\right) ^{2},\quad \eta \equiv M_{Pl}^{2}\frac{V^{\prime \prime }}{V}, \quad \xi _2 \equiv M_{Pl}^{4}\frac{V^{\prime }V^{\prime \prime \prime }}{V^{2}}, \end{aligned}$$
(2.6)

where primes on V denote derivatives with respect to the inflaton \(\phi \).

Expanding on earlier work [17,18,19], it is possible to derive an equation for the number of e-folds during reheating [20, 21] by relating the comoving Hubble scale wavenumber k at horizon crossing to the present scale wavenumber \(k_0=a_0 H_0\) as follows (also see [22, 23] for further details)

$$\begin{aligned} N_{re}= \frac{4}{1-3\, \omega _{re}}\left( -N_{k}-\frac{1}{3} \ln \left[ \frac{11 g_{s,re}}{43}\right] -\frac{1}{4} \ln \left[ \frac{30}{\pi ^2 g_{re} } \right] -\ln \left[ \frac{ k}{a_0 T_0} \right] -\frac{1}{4}\ln \left[ \frac{\rho _e }{H_k^4} \right] \right) . \end{aligned}$$
(2.7)

In the above equation, the number of degrees of freedom of species at the end of reheating is denoted by \(g_{re}\), while \(g_{s,re}\) represents the entropy number of degrees of freedom after reheating. The energy density at the end of inflation is denoted by \(\rho _{e}\), with \(a_0\) and \(T_0\) representing the scale factor and temperature today, respectively. The energy density above is model-dependent and can be expressed as \(\rho _e=\frac{3}{2}V_e\). Here, \(V_e\) represents the potential of the model at the end of inflation, while \(H_k\) is the Hubble function at the comoving Hubble scale wavenumber k.

An expression for the number of e-folds during reheating, in terms of energy densities, can be obtained [20] by solving the fluid equation assuming a constant EoS \(\omega _{re}\)

$$\begin{aligned} N_{re}= \frac{1}{3(1+\omega _{re})}\ln \left[ \frac{\rho _e}{\rho _{re}}\right] = \frac{1}{3(1+\omega _{re})}\ln \left[ \frac{\frac{3}{2}V_e}{\frac{\pi ^2g_{re}}{30}T_{re}^4}\right] , \end{aligned}$$
(2.8)

where \(T_{re}\) is the reheating temperature. From Eqs. (2.7) and (2.8) we get

$$\begin{aligned} N_{k}=- \frac{1-3\, \omega _{re}}{12(1+\omega _{re})}\ln \left[\frac{45V_e}{\pi ^2g_{re}T_{re}^4}\right]-\frac{1}{3} \ln \left[ \frac{11 g_{s,re}}{43}\right] -\frac{1}{4} \ln \left[ \frac{30}{\pi ^2 g_{re} } \right] -\ln \left[ \frac{ k}{a_0 T_0} \right] -\frac{1}{4}\ln \left[ \frac{27V_e M_{Pl}^4 }{2V_k^2} \right] , \end{aligned}$$
(2.9)

where \(V_k\equiv V(\phi _k)\) is the potential at the comoving Hubble scale wavenumber k. We can express the potential as \(V(\phi )=V_0 f(\phi )\), where \(V_0\) represents the overall scale and \(f(\phi )\) contains all the terms of the potential that depend on \(\phi \). This choice does not introduce any loss of generality. Equation (2.9) simplifies as follows

$$\begin{aligned} N_{k}=\ln \left( \frac{\left( \frac{43}{11g_{s,re}}\right) ^{1/3}g_{re}^{1/4}\sqrt{\pi }a_0T_0}{3\times 5^{1/4}k_p}\left( \frac{\pi ^2g_{re}T_{re}^4}{45M_{Pl}^4}\right) ^{\frac{1-3\omega _{re}}{12(1+\omega _{re})}}\left( \frac{V_0}{M_{Pl}^4}\right) ^{\frac{1+3\omega _{re}}{6(1+\omega _{re})}}\frac{f^{1/2}(\phi _k)}{f^{\frac{1}{3(1+\omega _{re})}}(\phi _e)}\right) . \end{aligned}$$
(2.10)

We can further simplify Eq. (2.10) eliminating \(V_0\). By using Eq. (2.5) we can solve for \(V_0\) in terms of \(\phi _k\) and the amplitude of scalar perturbations \(A_s\) at horizon crossing

$$\begin{aligned} V_0=\frac{12A_s\pi ^2\left( f^{\prime }(\phi _k)M_{Pl}\right) ^2}{f(\phi _k)^3}M_{Pl}^4, \end{aligned}$$
(2.11)

where \(f^{\prime }(\phi _k)\) is the derivative of \(f(\phi )\) with respect to \(\phi \) evaluated at \(\phi =\phi _k\). Finally, Eq. (2.10) can be written as

$$\begin{aligned} N_{k}=T_1+\frac{1}{3(1+\omega _{re})}\ln \left( \frac{\left( f^{\prime }(\phi _k)M_{Pl}\right) ^{1+3\omega _{re}}}{f(\phi _k)^{3 \omega _{re}}f(\phi _e)}\right) , \end{aligned}$$
(2.12)

where the term \(T_1\) is given by

$$\begin{aligned} T_1=\ln \left( \frac{\left( \frac{43}{11g_{s,re}}\right) ^{1/3}g_{re}^{1/4}\sqrt{\pi }a_0T_0}{3\times 5^{1/4}k_p}\left( \frac{\pi ^2g_{re}}{45}\right) ^{\frac{1-3\omega _{re}}{12(1+\omega _{re})}}\left( 12A_s\pi ^2\right) ^{\frac{1+3\omega _{re}}{6(1+\omega _{re})}}\right) +\frac{1-3\omega _{re}}{3(1+\omega _{re})}\ln \frac{T_{re}}{M_{Pl}}. \end{aligned}$$
(2.13)

We immediately notice that, for \(\omega _{re}=1/3\), \(T_1\) (and \(N_k\)) is independent of \(T_{re}\). Equation (2.12) is a general equation connecting reheating and inflation, valid for any single field potential and any \(\omega _{re}\). In the following section we show how to effectively use this approach by applying it to the Starobinsky potential.

3 The Starobinsky model, observables and the number of e-folds

The potential for the Starobinsky model in the Einstein frame is given by

$$\begin{aligned} V(\phi )=V_{0}\left( 1-e^{-\sqrt{\frac{2}{3}}\frac{\phi }{M_{Pl}}}\right) ^{2}. \end{aligned}$$
(3.1)

By expressing the model in the Einstein frame and identifying a scalar potential \(V(\phi )\), we can use standard expressions for slow-roll single-field inflation. This enables us to establish relationships between cosmological observables such as r, \(n_{sk}\), \(n_{tk}\), and the scalar spectral index \(n_s\). Here, \(n_{sk}\) represents the running of the spectral index, while \(n_{tk}\) denotes the tensor running.

To find relations between \(n_{s}\) and the rest of the quantities of interest, it is convenient to derive a closed-form expression for the inflaton field at horizon crossing, \(\phi =\phi _k\). This can be accomplished by solving Eq. (2.2), which results in

$$\begin{aligned} \phi _{k}=\sqrt{\frac{3}{2}}M_{Pl}\ln \left( \frac{4+3\delta _{n_s}+4\sqrt{1+3\delta _{n_s}}}{3\delta _{n_s}}\right) , \end{aligned}$$
(3.2)

where \( \delta _{n_s}\equiv 1-n_s.\) We can determine the number of e-folds during inflation \(N_k\), after the pivot scale of wavenumber \(k\equiv a_k H_k\) left the horizon, using the SR approximation

$$\begin{aligned} N_{k}=-\frac{1}{M_{Pl}^{2}}\int _{\phi _{k}}^{\phi _{e}} \frac{V}{V'} d\phi = \frac{3}{4}\left( e^{\sqrt{\frac{2}{3}}\frac{\phi _{k}}{M_{Pl}}}-e^{\sqrt{\frac{2}{3}}\frac{\phi _{e}}{M_{Pl}}}-\sqrt{\frac{2}{3}}\left( \frac{\phi _{k}}{M_{Pl}}-\frac{\phi _{e}}{M_{Pl}}\right) \right) . \end{aligned}$$
(3.3)

Where \(\phi _{e}\) is the field evaluated at the end of inflation which, following [24], we approximate as \(\phi _e \approx 0.615 M_{Pl}\). Given the horizon exit value for \(\phi _k\) by Eq. (3.2), it is possible to express \(N_k\) in terms of the spectral index \(n_s\).

At the origen the Starobinsky model is well approximated by a quadratic potential, in this case \(\omega _{re}=0\). Also, for the Starobinsky model \(T_{re}\) has been determined to be \(3.1\times 10^9 \textrm{GeV}\) [11]. In this case \(T_1\approx 58.7261+ \frac{1}{3}\ln \frac{T_{re}}{M_{Pl}}\approx 51.8988\). From Eq. (2.12) we get

$$\begin{aligned} N_{k}=T_1+\frac{1}{3}\ln \left( \frac{f^{\prime }(\phi _k)M_{Pl}}{f(\phi _e)}\right) , \end{aligned}$$
(3.4)

where \(f^{\prime }(\phi _k)=\frac{2}{M_{Pl}}\sqrt{\frac{2}{3}}e^{-\sqrt{\frac{2}{3}}\frac{\phi _k}{M_{Pl}}}\left( 1-e^{-\sqrt{\frac{2}{3}}\frac{\phi _k}{M_{Pl}}}\right) \) and, at the end of inflation, \(f(\phi _e)= \left( 1-e^{-\sqrt{\frac{2}{3}}\frac{\phi _e}{M_{Pl}}}\right) ^{2}\). Having obtained \(N_k\) in terms of \(\phi _k\), and \(\phi _k\) in terms on \(n_s\), we can solve Eq.(3.4) directly for \(n_s\) (see Fig. 1). This yields the following result

$$\begin{aligned} n_s= 0.96235, \end{aligned}$$
(3.5)

with the last digit rounded off.

Fig. 1
figure 1

The solid (blue) curve represents the lhs of Eq. (3.4) while the dashed (orange) curve represents the slowly varying function of the rhs of the same equation. For the Starobinsky model in the Einstein frame defined by Eq. (3.1) we approximate the EoS by \(\omega _{re}=0\), and \(T_{re}=3.1\times 10^9\,\textrm{GeV}\) [11]. The intersection point gives us the value \(n_s=0.96235\) for the scalar spectral index and, from the consistency relations for the Starobinsky model given by Eqs. (2.1), (3.6), (3.7), and (3.8), we obtain the other observables

Full size image

The consistency relations for the Starobinsky model [25] provide the following expressions for the tensor-to-scalar ratio r, the running of the spectral index \(n_{sk}\), and the running of the tensor \(n_{tk}\), in terms of \(n_s\)

$$\begin{aligned} r=\frac{4}{3}\left( 2-2\sqrt{1+3\delta _{n_s}}+3\delta _{n_s}\right) , \end{aligned}$$
(3.6)
$$\begin{aligned} n_{sk}=\frac{1}{18}\left( 4+9\delta _{n_s}-9\delta _{n_s}^2-\sqrt{1+3\delta _{n_s}}\left( 4+3\delta _{n_s}\right) \right) , \end{aligned}$$
(3.7)
$$\begin{aligned} n_{tk}=\frac{1}{36}\left( 8+3(4-3\delta _{n_s})\delta _{n_s}-8\sqrt{1+3\delta _{n_s}}\right) , \end{aligned}$$
(3.8)

where \(\delta _{n_s}\equiv 1-n_s\), as before. Using these equations, we can compute the values for the other observables. We can also determine the number of e-folds of inflation, reheating and, from entropy conservation after reheating [26], the number of e-folds of radiation

$$\begin{aligned} N_{rd}\equiv \ln \left( \frac{a_{eq}}{a_r}\right) =\ln \left( \frac{a_{eq} T_{re}}{\left( \frac{43}{11 g_{s,re}} \right) ^{1/3}a_0T_0}\right) , \end{aligned}$$
(3.9)

where \(a_{r}\) denotes the scale factor at the end of reheating or, equivalently, at the beginning of the radiation epoch and \(a_{eq}\) is the scale factor at radiation-matter equality. Values of observables as well as number of e-folds are given in the Table 1.

Finally, we can use the expression

$$\begin{aligned} N_{keq}\equiv \ln \left( \frac{a_{eq}}{a_k}\right) = \ln \left( \frac{a_{e}}{a_k}\right) + \ln \left( \frac{a_r}{a_e}\right) + \ln \left( \frac{a_{eq}}{a_r}\right) =N_{k}+N_{re}+N_{rd}, \end{aligned}$$
(3.10)

as a consistency check. This equation can be written more concisely as \(N_{keq}=\ln \left( \frac{a_{eq} H_k}{k_p}\right) =\ln \left( \frac{a_{eq}\pi \sqrt{A_s r}}{\sqrt{2}k_p}\right) \) [23]. Thus, the number of e-folds from the time scales of wavenumber \(k=a_k H_k\) leave the horizon at \(a_k\) to the time of radiation-matter equality at \(a_{eq}\) is essentially given by the parameter r, equivalently, by the value of the scale factor at horizon crossing \(a_k\) [27]. We find that \(N_{k}+N_{re}+N_{rd}\approx 113.182\), and the same value using the formula \(N_{keq}=\ln \left( \frac{a_{eq}\pi \sqrt{A_s r}}{\sqrt{2}k_p}\right) \).

Table 1 Above, we have listed the parameter values used in the calculations, while below, we present the values obtained for the observables and number of e-folds. Note that for \(g_{re}\approx g_{s,re}\), Eq. (3.4) is practically independent of both \(g_{re}\) and \(g_{s,re}\)
Full size table

4 Conclusions

Our approach allows us to find solutions for the cosmological observables, including the number of e-folds during inflation, reheating, and radiation, with minimal assumptions. By imposing the reheating conditions, we establish a connection between inflation and reheating and derive Eq. (2.12) which, for the Starobinsky model and \(\omega _{re}=0\), is solved for the spectral index \(n_s\). We use the consistency relations of the model to determine the values for the other observables. The number of e-folds during inflation \(N_k\) and the number of e-folds during reheating \(N_{re}\) are also determined by their respective formulas involving \(n_s\), while the number of e-folds during radiation \(N_{rd}\) is determined by the reheating temperature \(T_{re}\). The results show remarkable agreement between the Starobinsky model and current measurements of the power spectrum of primordial curvature fluctuations and the present bounds on the spectrum of primordial gravitational waves.