1. Introduction
In solving non-linear differential equations, initial conditions play a crucial role in determining the dynamical behavior of the system as different initial conditions may lead to different trajectories that may behave completely differently. Since gravity is highly non-linear in nature, solving early Universe evolution depends highly on the initial conditions as well. However, the already highly acceptable paradigm of the early Universe—the inflationary paradigm [
1,
2,
3,
4,
5,
6,
7,
8,
9,
10,
11,
12,
13], perhaps succeeds based on the fact that, not only it is in line with the recent tighter observational constraints [
14,
15] but it also provides an attractor solution [
16,
17]. This implies that the early Universe inflationary solution does not depend on the initial conditions and because of this reason, even if the Universe contains additional matter(s), inflation makes them decay exponentially and the Universe solely dominates by the inflaton field responsible for inflationary evolution.
However, the greatest achievement comes from the fact that the inflationary solution can also evade Belinsky–Khalatnikov–Lifshitz (BKL) instability: the anisotropic energy density can grow faster than the energy density responsible for the evolution of the Universe [
18] which can cause a highly unstable system. However, even with this much success, the inflationary paradigm also suffers many crucial difficulties. First, despite the ever-tightening observational constraints, there seem to exist many inflationary models that continue to remain consistent with the data [
19,
20,
21,
22], even leading to the concern whether inflation can be falsified at all [
23]. Inflation being insensitive to initial conditions is still debatable. For instance, in Ref. [
24], using non-perturbative simulations, authors showed that the inflationary expansion starts under very specific circumstances. In Ref. [
25], authors pointed out the fact that small field potential fails to start the inflationary expansion under most initial conditions. It suffers from the trans-Planckian problem: the cosmological scales that we nowadays observe may correspond to length scales smaller than the Planck length at the onset of inflation, which in contradiction to the inherent assumption that the scales are classical even at the initial stage of inflation [
26].
These issues lead to in search for alternatives to the inflationary paradigm and the most popular one is the classical non-singular bouncing scenario where, the Universe undergoes a phase of contraction until the scale factor reaches a minimum value, before it enters the expanding phase [
27,
28,
29,
30,
31,
32]. However, the problem with most bouncing models is that, it is extremely difficult to construct. While only known attractors are the ekpyrotic models [
33], in fact, in Ref. [
34], authors showed that the ekpyrotic contraction is a ‘super-smoother’, i.e., it is robust to a very wide range of initial conditions and avoid Kasner/mixmaster chaos and a similar statement could never be proven about any other primordial scenario, it fails to be in line with the observations. On the other hand, while matter bounce models may provide correct theoretical predictions consistent with the observations, the solution fails to be stable and thus, any tiny initial deviation can grow exponentially to cause an instability in the system. Another difficulty is, while constructing the contracting phase is rather easy to achieve, obtaining the non-singular bouncing phase is extremely difficult as it requires to violate the null energy condition: often called as the theoretical no-go theorem [
35,
36,
37,
38,
39,
40,
41,
42,
43,
44]. Most importantly, even if one may construct a model evading the instability and theoretical no-go theorem, perturbations suffer from many difficulties such as gradient instability, and these models fail to be in line with the observational constraints: a small tensor-to-scalar ratio
and simultaneously, very small scalar non-Gaussianity parameter
: referred to as the observational no-go theorem [
45,
46] (This has also been addressed in the past in the context of the so-called pre-big bang scenario based on the string cosmology equations, as reported for instance in a recent review in Ref. [
47]). Apart from these main three issues, in general, bouncing models possess another, rather weak, difficulty as well: a natural exit mechanism from the bouncing phase to enter into the conventional reheating phase.
In solving these problems, it is soon realized that one needs to go beyond the canonical theories and consider non-minimal couplings, viz. the Horndeski theories or even beyond Horndeski theories [
35,
36,
37,
38,
39,
40,
41,
42,
45,
46,
48,
49,
50,
51,
52,
53,
54,
54]. However, it remains an open problem until in Ref. [
55], we showed that a simple non-minimal coupling can solve all those issues.
In Refs. [
56,
57], we have shown that, the issue of stability can easily be resolved by using a simple non-minimal coupling. This has been achieved with the help of conformal transformation and we also have shown that, it may lead to viable tensor spectra as well [
58]. Recently, we extended the work to construct the first viable non-minimal bouncing model which evades all the issues and is consistent with the constraints [
55]. This model is constructed in such a way that the scalar sector of the action is conformal to that of the slow-roll inflationary model. Since, under conformal transformation, perturbations remain invariant and the constraints are obtained from the correlation functions of the scalar and tensor perturbations, we achieved to bypass observational no-go theorem. The choice of the coupling function also helps us to achieve the non-singular bouncing phase without any instability and therefore, it evades the theoretical no-go theorem and leads to the conventional reheating scenario. It has been pointed out that since the work is similar to works in Refs. [
56,
57], it can also resolve the issue of stability, mainly the BKL instability as well. However, the model contains a free model parameter
(which needs to be greater than zero to achieve the bouncing solution) and there is no bound on it. It is not clear how the system behaves in the presence of an additional matter.
At this moment, we need to stress that, when a model is studied for perturbations and is compared with the observations, it is assumed that, either the initial conditions are chosen with extreme measure, or the model solution is stable enough so that there is no need to fine-tune for choosing the initial conditions. Otherwise, even without the external matter, the perturbations themselves can cause instability as it also possesses a small amount of energy. This is the main reason why a stable solution is always preferred (one may even argue that only stable solutions are allowed). Therefore, even our model is conformal to inflationary solutions and obeys all relevant observational constraints, one needs to verify the stability of the model, especially in the presence of an external matter. This leads to the main motivation of this article.
In this work, we analyze the stability of the model in presence of an additional barotropic matter. The aim is to find the general attractor behavior of the bouncing solution that is consistent with the observations. There is mainly one free model parameter
(for minimal theory,
, since for this value, the coupling function becomes unity) along with the equation of state of the additional barotropic fluid
. In addition to that, since we consider conformal single field model that leads to slow-roll inflation, we also consider the potential slow-roll parameter in the minimal Einstein frame to be (nearly) constant (in case of exponential potential is strictly constant, as examined in Refs. [
16,
17]). With the help of these model parameters, we find the condition for stability for the bouncing solution and show that, to avoid the BKL instability leads to
. Since bouncing solution requires
, the bound becomes
. To understand the behavior of the stability, we study the phase space in the vicinity of the fixed point corresponding to the bouncing solution and find that the system acts like the Universe contains three different types of matter and due to the attractor behavior, only the leading bouncing matter component remains intact, whereas, other two decay very fast. This helps us realize how the external barotropic fluid influences the system and how the attracting nature of the bouncing solution gets rid of the initial deviations.
A few words on our conventions and notations are in order at this stage of our discussion. In this work, we work with the natural units such that , and we define the Planck mass to be . We adopt the metric signature of . We should mention that, while the Greek indices are contracted with the metric tensor , the Latin indices contracted with the Kronecker delta . Moreover, we shall denote the partial and the covariant derivatives as ∂ and ∇. The overdots and overprimes, as usual, denote derivatives with respect to the cosmic time t and the conformal time associated with the Friedmann-Lemaître-Robertson-Walker (FLRW) line-element, respectively. The sub(super)script ‘I’ and b denote the quantity in the minimal Einstein theory and non-minimal theory, respectively.
2. A Brief Introduction of the Bouncing Model
The non-minimal bouncing model is constructed in such a way that the scalar sector of the non-minimal action transforms conformally to that of a minimal canonical Einstein model that leads to slow-roll inflation. The slow-roll inflationary model can easily be constructed by using a single canonical scalar field
minimally coupled to the gravity, i.e., the Einstein gravity as
is the Ricci scalar for the metric
, and
is the potential responsible for the inflationary solution. The part of the action
indicates the presence of an additional fluid. As mentioned above, since inflation is an attractor, the evolution of the field does not depend on the additional fluid and is solely governed by the scalar sector of the system. Assuming the above theory is responsible for the slow-roll inflationary dynamics, the scale factor solution, during inflation, can approximately be written as a function of the scalar field
as
with
. Now we shall construct a model which is conformal to the scalar part of the above inflationary action (
1) in such a way that the new scale factor behaves as a bouncing solution. The transformation can be written as
and
are the required bouncing metric and scale factor solutions, respectively, along with
being the coupling function. Using the above conformal transformation along with the action (
1), we can construct the action responsible for such bouncing solution as
and the potential
depend on the coupling function
and the inflationary potential
as
where
.
is the external fluid present in the system in the non-minimal theory. Please note that, the above action resembles very similar to the gravi-dilaton string effective action, with the inclusion of non-perturbative string-coupling corrections [
47]. Therefore, one can explain the existence of such kind of coupling functions in the domain of string theory. Note that, while the scalar part of the above action (
4) is conformal to that of (
1), these two actions in reality are not conformal due to the presence of the additional fluid. However, if the scalar field dominated solution in the non-minimal theory also becomes stable, similar to minimal theory, then, for a given inflationary model (
1) with its solution (
2), one can obtain the bouncing scale factor solution
in the non-minimal theory by setting
in a certain manner. For simplicity, if we desire the bouncing solution of the form
where
is the comoving time, then by using Equations (
2) and (
3), one can obtain the required solution of the coupling function
as
is normalized in such a way that at the minima of the potential,
becomes unity. The non-minimal theory (
4) is the desired bouncing model: given an inflationary potential
, we can completely determine the coupling function
, which in turn, provides non-minimal theory (
4) responsible for the bouncing scale factor solution (
6). Note that, since Equation (
2) is an approximated form, by using the above expressions (
7) and (
3), the obtained scale factor solution
is also not exact, but close to (
6) (we will evaluate it later). We need to stress that the form of the scale factor (
6) is valid only in the contracting phase, and during and after the bounce, it is extremely difficult to obtain an analytical solution. Therefore, the bouncing model, in our case, is also asymmetrical. The underlying assumption of the theory is that, either the bouncing solution is stable or we must choose a precise initial condition that leads to the bounce. Last of all, while
provides a bouncing solution, there are actually no bound on
as it can take any value and negative values of it will lead to different Universe with different scale factor solution (
6). For instance,
brings us the conventional minimal Einstein scenario as, for this value, the coupling function becomes unity. In this work, we will concentrate only on the bouncing solutions where
and obtain the condition for stability and the behavior of the system in the vicinity of the contracting solution, which we will explore in the next section.
Few things to note before we move in the next section. In case of scalar perturbation, in the Einstein frame (
1), the perturbed action is
where,
is the curvature perturbation and ∇ is the divergence operator. In the non-minimal frame (
4), the expression of the action for the scalar perturbation changes by replacing
with
. However,
is simply the conformally transformed
, i.e.,
This implies that,
which leads to the action being identical in both frames and hence,
at linear order. As the speed of sound is conformally invariant,
is unity in both minimal inflationary as well as in non-minimal bouncing scenarios. Similarly, one can evaluate the perturbed interaction Hamiltonian (for detailed evaluation, see Refs. [
59,
60,
61,
62,
63]) at any order and show that it also remains invariant in all conformal frames. This implies that the curvature perturbations at any order in both the frames are identical, which is not surprising, as we know, under conformal transformation, curvature perturbation remains invariant. The argument extends to tensor perturbation as well. Therefore, given a viable inflationary model which is consistent with the observations, in our bouncing model, there appear no divergences or instabilities (e.g., gradient instability) in the perturbations anywhere, even at the bounce and, in addition to that, the model also satisfies all observational constraints: thus evading the no-go theorem. The model leads to stable non-singular bounce as well as, shortly after the bounce, it leads to conventional reheating phase (for detailed information, see
Figure 1 and
Figure 2, along with Ref. [
55] as well as [
64], to compare reheating in minimal theory with the same in non-minimal theory).
In short, for
, the model (
4) leads to a stable, asymmetrical bouncing Universe, where the scale factor in the contracting region behaves as
If the slow-roll inflationary potential
satisfies the observational constraints (e.g., the Starobinsky potential), the corresponding bouncing model, irrespective of the value of
, also satisfied all observational constraints, as the perturbations remain invariant under conformal transformation, thus evading the observational no-go theorem.
3. Stability Analysis
Let us now analyze the stability condition for the above system (
4). Assuming the additional fluid is barotropic, in nature, with the equation of state
, the governing equations become
where,
is the energy-momentum tensor of the additional barotropic fluid in the non-minimal theory (
4), satisfying
We have not used any subscript
b since, from now on, all concerning equations will correspond only to the non-minimal theory. These equations, using Friedmann–Robertson–Walker (FRW) homogeneous and isotropic Universe with the line element
turn into
where
and
are the energy density and the pressure of the barotropic fluid with
being the equation of state. To study the stability analysis, instead of using quantities with dimensions, we can express and re-write these field equations in terms of dimensionless quantities. These are defined as
Using these dimensionless variables, one can quickly express the evolution equation of them as
along with the constraint equation
N is defined as the e-folding number in the non-minimal frame
is the initial value of the scale factor and
is the potential in the minimal Einstein frame and is related to
through Equation (
5). Note that, the model parameter
actually represents a parameter in the Einstein frame which is completely arbitrary and can be fixed from the observations. For example, while near the pivot scale, for chaotic inflation,
, for Starobinsky model, at the same scale, it becomes
. However, chaotic inflation has already been ruled out from the observations while, the Starobinsky model remains consistent with it and one of the most important models in the paradigm. We should stress that the phase during the stability is analyzed is the contraction in a bouncing model and therefore, instead of the e-N-fold time convention, here we are using the e-fold time convention.
One can also express the slow-roll parameter
in terms of these dimensionless parameters as
Using the above expression, one can define the effective equation of state as
Although is a function of and thus, it is dynamical, the minimal theory that we are considering leads to slow-roll inflation and thus one can consider it to be nearly constant, similar to exponential potential.
In this system, as one can see, there exist three model parameters: . is arbitrary as it describes the nature of the additional field. For example, in the case of anisotropic stress fluid that corresponds to the stiff matter, . This field results in the BKL instability that, later, we will investigate. As mentioned above, is nearly constant and hence, we are left with one model parameter that can be constrained. In the next section, we will try to constrain it using the Lyapunov exponents.
3.1. Fixed Points and the Bouncing Solution
Using the evolution equation, one can find the fixed/critical point of the system and these points correspond to the solutions of the system, in this case, the Universe depicted by the non-minimal theory (
4). These can be found by setting Equations (
18) and (
19) to be equal to zero, i.e., the velocities of
x and
y vanishes at these points. There are seven such fixed points:
The sign of the
signifies whether the Universe is expanding or contracting. Since,
H appears in the denominator in the expression of
y in (
17), the negative sign of
y leads to contraction (bouncing), while positive
y corresponds to expanding Universe. Therefore,
leads to contraction, and
implies the expanding Universe. In fact,
is our desired solution and the scalar field dominated solution as the fractional energy density of the additional fluid and the effective equation of state (see Equation (
22)) of corresponds to this point are
which corresponds to the the scale factor solution as
As you can see, the exponent is not exactly
but
, which is contrary to the Equation (
6). However, as mentioned before,
for slow-roll models and only by using such model, we achieved such non-minimal bounce, and therefore,
. Therefore,
corresponds to our bouncing solution that we will investigate thoroughly in the next section. Other fixed points are
Since our desired fixed point is , in this work, we will focus only on this point. In the next subsection, we will explore the Lyapunov exponents that determine the stability of the solution, and then later, we will focus on the behavior of the Universe close to the fixed point.
3.2. Lyapunov Exponent and the Stability of the Bouncing Solution
Now, in order to study the stability of these fixed points, we need to linearize the Equations (
18) and (
19) as
where,
and
are the right-hand side of (
18) and (
19), respectively.
denotes the value at the fixed point.
and
are the deviations of the fixed points. By linearizing equations, we assume that we are studying the stability condition in the vicinity of the fixed points, i.e., the deviation from the fixed points are very small. We are considering homogeneous and isotropic background and assume that the anisotropies and inhomogeneities caused by the deviations do not impact on the background. This can be achieved by assuming that the energy fractions due to the deviations are extremely small.
In order to check the stability of the fixed points, we need to evaluate the eigenvalues of this matrix. These eigenvalues are the Lyapunov exponents that arise in the evolution of both
and
and thus, only these values determine the stability of the fixed points as
As defined earlier, N is the e-folding number; in an expanding Universe, while it increases, it decreases for contracting solution, and therefore, if both the eigenvalues are negative (positive) in an expanding (contracting) Universe, then and approach zero as N approaches . This implies that the deviation from the actual trajectory reduces with time and the new trajectory returns to the original trajectory asymptotically in time. In other words, both eigenvalues being negative (positive) in an expanding (contracting) Universe implies that the fixed point is stable and the corresponding solution is an attractor solution.
Let us now calculate the Lyapunov exponents for our desired fixed point (
23). It takes the form
The above expression bears one of the main results of this work. As one can see properly since the desired fixed point is the bouncing solution, we require both of the
’s to be positive. Few things to note: since
represents the minimal solution as mentioned before (coupling function (
7) becomes unity), one can easily obtain the results in Refs. [
16,
17] by substituting
. Since
and
are positive, increasing values of
and
narrow the permissible regime of
for being an attractor. Anisotropic fluid or the stress fluid resembles similar to a stiff matter with the equation of state
and the energy density
, where
a is the scale factor. Therefore, since BKL instability is caused by it, the condition to evade the BKL instability can easily be obtained by substituting
:
Therefore, for exponential potential,
is constant and positive domain of
can only be obtained if
. In the case of slow-roll,
and therefore, for slow-roll inflationary models which are conformal to the non-minimal bounce, the condition becomes
as the bouncing solution can only be obtained for
. This is the main result of this article as, now,
cannot possess arbitrary positive values for a stable viable solution. In the next section, in order to understand how the deviations impact the system near the fixed point, we will study the phase space in the concerned region.
4. Evolution of the Deviations
Using the eigenvalues and the eigenvectors, one can easily solve (
33). For the fixed point (
23), it becomes
where,
and
are given in (
34). The
A’s can be evaluated as
These solutions are obtained by using the initial conditions
. The slow-roll parameter can also be obtained in the vicinity of the fixed point (
23) by perturbing the Equation (
21) and this becomes
where,
Once we obtain the solution expression for the slow-roll parameter, we can solve the energy density equation as
and by performing the integral, one can easily solve for
H. It becomes,
is the initial value of the Hubble parameter. Using (
37) and (
43), and by using the fact that,
and
are tiny and act like perturbations, we obtain the evolution of the Hubble parameter as
where,
One can easily obtain the the above expressions in terms of model parameters
by using relations (
39)–(42), (
44) and (45). It is clear from the above equation that the system behaves as the Universe contains three fluids, where the leading solution corresponds the scalar field dominated solution, in the case of additional fluids having the initial energy density fraction
and
, the energy density of them grow/decay as
and
. It is apparent that, relative to the leading order term with
, the other two decay/grow with the factor
and
, respectively. Therefore, for a contracting solution,
’s being positive signifies that the those additional fluids represented by the
and
decay. This can even be seen from the expression of the energy fraction of the additional fluid (
20). This becomes
From the above expression, it becomes obvious that, for a bouncing solution, if ’s are positive, the energy fraction vanishes after a while. To understand the correct nature of it, below, we presented three different scenarios with different choices of model parameters that try to help us realize the system.
4.1. Radiation Bounce
As we know from (
36),
must be positive and less than 2 in order to satisfy bouncing solution as well as to evade the BKL instability. Therefore, in the first example, we consider radiation bounce where
, i.e., the scale factor solution approximately becomes
. We consider the additional fluid as the anisotropic stress fluid with stiff equation of state
. Assuming typical conformal slow-roll model (e.g., chaotic inflation with
) with
, the Hubble solution takes the form
with
The equation represents the attractor nature of the solution as the leading order solution, i.e., the effective radiation matter solution dominates over the additional fluid evolution. The actual energy density fraction of the additional fluid (
49) takes the form
For bouncing solution, as it is obvious, it quickly vanishes after a while.
4.2. Comparing with the Minimal Scenario
As it has been mentioned before, negative
represents the expanding solution, and
reduces to the minimal Einstein inflationary scenario as for this value, the coupling function becomes unity. In this subsection, we will compare the above result with the minimal inflationary case. Again using the similar value of the additional fluid equation of state as well the
, i.e.,
, the Hubble solution becomes
with
Notice that, in this case, the additional fluids with
and
initial mass fractions decay faster compared to the radiation bounce scenario as the Lyapunov exponents, i.e., the
’s are having high values. This can even be seen from the expression of energy fraction of the additional fluid (
49):
It vanishes faster than that of the radiation bounce case. Therefore, one may consider, by judging the capability of stability solution, minimal inflationary models are preferred compared to the non-minimal bouncing models. In the next subsection, we will show that, this is clearly not the case.
4.3. Comparing with the Viable Ekpyrosis
In Ref. [
34], it has been studied extensively and showed that the ekpyrosis is a super-attractor (although it has been studied for minimal models, one can easily extend it for non-minimal solution as well. See Ref. [
56] in this context). It can also been seen from (
34) as
appears in the denominator in the expressions. This implies that the smaller values of
leads to higher values of Lyapunov exponents and this, in return, implies the solutions are, in comparison, more stable. Consider the example of
: the solution leads to the ekpyrosis solution. By using similar values of
and
, one can obtain the Hubble solution for the ekpyrosis as
with
Again, compared to the leading solution, the externals fluids decay much faster than the minimal case, shown above. This is because the Lyapunov exponents are, in this case, 38 and 30, respectively. For instance, the energy fraction of the additional fluid becomes
Therefore, even if all the examples, that are presented in this work, are in line with the observations, there is a difference in the behavior of the stability and among them, ekpyrosis is the most stable, and in this sense, the ekpyrosis is the most preferred, which again, has been re-established.
5. Conclusions
In this article, we construct a non-minimal bouncing model that is conformal to the scalar sector of the minimal inflationary model. However, along with the presence of the additional matter in the system, the minimal and non-minimal models do not conform and therefore, the stability of the two systems is not related by the conformal factor. is the main model parameter and by setting it to be equal to , one can arrive at the stability conditions of the minimal model as well as it had been studied before in the literature. While, the scalar part of the non-minimal bounce is conformal to that of the minimal inflationary model, choosing a viable inflationary model that satisfies the observational constraints, the non-minimal bouncing model also satisfies the observation constraints with the inherent assumption that the model solution is an attractor. If the solution is not stable, then even the quantum fluctuations can grow and diverge and the system becomes unstable. In this article, therefore, we study the stability of the non-minimal bouncing model with the presence of an additional fluid with the equation of state . We show that, while being positive leads to a viable bouncing solution, not all positive values of it leads to a stable solution, and to avoid BKL instability, we obtained the constraint in as the domain for stability as well as bounce as We also showed that the ekpyrosis is the best attractor among them and thus is always preferred.
While we studied the system for homogeneous and isotropic background, one can easily extend the work for homogeneous and anisotropic systems as well and may arrive at a similar condition. One can even go beyond the homogeneous system and study the stability using numerical relativity. However, this is beyond the scope of this article. We did not consider the negative values of
as well as
case. The last one signifies the emergent gravity scenario. Apart from that, we also did not study the behavior of the other fixed points (
28)–(32). While in the minimal scenario, most of them are unstable, in our case, it may not be the same and may impact the system heavily. We reserve these works for our future prospects.