Starobinsky-like two-field inflation

We consider an extension of the Starobinsky model, whose parameters are functions of an extra scalar field. Our motivation is to test the robustness (or sensitivity) of the Starobinsky inflation against mixing scalaron with another (matter) scalar field. We find that the extended Starobinsky model is (classically) equivalent to the two-field inflation, with the scalar potential having a flat direction. For the sake of fully explicit calculations, we perform a numerical scan of the parameter space. Our numerical calculations support the viability of the Starobinsky-like two-field inflation for the certain range of its parameters, which are characterized by the scalar index $n_s=0.96\pm 0.01$, the tensor-to-scalar ratio $r<0.06$, and small running of the scalar index at $|\alpha_s|<0.05$.


Introduction
Cosmological inflation in the early Universe is practically well established both theoretically and experimentally. It gives the universal solution to many problems of the Standard Cosmology, because it predicts homogeneity of our Universe at large scales, its spatial flatness, its large size and entropy, as well as the almost scale-invariant spectrum of cosmological perturbations, in remarkable agreement with the COBE, WMAP, PLANCK and BICEP2 measurements of the Cosmic Microwave Background (CMB) radiation spectrum. Inflation is also thought of as the amplifier of microscopic quantum field fluctuations in vacuum, and it is the only known mechanism for seeds of the macroscopic structure formation.
The standard mechanism of inflation in field theory uses a scalar field (called inflaton), whose potential energy drives inflation. The inflaton scalar potential should be flat enough to meet the slow-roll conditions during the inflationary stage. Physical nature and fundamental origin of inflaton and its interactions to the Standard Model (SM) elementary particles are unknown.
Starobinsky inflation [1,2,3,4,5] offers the gravitational origin of inflaton by identifying it with the spin-0 part of space-time metric. In the Higgs inflation [6,7,8] inflaton is identified with the Higgs field of the SM. Both those single-field inflationary models offer the very economic and viable descriptions of chaotic inflation together with the clear origin of the inflaton field either from gravitational theory or from particle theory, respectively. As regards slow-roll inflation, the predictions of the Starobinsky and Higgs inflationary models are essentially the same (see below in this Section).
The simplest Starobinsky model of inflation is based on the modified gravity action [1] (1. 2) The inflationary model (1.1) has a single mass parameter M whose value is fixed by the observational Cosmic Microwave Background (CMB) data as M = (3.0 × 10 −6 )( 50 Ne ) where N e is the e-foldings number. The predictions of the Starobinsky model (1.1) for the spectral indices n s ≈ 1 − 2/N e ≈ 0.964, r ≈ 12/N 2 e ≈ 0.004 and low non-Gaussianity are in agreement with the WMAP and PLANCK 2013 data (r < 0.13 and r < 0.11, respectively, at 95% CL) [9], though are in disagreement with the BICEP2 measurements (r = 0.2 + 0.07, −0.05) [10]. The enhancement of the tensor-to-scalar-ratio r of the Starobinsky model to the higher values can be achieved via modification of the simplest Ansatz (1.1) by (matter) quantum corrections (beyond one loop) [11,12]. However, the Planck 2015 data [13] excludes a significant enhancement of r beyond r = 0.08. Therefore, the Starobinsky model (1.1) still perfectly fits the current observational data.
It raises the natural question on the theory side about robustness of the simplest Starobinsky model (1.1) against mixing scalaron with other (matter) scalars. Though the current observational data favors a single-field inflation, it is very unlikely that any single-field inflationary model is capable to provide the ultimate description of inflation. As regards a more fundamental description of inflation in supergravity and string theories, multi-field inflation is a must, see e.g., Refs. [14,15,16]. The direct observational evidence for multi-field inflation would be a detection of primordial isocurvature perturbations beyond the adiabatic spectrum (see Subsec. 3.3 for details).
In this paper, we study the two-field extensions of the Starobinsky model by nonminimal couplings, motivated by generic supergravity extensions of Eq. (1.1) in Ref. [17].
The action (1.1) can be dualized by the Legendre-Weyl transform [18,19] to the standard (quintessence) action of the Einstein gravity coupled to a single physical scalar (canonically normalized inflaton) φ having the scalar potential During slow-roll inflation the R 2 term dominates in the action (1.1), whereas the coupling constant in front of it is dimensionless. It implies the (approximate) rigid scaling invariance of the Starobinsky inflation in the high curvature R (or in the large field φ → +∞) limit [8]. The scaling invariance is not exact for finite values of R, and its violation is exactly measured by the slow-roll parameters, in full correspondence to the observed (nearly conformal) spectrum of the CMB perturbations. The (approximate) flatness of the inflaton scalar potential implies the (approximate) shift symmetry of the inflaton field. It also implies the alternative physical interpretation of the inflaton field as the pseudo-Nambu-Goldstone boson associated with spontaneous breaking of the scale invariance [20,21,22].
Similar observations apply to the Higgs inflation in the presence of a non-minimal coupling of the Higgs field to the space-time scalar curvature [6]. It also has the approximate (rigid) scale invariance and, actually, the same scalar potential (1.3) during slow roll inflation [8]. The Higgs inflation is based on the Lagrangian (in the Jordan frame) [6] having the real scalar field φ J (x) non-minimally coupled to gravity with the coupling constant ξ, and the Higgs scalar potential The action (1.4) can be rewritten to the Einstein frame by the Weyl transformation . (1.6) It gives rise to the standard Einstein-Hilbert term (− 1 2 R) for gravity in the Lagrangian. However, it also leads to a non-minimal (or non-canonical) kinetic term of the scalar field φ J . To get the canonical kinetic term, a scalar field redefinition is needed, φ J → ϕ(φ J ), subject to the condition As a result, the non-minimal theory (1.4) is classically equivalent to the standard (canonical) theory of the scalar field ϕ(x) minimally coupled to gravity, and having the scalar potential [6] V (1.9) Given a large positive ξ ≫ 1, one easily finds in the large field limit ϕ ≫ 2 The LHC and TEVATRON measurements of the masses of Higgs and t-quark, however, imply (via the renormalization group and the Standard Model particle content) that the effective Higgs potential coupling constant λ becomes negative at around 10 11 GeV [23], which is lower than the expected scale of inflation. It means that the SM has to be extended by new particles and new physics.
It is still possible that inflaton is neither Starobinsky scalaron nor Higgs field, but a mixture of them. This possibility leads to a two-field inflation also. Another motivation to study the Starobinsky-like two-field inflation comes from 4D, N=1 supergravity with chiral matter superfields, where inflaton is automatically extended to a complex field as the leading bosonic field component of an N=1 scalar supermultiplet. For example, as was demonstrated in Ref. [17], a generic N=1 supergravity extension of the simplest Starobinsky model (1.1) leads to the non-minimal couplings of the Higgs field to both R and R 2 gravity terms. And it is commonplace in string cosmology that inflaton is mixed with other scalars (moduli), so that a stabilization of the latter is required for inflation.
Our paper is organized as follows. In Sec. 2 we define the new class of two-field inflationary models as a combination (and a generalization) of Eqs. (1.1) and (1.4), and rewrite them to the more standard (dual) form. Those inflationary models interpolate between the Starobinsky and Higgs (single-field) inflationary models, and can accommodate a broader range of values for the tensor-to-scalar ratio. In Sec. 3 we set up the equations of motion, and classify our model against the other two-field inflationary models studied in the literature. In Sec. 4 we focus on the particular case by dropping the Higgs part of the scalar potential. In Sec. 5 we summarize our numerical findings in the special two-field model of the Starobinsky-like inflation. Sec. 6 is our Conclusion. The technical details about linear perturbations, their spectra and evolution are collected in Appendix A.

Starobinsky-and Higgs-inspired two-field inflation models
Our new inflationary model (in Jordan frame) of a real scalar field φ non-minimally coupled to the Starobinsky (R + R 2 ) gravity is given by Both functions enter the Lagrangian (2.1) via their squares, in order to avoid ghosts. It is worth mentioning that both non-minimal couplings are required by renormalization of the (R + R 2 ) gravity coupled to matter. In other words, we just replaced the parameters of the Starobinsky (R + R 2 ) gravity by functions of a (Higgs) scalar field φ.
Should the scalar field φ be stabilized by its scalar potential V (φ) to some vacuum expectation value φ 0 , our model reduces to the standard Starobinsky model (Sec. 1). Should the M 2 (φ) be sent to infinity, the Higgs inflationary setup is recovered.
In the case of the truly Higgs field φ, its scalar potential takes the form in terms of the coupling constants λ > 0 and v = φ 0 ≡ φ 0 . Thus, the model (2.1) describes all the quintessence models with a non-minimal coupling to R (like the Higgs inflation) and the R + R 2 gravity model of Starobinsky (1.1) as the particular cases. A non-minimal coupling to the R 2 term is our new feature when M (φ) is truly field-dependent.
In order to understand the physical significance of our model, and put it under the standard treatment in theoretical cosmology (in Einstein frame), let us replace the R 2 by the 2χR − χ 2 in Eq. (2.1), where the new scalar χ has been introduced as follows: It is easy to check that it is classically equivalent to the original model (2.1) because the equation of motion of the χ field is algebraic, and its solution reads χ = R.
Introducing the notation allows us to rewrite Eq. (2.3) to the Brans-Dicke-type form When assuming positivity of A (to avoid ghosts), i.e.
the Weyl transformation of metric, g µν → A(φ, χ)g µν , gives rise to the standard (Einstein-Hilbert) term for gravity in the classically equivalent (dual) Lagrangian, where we have used the Einstein-frame metric and have ignored an additive total derivative. Hence, our model takes the form of the Non-Linear Sigma-Model (NLSM) [25] having the NLSM metric (i, j = 1, 2, the primes denote the derivatives with respect to φ) in terms of the two scalars φ 1,2 = (φ, χ), minimally coupled to gravity in the physical (Einstein) frame, and having the scalar potential (2.10) The NLSM kinetic terms have no ghosts under the condition (2.6) because (2.11) The NLSM in Eqs. (2.8) and (2.9) can be further simplified by considering A as the new independent scalar field (instead of χ) and doing the field redefinition which leads to the canonical kinetic term of the scalar field ρ. We find (2.13) Hence, the full scalar potential V (φ, ρ) is given by 14) The first term of the scalar potential has a flat direction along while the second term in the case (2.2) has the absolute minimum at φ = v. Therefore, along the flat direction (2.15), the full scalar potential (2.14) has the absolute minimum at φ = v in the Minkowski vacuum. Moreover, along the flat direction (2.15), our model reduces to a single-field model having the scalar potential V (φ) in terms of the non-canonical scalar φ. The latter can be traded for a canonical scalar by a field redefinition, similarly to that of Eq. (1.7). We would like to emphasize that the discovered existence of a flat direction is automatic in the class of models under consideration, and it does not require supersymmetry. In a generic solution, two scalar fields φ and ρ are going to evolve towards the flat direction.
In order to make our model to be more specific and treatable for numerical calculations, in what follows we eliminate the functional freedom above by choosing with some (non-negative) coupling constants α, β and M . The choice (2.16) is also motivated by renormalizability. Though each of the quantum field theories (2.1) and (2.13) is not renormalizable as a theory of quantum gravity, it still makes sense to demand the (limited) renormalizability of the quantized scalar sector in a classical (curved) gravitational background. Then the non-minimal couplings (2.16) naturally arise with the renormalization counterterms [26,27]. The Higgs potential (2.2) also fits the limited renormalizability requirement.
The field theory (2.13) has two real scalars (ρ, φ) minimally coupled to the Einstein gravity and having the scalar potential (2.14). The kinetic term of the ρ scalar is canonically normalized, whereas the canonical term of the φ scalar has the ρ-dependent factor. In the next Sections we study two-field inflation in those models.
3 Classification of our model against the literature Though our model (2.13) has the non-canonical kinetic term of φ, it falls into the class of the two-field inflationary models studied, e.g., in Ref. [28] in the slow roll approximation and having the action As regards cosmological perturbations, their spectra and evolution in the theory (3.1), we employ the results of Refs. [28,29,30,31,32] in the linear approximation with respect to the slow-roll parameters. To make our paper self-contained and complete, a derivation of the relevant equations is summarized in Appendix A. Those results are used in Sec. 5 for our numerical analysis.

Equations of motion
When assuming a spatially flat Friedmann-Lemaitre-Robertson-Walker (FLRW) universe with the metric the field equations in the theory (3.1) take the form

4)
where the dots stand for the time derivatives, the subscripts after the comma denote the derivatives with respect to the fields, and H = • a /a is the Hubble function. In two-field inflation the slow-roll parameters form 2 × 2 matrices, (I, J = (ρ, φ) = 1, 2), whose entries are given by In the case (3.1) one finds (3.11) and

Correspondence of our model to the literature
In our case (2.13), we have to specify above that as well as 14) and where we have used the scalar potential (2.14). 4 Special case with V H = 0 A simple two-field inflationary model of the same type defined in our Eq. (2.1), though with a mass term instead of the Higgs scalar potential and without non-minimal interactions to R or R 2 , was considered in Ref. [32]. It was found by numerical calculations in Ref. [32] that the Starobinsky inflation is robust against that extension for a certain range of the ratio of two scalar masses. The model of Ref. [32] is in good agreement with the Planck data [13]. Multi-field dynamics of Higgs inflation in the presence of non-minimal couplings was analyzed in Ref. [33] where it was found that it is also in very good agreement with the Planck measurements.
In what follows, we restrict ourselves to the different case with V H = 0 (or in the limit λ → 0), when the functions f (φ) and M (φ) are given by Eq. (2.16), for simplicity. Then the scalar potential reads (4.1) The slow-roll matrices take the form whereas the η IJ is the same as that in Eq. (3.12) with (4.5) In particular, the Hubble function squared of Eq. (3.5) reads as (4.6)

Numerical results
The two-field scalar potential (4.1) is semi-positively definite. It reduces to the Starobinsky scalar potential (1.3) at a fixed (or stabilized) φ, and to a power-law scalar potential V /M 2 (as a sum of the α 2 φ 4 and the α 2 βφ 6 terms) at a fixed (or stabilized) ρ. The shape of the Starobinsky scalar potential at α = β = 0 is given in Fig. 1. In that case the scalar potential does not depend upon φ at all. The one-dimensional deformations of the scalar potential in the α-and β-directions are given by Fig. 2 (a) and (b), respectively.
The deformed scalar potential in the β-direction (Fig. 2(a)) at α = 0 essentially amounts to rescaling the M 2 to the M 2 (1 + βφ 2 ), though the effect of stabilization of φ is already visible in Fig. 2(a). A change of the parameter β at α = 0 merely affects the amplitude of CMB fluctuations that fixes the effective scalaron mass, M eff = (3.0 × 10 −6 )( 50 Ne ). Hence, the Starobinsky inflationary pattern is very robust against changes of β as long as βφ 2 ≪ 1 or, simply, when β is much less than 1. The situation does not significantly change under small finite values of the parameter 0 < α ≪ 1, as is illustrated by Fig. 2(b). When either α, or β, or both, grow well beyond 1, the field φ is quickly stabilized, as is illustrated by Fig. 3, (a) and (b). However, our numerical calculations show that the inflation becomes not viable by failing to get the observed value of the spectral (scalar) index n s . Therefore, in what follows, we assume that both α and β are well below 1. As our primary example, we investigate the viable inflationary model specified by the parameters α = 0.01 and β = 0.001 in more detail below. The profile of its scalar potential is given in Fig. 4.
The (time) running of the slow-roll parameters ǫ and and η in our special example is given by Fig. 5(a) and (b), respectively. The spectral scalar index at the pivot scale is given by n s = 0.96 ± 0.01. As to the tensor-to-scalar ratio r, we get r = 0.056 ± 0.003. The spectral scalar index running α s ≡ dn s /d ln k is |α s | < 0.05 in all our models.  To get those results, we used numerical solutions to the background equations of motion, whose graphs are given by Fig. 6, (a) and (b), for the fields ρ and φ, respectively.
Our numerical calculations in this Section support the qualitative conclusion that the Starobinsky inflation is robust against the field dependence in the non-minimal functions f (φ) and M (φ), as long as the non-minimal coefficients α and β are much less than 1. In other words, the Starobinsky inflation is stable against small deformations of the nonminimal couplings as long as those deformations are much less than of the order 1 (in Planck units). In the case of large deformations, inflation persists but is not viable.
We also found that at the end of inflation the scalaron field ρ oscillates near its minimum and thus contributes to (pre)heating, whereas the (matter) scalar field φ does not, approaching a constant value. It can be already seen in Fig. 6, (a) and (b), but is much better illustrated by our numerical findings in Fig. 7(a) and (b). Actually, the field φ starts oscillating and thus contributing to the reheating only when the parameter α is much larger than 1, however, it does not lead to a viable inflation.
The profile of the ρ-solution does not significantly change under varying the parameters

Conclusion
We found that the Starobinsky inflation is robust against mixing scalaron with another (matter) scalar via non-minimal interactions of the latter with both R and R 2 terms in the original (Jordan) frame, as long as the non-minimal field couplings are much smaller than one (in the Planck units). The non-minimal couplings were introduced by promoting We confirmed numerically that the inflationary trajectory in our two-field inflationary models remains close to the single-field attractor solution in the original Starobinsky model [1] under adding small non-minimal couplings to the R and R 2 terms in Eq. (2.1). Our main statement is reflected in the title of our paper by calling our two-field inflationary models the Starobinsky-like ones. Though our numerical solutions to the dynamical equations (Sec. 5) were obtained by using some initial conditions for inflation, we found that the dependence of our solutions upon small changes in the initial conditions is weak and rather unimportant. It is related to the facts that (i) our numerical solutions also exhibit an attractor-type behavior (see e.g., Refs. [36,37] for more), and (ii) our scalar potentials do not have ridges that are generically present in multi-field inflation caused by non-minimal couplings and whose presence leads to strong dependence upon the initial conditions at the onset of inflation [38].
The two-field Starobinsky-like inflation becomes not viable when any of the non-minimal parameters is of the order one or larger. Our results are complementary to the findings of Ref. [32] where the robustness of the Starobinsky inflation was established in another two-field Starobinsky-like limit with α = β = 0 and a non-vanishing mass term of the matter scalar.
The main difference of our two-field inflationary models against the single-field Starobinsky model is the presence of isocurvature perturbations. However, those perturbations turn out to be very small and (currently) undetectable. As was argued in Ref. [39], significant isocurvature perturbations in generic multi-field inflationary models with non-minimal couplings may account for the observed low power in the CMB angular power spectrum of temperature anisotropies at low multipoles [40]. However, in our models the isocurvature perturbations are not amplified enough to be the reason for that observation.
Though we did not investigate primordial non-Gaussianities in our Starobinsky-like two-field inflationary models, we expect them to be negligible, like the original (singlefield) Starobinsky model. The field-dependent couplings are quite natural from the viewpoint of string theory where all coupling constants are given by expectation values of scalar fields. As regards the physical meaning of our two scalars from the viewpoint of string theory, it is conceivable that scalaron is related to string theory dilaton, whereas another (matter) scalar is given by one of the moduli arising from superstring compactification. A detailed investigation of the possible connection of our models to string theory is beyond the scope of this paper.
The Einstein equations lead to the energy and momentum constraints as follows: In terms of the Mukhanov-Sasaki (MS) variables [34] the perturbed Klein-Gordon equations take the form where the background equations and the energy-momentum constraints above have been used in the notation [28], Adiabatic and entropy perturbations: It is common in the literature to decompose cosmological linear scalar perturbations into two directions that are either parallel or orthogonal to the trajectory in the field space [34]. The first type of perturbations is called curvature (or adiabatic) perturbations, whereas the second type is called isocurvature (or entropy) perturbations, respectively. It can be done by introducing the linear combinations (they do not refer to new scalar fields) δσ ≡ cos θδρ + sin θe b δφ and δs ≡ − sin θδρ + cos θe b δφ, (A. 15) where we have used the notation cos θ ≡ρ σ , sin θ ≡φ e ḃ σ ,σ ≡ ρ 2 + e 2bφ2 .
(A. 16) The corresponding MS variables Q σ ≡ δσ +σ H Φ are given by the linear combinations The gauge-invariant quantity, known as the co-moving curvature perturbation [34], reads in terms of Q σ as follows: while the renormalized entropy (iso-curvature) perturbation [34] is given by In terms of the adiabatic and entropy "vectors" in field space, defined by [28] E I σ = (cos θ, e −b sin θ), the corresponding first order derivatives are and the second order derivatives are Given the notation above, the background equations of motion in the adiabatic and entropy directions areσ respectively. The equations of motion for perturbations are given bÿ where we have used the notation [28] again: with the s θ and c θ standing for the sin θ and cos θ, respectively.
Perturbation spectra: We are now in a position to study the perturbation spectra of the two-field inflation and evolution of its perturbations. The power spectra of the adiabatic and entropy perturbations are given by the correlation functions [28] The cosmologically important scales are given by (i) the horizon crossing (inside the Hubble scale) and (ii) the scales over the Hubble scale, so that it is natural to evaluate the correlation functions at those scales, along the standard procedure [28].

Evolution of perturbations at the horizon crossing:
In terms of the conformal time τ = where the primes denote the derivative with respect to the conformal time τ . In the slow-roll approximation these equations can be further simplified as with the notation In Eqs. (A.38) and (A.39) we kept only the linear terms with respect to b ,ρ because it is suppressed by the M Pl and b ,ρρ = 0 in our case (2.13). The terms proportional to ξ are written down separately, in order to emphasize the difference between the canonical and non-canonical kinetic terms. According to Ref. [28], it is convenient to introduce L and Q as and then rewrite Eq. (A.37) to the standard (in mathematical physics) form As the next step, again following Ref. [28], let us introduce the time-dependent matrix R which satisfies R ′ = −LR, and the new vector v defined by u = Rv. Then the equation above can be resolved for v as in the linear order with respect to the slow-roll parameters. It follows where the second term reads .

(A.46)
It is usually assumed in the literature that the slow-roll parameters vary slowly enough during the few e-folds when the inflationary scale crosses the Hubble radius. In that case, one can replace the time-dependent matrix with its value at the Hubble crossing. In other words, the matrix on the r.h.s. of Eq. (A.46) is supposed to be evaluated at k = aH, so that the remaining time dependence only exist in the overall coefficient 3/τ 2 . Then one can always diagonalize this matrix by using a time-independent rotation matrix where the star subscript refers to the horizon crossing.
Since R varies slowly around the Hubble crossing, one can replace R by R * . When using the notation [28]λ and one can also rewrite Eq. (A.42) to The solution to Eq. (A.53) with the proper asymptotic behavior reads [28] in terms of the Hankel function H (1) µ of the first kind and of the order µ A , where [28] µ A = 9 4 + 3λ A , (A.56) and e A , A = 1, 2, are the independent orthonormal (Gaussian) random variables, Because of independence of w 1 and w 2 , the correlations of u σ and u s around the Hubble crossing can be expressed in terms of Q σ and δs as where we have Therefore, after taking into account that one finds [28] Given λ A ≪ 1, i.e. µ A ≃ 3 2 +λ A , one can further simplify the result above, by expanding G A (x) as follows: (1 + 2λ A g(x)) = (1 + x 2 )(1 + 2λ A g(x)) , (A.66) where the new function has been introduced, Evolution of perturbations on super-Hubble scales: A two-field inflationary model is reduced to a single field inflationary model when the isocurvature perturbations are suppressed. Then the adiabatic spectrum takes the form However, it does not apply to our model (Sec. 2) in a generic case where it does not reduce to the Higgs or the Starobinsky (single field) inflationary model. The existence of the isocurvature modes is a generic feature of two-field inflationary models, and it is going to affect adiabatic perturbations also during the super-Hubble scale evolution, so that Eq. (A.71) does not apply, in general. To get the power spectra and the correlation functions in that case, one should solve the coupled system of Eqs. (A.25) and (A.26). A numerical approach is the only way in most cases.
However, in some spacial cases, when the slow-roll approximation is at work, one may analytically solve the equations of motion on the super-Hubble scales too. The example considered in Ref. [28]  Unfortunately, as was already noticed in Ref. [28], the constant slow-roll approximation does not hold for many e-folds, and breaks down long before the exit from inflation. In another analytically treatable case, with the mass terms as the scalar potential and the canonical kinetic terms for scalars, the curvature and iso-curvature perturbations were computed in Ref. [35].