<figure><img src="image_1.png" /></figure>

# Tracing primordial black holes in nonsingular bouncing cosmology

Jie-Wen Chen, Hao-Lan Xu, Yi-Fu Cai

► To cite this version:

Jie-Wen Chen, Hao-Lan Xu, Yi-Fu Cai. Tracing primordial black holes in nonsingular bouncing cosmology. Physics Letters B, 2017, 111, 10.1016/j.physletb.2017.03.036. hal-01494881

HAL Id: hal-01494881

https://hal.sorbonne-universite.fr/hal-01494881

Submitted on 24 Mar 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /></figure>

ELSEVIER

Contents lists available at [ScienceDirect](https://www.sciencedirect.com)

Physics Letters B

[www.elsevier.com/locate/physletb](https://www.elsevier.com/locate/physletb)

<figure><img src="image_4.png" /></figure>

# Tracing primordial black holes in nonsingular bouncing cosmology

Jie-Wen Chen<sup>a</sup>, Hao-Lan Xu<sup>a,b</sup>, Yi-Fu Cai<sup>a</sup>

## ARTICLE INFO

**Article history:**

Received 18 September 2016

Received in revised form 18 March 2017

Accepted 18 March 2017

Available online xxxx

Editor: M. Trodden

## ABSTRACT

We in this paper investigate the formation and evolution of primordial black holes (PBHs) in nonsingular bouncing cosmologies. We discuss the formation of PBH in the contracting phase and calculate the PBH abundance as a function of the sound speed and Hubble parameter. Afterwards, by taking into account the subsequent PBH evolution during the bouncing phase, we derive the density of PBHs and their Hawking radiation. Our analysis shows that nonsingular bounce models can be constrained from the backreaction of PBHs.

© 2017 Published by Elsevier B.V. This is an open access article under the CC BY license ([http://creativecommons.org/licenses/by/4.0/](http://creativecommons.org/licenses/by/4.0/)). Funded by SCOAP³.

## 1. Introduction

The matter bounce scenario [1–3] is one type of nonsingular bounce cosmology [4–9], which is often viewed as an important alternative to the standard inflationary paradigm [10–13]. By suggesting that the universe was initially in a contracting phase dominated by dust-like fluid (with a vanishing equation-of-state parameter $w = 0$), then experienced a phase of nonsingular bounce, and afterwards entered a regular phase of thermal expansion. The matter bounce cosmology can solve the horizon problem as successfully as inflation and match with the observed hot big bang history smoothly. Based on primordial fluctuations generated during matter contracting and their evolution through the nonsingular bounce, one can obtain a scale invariant power spectrum of cosmological perturbations. Unlike inflation, the matter bounce scenario does not need a strong constrain on the flatness of the potential of the primordial scalar field that drives the evolution of the background spacetime [14,15]. Also, this scenario can avoid the initial singularity problem and the trans-Planckian problem, which exists in inflationary and hot big bang cosmologies [16,17].

The aforementioned scenario has been extensively studied in the literature, such as the quintom bounce [18,19], the Lee–Wick bounce [20], the Horava–Lifshitz gravity bounce [21–23], the $f(T)$ teleparallel bounce [24–26], the ghost condensate bounce [27], the Galileon bounce [28,29], the matter-ekpyrotic bounce [30–32], the fermionic bounce [33,34], etc. (see, e.g. Refs. [1,35] for recent re-

views). In general, it was demonstrated that on length scales larger than the time scale of the bouncing phase, both the amplitude and the shape of the power spectrum of primordial curvature perturbations can remain unchanged through the bouncing point due to a no-go theorem [36,37]. A challenge that the matter bounce cosmology has to address is how to obtain a slightly red tilt on the nearly scale invariant primordial power spectrum. To address this issue, a generalized matter bounce scenario, which is dubbed as the Λ-Cold-Dark-Matter (ΛCDM) bounce, was proposed in [38] and predicted an observational signature of a positive running of the scalar spectral index [39,40].

As a candidate describing the very early universe, the matter bounce scenario is expected to be consistent with current cosmological observations and to be distinguishable from the experimental predictions of cosmic inflation as well as other paradigms [7, 41]. Meanwhile, a possible probe of primordial black holes (PBHs) may offer a promising observational approach to distinguish various paradigms of the very early universe [42,43]. PBHs could form at very early times of the universe, where a large amplitude of density perturbations would have obtained. Correspondingly, the formation process and the abundance of PBHs strongly depend on those early universe models, in which fluctuations of matter fields are responsible for such large amplitudes of density perturbations [44].

In the literature, most of attentions were paid on the computation of PBH predictions from the inflationary paradigm (for instance see [45–49]), while so far, only a few works addressed the PBH formation in a bouncing scenario [50,51]. Furthermore, those studies of PBHs in a bouncing scenario have not yet been discussed in detail, for specific cosmological paradigms or been applied to

<footnote><sup>a</sup> CAS Key Laboratory for Researches in Galaxies and Cosmology, Department of Astronomy, University of Science and Technology of China, Hefei, Anhui 230026, China</footnote>

<footnote><sup>b</sup> Institut d'Astrophysique de Paris, UMR 7095-CNRS, Université Pierre et Marie Curie, 98 bis boulevard Arago, 75014 Paris, France</footnote>

<footnote>E-mail addresses: chjw@mail.ustc.edu.cn (J.-W. Chen), xhl1995@mail.ustc.edu.cn (H.-L. Xu), yifucai@ustc.edu.cn (Y.-F. Cai).</footnote>

<footnote>http://dx.doi.org/10.1016/j.physletb.2017.03.036</footnote>

<footnote>0370-2693/© 2017 Published by Elsevier B.V. This is an open access article under the CC BY license ([http://creativecommons.org/licenses/by/4.0/](http://creativecommons.org/licenses/by/4.0/)). Funded by SCOAP³.</footnote>
falsify various early universe cosmologies, especially the matter bounce scenario. In the context of matter bounce cosmology, there are several differences on the computation of the PBH abundance comparing with that in an expanding universe. First, comparing with inflation where the primordial fluctuations become frozen at the moment of the Hubble exit, those primordial fluctuations on matter fields in bounce cosmology would continue to increase after the Hubble exit during the contracting phase until the universe arrive at the bouncing phase [9,20], and the contracting phase would yield a different initial condition for the PBH formation and evolution. Second, once these PBHs have formed, the contraction of spacetime could also compress and enlarge the primordial matter density, thus change the PBH horizon radius which then can lead to effects on their evolution.

In this paper, we perform a detailed survey on the PBH formation and evolution in the background of the matter bounce cosmology. In Section 2, we briefly introduce the matter bounce scenario and describe the formation of the power spectrum of primordial curvature perturbation in an almost model-independent framework. In Section 3, a physical picture of the PBH formation in the contracting background is presented. After a process of detailed calculations, the threshold for forming PBHs and the corresponding mass fraction are provided. In Section 4, we discuss the evolution of PBHs in the bouncing phase by taking into account the effects arisen from the contraction of the background and the Hawking radiation. In Section 5, we summarize our results and discuss on some outlook of the PBH physics within the nonsingular bouncing cosmology.

## 2. Nonsingular bounce

Nonsingular bounce can be achieved in various theoretical models, namely, to modify the gravitational sector beyond Einstein, to utilize matter fields violating the Null Energy Condition (NEC), or in the background of non-flat geometries (see e.g. [52,53]). It is interesting to notice that, in general, on length scales larger than the time scale of the nonsingular bouncing phase, primordial cosmological perturbations remain almost unchanged throughout the bounce [36,37]. In this regard, one expects that the effective field theory approach should be efficient to describe the information of a nonsingular bounce model at background and perturbation level. Recently, it was found in [30] that a nonsingular bounce model can be achieved under the help of scalar field with a Horndeski-type, non-standard kinetic term and a negative exponential potential. Within this model construction, the matter contracting phase can be obtained directly by including the dust-like fluid or involving a second matter field [31]. Note that, in the present study we assume that the effective field approach of bouncing cosmology is valid through the whole evolution without modifications to General Relativity.

### 2.1. The model

It turns out that, under the description of the effective field theory approach, the background dynamics of the nonsingular bouncing cosmology can be roughly separated into three phases: the matter-dominated contraction, the non-singular bounce, and the thermal expansion. We consider a simple model starting with a matter contracting phase ($t < t_-$) from an initial time $t_{initial} \to -\infty$, and then entering into a nonsingular bouncing phase at $t_-$, which lasts till $t_+$. After the bounce ends at $t_+$, the universe begins the hot big bang expansion, which is in accordance to the current observations.

The evolution of the matter bounce cosmology in each stage can be approximately described as follows.

(i) In the matter contracting phase, the scale factor of the universe shrinks as

$$a(t) = a_{-} \left( \frac{t - \tilde{t}_{-}}{t_{-} - \tilde{t}_{-}} \right)^{2/3}, \qquad (2.1)$$

where $a_{-}$ is the scale factor at time $t_{-}$, and $\tilde{t}_{-}$ is related to the Hubble parameter at $t_{-}$ via the relation $t_{-} - \tilde{t}_{-} = \frac{2}{3H_{-}}$. The equation-of-state parameter during this phase is $w = 0$, which can be realized in many ways, such as by cold dust, by massive field or by the gravity sector involving non-minimal couplings. We parameterize these different mechanisms by introducing the sound speed $c_s$, which can affect the propagation of primordial perturbations in the gradient terms.

(ii) In the nonsingular bouncing phase, the scale factor of the universe can be approximately described as [30]

$$a(t) = a_B e^{-\frac{\gamma t^2}{2}}, \qquad (2.2)$$

where the coefficient $a_B$ is the scale factor exactly at the bouncing point, and from Eq. (2.2) one obtains $a_{-} = a_B \exp[\gamma t^2/2]$. $\gamma$ is a model parameter describing the slope of Hubble parameter to time during the bouncing phase, as:

$$H(t) = \gamma t. \qquad (2.3)$$

It can be seen that $t=0$ corresponds to the bouncing point when the universe stops the contraction and starts the expansion. Thus, it can be found that $t_{-} = H_{-}/\gamma$, and the value of $H$ vanishes at $t=0$ which is at the bouncing point. Inserting Eqs. (2.2) and (2.3) into the Friedmann equation, one can see that the null energy condition $\rho + p > 0$ is violated around the bouncing point. It is the negative pressure that avoids the singularity and drives the universe to evolve from a contracting phase to an expanding phase. (iii) In the era of radiation-dominated expansion, we have

$$a(t) = a_{+} \left( \frac{t - \tilde{t}_{+}}{t_{+} - \tilde{t}_{+}} \right)^{1/2}, \qquad (2.4)$$

where $t_{+} = H_{+}/\gamma$, $t_{+} - \tilde{t}_{+} = \frac{1}{2H_{+}}$ and $a_{+} = a_B e^{\frac{\gamma t_{+}^2}{2}}$. In present analysis we have adopted the assumption that the heating process happens instantly after the bounce (see [54,55] for relevant analyses).

From such parameterizations, the matter bounce cosmology can be approximately described by model parameters $H_{-}$, $H_{+}$, $\gamma$, and also $c_s$ if perturbations are taken into account. In Fig. 1 we depict the evolution of comoving Hubble length $|H^{-1}| = |aH|^{-1}$, and one can read that one Fourier mode of cosmological perturbation in the matter bounce cosmology could exit the Hubble radius during the contracting phase, and then enter and re-exit the Hubble radius during the bouncing phase, and eventually re-enter the Hubble radius again in the classical Big Bang era. Note that the change of the Hubble radius in the vicinity of the bouncing point can be very large. In the literature, observational constraints upon the bounce cosmology can be derived from various cosmological experiments such as the cosmic microwave background [41] and primordial magnetic fields [56]. In this work we will provide the independent constraints on the model parameters from PBHs.

### 2.2. Curvature perturbation during matter contracting phase

During the matter contracting phase, equation of motion for the curvature perturbation can be expressed as

$$v''_{k} + (c_{s}^{2} k^{2} - \frac{z''}{z}) v_{k} = 0, \qquad (2.5)$$
<figure><img src="image_5.png" /><figcaption>**Fig. 1.** Sketches of the evolutions of the comoving Hubble radius (green curve) $|\mathcal{H}^{-1}| = |aH|^{-1}$ and an effective length scale from the coefficient $|z''/z|^{-1/2}$ (red dashed) where $z = \frac{a\sqrt{\rho+\bar{p}}}{H}$, in the nonsingular bounce cosmology. The black curve is the comoving length for wavenumber $k$. The horizontal axis corresponds to the comoving time which is defined by $\eta \equiv \int \frac{dt}{H}$. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

where $v_k = zR_k$ is the Mukhanov-Sasaki variable with $R_k$ being the comoving curvature perturbation, $z = \frac{a\sqrt{\rho+\bar{p}}}{c_sH}$ relies on the detailed evolution of the background dynamics, and the prime represents the derivative with respect to the comoving time $\eta$. Assuming that primordial perturbations originated from vacuum fluctuations at initial times,<fnref n="1" /> one can derive the solution of Eq. (2.5)

$$v_k = \frac{\sqrt{\pi(-\eta)}}{2} H_{3/2}^{(1)}[c_s k(-\eta)], \quad (2.6)$$

where $H_{3/2}^{(1)}$ is the $\frac{3}{2}$-order Hankel function of the first kind. From Eq. (2.6), the power spectrum of $R_k$ is then given by

$$\Delta_R^2 \equiv \frac{k^3}{2\pi^2} \left|\frac{v_k}{z}\right|^2 = \frac{c_s^2 k^3 (-\eta)}{24\pi M_P^2 a^2} \left|H_{3/2}^{(1)}(-c_s k \eta)\right|^2, \quad (2.7)$$

during the matter contracting phase. At large scales $c_s k \ll |\mathcal{H}|$, $\Delta_R^2 \simeq \frac{a^2 H^2}{48\pi^2 c_s M_P^2 a^3}$ is almost scale independent, and $R_k \simeq \sqrt{\frac{a^2 H^2}{24k^3 c_s M_P^2} a^{-3/2}}$. This result is different from that in an expanding spacetime, in which $R_k$ is time independent at large scales. Moreover, at small scales $c_s k \gg |\mathcal{H}|$, curvature perturbation is of quantum fluctuation with $\Delta_R^2 \simeq \frac{a^2 H^2}{12\pi^2 c_s M_P^2 a^3} \left(\frac{c_s k}{H}\right)^2$.

### 3. PBH formation during matter contracting phase

PBHs originate from the collapsed over-dense regions seeded by cosmological perturbations in early universe. When an over-dense region starts to form a black hole (BH), its radius $R$ is expected to be larger than the Jeans radius $R_J \equiv \frac{c_s}{2} \sqrt{\frac{\pi}{\bar{\rho}}} = \sqrt{\frac{2}{3}} \pi c_s |H^{-1}|$ and smaller than the Hubble length $|H^{-1}|$ [42,44,48,58,59], which implies,

$$\sqrt{\frac{2}{3}} \pi c_s |H^{-1}| \le R \le |H^{-1}|. \quad (3.1)$$

For scales smaller than $R_J$, the pressure gradient force would prevent the collapsing of the over-dense region. It is clear that $R_J \ge$

$|H^{-1}|$ when $c_s \ge 0.39$, and hence, there would be no PBH in this case.

We note that the PBH formation in the contracting phase is different from that in a purely expanding universe. In the expanding universe, the PBH formation would begin only when the size of a region reenters the Hubble scale. However, in the contracting phase, when the size of the region exits the Jeans scale, the PBH formation would have started if this region is dense enough. Thus, we would like to focus on the fluctuations within the Jeans scale instead of the Hubble scale in the contracting universe.

#### 3.1. PBH mass fraction

The background dust-like fluid inside the over-dense regions may collapse into PBHs, and the abundance of PBHs is depicted by the mass fraction, which in the matter contracting phase is

$$\beta = \frac{\rho_{\text{PBH}}}{\rho_{\text{PBH}} + \rho_{\text{bg}}}, \quad (3.2)$$

where $\rho_{\text{PBH}}$ denotes the density of PBHs and $\rho_{\text{bg}}$ the density of the remnant dust-like fluid. In general, PBHs and remnant dust-like fluid together drive the evolution of background $\rho_{\text{PBH}} + \rho_{\text{bg}} = 3M_P^2 H^2$. Note that this relation only holds when the subsequent evolutions, including the accretion and Hawking radiation, of the formed PBHs are not considered.

The mass fraction can be obtained from the Press-Schechter theory [60] (see [42,44,48] also):

$$\begin{aligned} \beta(t) &= \int_{\delta_c}^{\delta_m} \frac{2}{\sqrt{2\pi}\sigma} \exp(-\frac{\delta^2}{2\sigma^2}) d\delta \\ &= \operatorname{erfc}\left(\frac{\delta_c}{\sqrt{2}\sigma(t)}\right) - \operatorname{erfc}\left(\frac{\delta_m}{\sqrt{2}\sigma(t)}\right), \end{aligned} \quad (3.3)$$

where $\delta \equiv \frac{\delta\rho}{\bar{\rho}}$ is the fractional density fluctuations, $\sigma \equiv \sqrt{\langle\delta^2\rangle}$ is the mean mass fluctuation at the Jeans scale, which can be derived from the power spectrum of curvature perturbation $\Delta_R^2$, $\delta_c$ is the threshold of the PBH formation, and the upper limit $\delta_m$ ensures that the PBHs are no larger than the Hubble scale [42,44,48]. Note that Eq. (3.3) is based on the assumption that $\delta$ obeys the Gaussian distribution $N[0, \sigma]$.

The value of threshold $\delta_c$ is determined as follows. Considering a spherical over-dense region with a radius $R$, the space inside the region satisfies Friedmann equation [61,62]

$$(H + \delta H)^2 = H^2 (1 + \delta) - \frac{\delta K}{a^2} = H^2 (1 + \tilde{\delta}), \quad (3.4)$$

where $H$ is the Hubble parameter of the background, $a$ is the scale factor inside the over-dense region, $\delta H$ and $\delta K$ are the perturbed Hubble parameter and curvature respectively, and $\tilde{\delta} \equiv \delta - \frac{\delta K}{a^2 H^2}$. It is noticed that $H < 0$ for a contracting phase, and $\delta H < 0$ in the over dense region. The outer region ($> R$) is thought to be unperturbed for simplicity. We assume that when the region collapses to a BH, its surface has an additional physical speed $v = 1$ with respect to the conformally static background, i.e. collapsing at the speed of light. In the comoving slicing, it is written as

$$\delta H \cdot R = -1. \quad (3.5)$$

Inserting Eq. (3.4) into the above equation, one obtains the threshold for an arbitrary scale $R$

$$\tilde{\delta}_c = \delta_c - \frac{\delta K_c}{a^2 H^2} = \frac{1}{H^2 R^2} + \frac{2}{|H|R}. \quad (3.6)$$

On the Jeans scale $R_J$, one has $\tilde{\delta}_c = \frac{\sqrt{6}}{\pi c_s} + \frac{3}{2\pi^2 c_s^2}$.

<footnote n="1">¹ Note that the vacuum state is not necessarily the only choice of the initial condition for primordial curvature perturbations, namely, they may also arise from fluctuations of a thermal ensemble [20,57].</footnote>
The mass of BH with a radius R, in the contracting background, can be determined by

$$M(R) = \frac{4\pi}{3} R^3 \bar{\rho} (1 + \delta_c) = \frac{R}{2G} (1 + |H|R)^2 . \quad (3.7)$$

Note that, both the density and curvature fluctuations could contribute to the BH mass, as shown in [48]. We mention that the above description of BH is a rough estimate, and hence, Eq. (3.7) slightly differs from the Misner–Sharp mass $M = \frac{R}{2G}$. The accurate description of the BH in the matter contracting phase should be based on the Tolman–Bondi–Lemaitre metric [63–65], which will be addressed in our follow-up study. In present analysis, however, we note that at small scale $R \ll |H^{-1}|$ the above model naturally recovers the solution of a Schwarzschild BH $M = \frac{R}{2G}$ in a flat spacetime, which indicates that this estimate remains reliable.

To derive the values of $\delta_c$ and $\sigma$, one needs to know the relations among the variables $\delta$, $\delta K$ and $R_k$. It is convenient to take the comoving gauge, in which

$$H^2 \delta = \frac{2}{3} \nabla^2 \Psi(x, t), \quad \frac{\delta K}{a^2} = -\frac{2}{3} \nabla^2 \mathcal{R}(x, t), \quad (3.8)$$

where $\Psi(x, t)$ is the Bardeen potential and $\mathcal{R}(x, t)$ is the curvature perturbation [45,66,67]. The relation between their Fourier components $\Psi_k$ and $R_k$ is given by [45,67]

$$-(1+w)R_k = \frac{5+3w}{3}\Psi_k + \frac{2\Psi_k}{3H}, \quad (3.9)$$

where the dot denotes the derivative with respect to the cosmic time $t$. The solution during the matter contracting era $R_k \propto a^{-3/2}$ is

$$\Psi_k = -\frac{3}{2}R_k, \quad (3.10)$$

differing from $\Psi_k = -\frac{3}{5}R_k$ in the matter dominant expanding era. One also obtains

$$\delta_k = \frac{k^2}{a^2 H^2} R_k, \quad (3.11)$$

where $\delta_k$ is the Fourier component of the density perturbation $\delta$. From (3.6) and (3.8), one also has $\delta = \frac{3\delta K}{2a^2 H^2} = 3\delta$. Therefore, the threshold is $\delta_c = \frac{3\sqrt{6}}{\pi c_s} + \frac{9}{2\pi^2 c_s^2}$.

The value of the upper limit $\delta_m$ for the Jeans scale is fixed as follows. If $\delta = \delta_m$, the fluid inside the Hubble radius would form a BH,<fnref n="2" /> which would have

$$\frac{4\pi}{3} R_J^3 \bar{\rho} (1 + \delta_m) + \frac{4\pi}{3} \bar{\rho} (|H^{-3}| - R_J^3) = \frac{2|H^{-1}|}{G}, \quad (3.12)$$

where the first term on the left side is the mass of the overdense region and the second term is the mass outside the $R_J$ but inside the Hubble radius. As a result, we obtain $\delta_m = 3\delta_m = 9(\sqrt{\frac{2}{3}}\pi c_s)^{-3}$.

From Eq. (3.11), one can get the power spectrum of the density perturbation

$$\Delta_{\delta}^{2} = \frac{k^3}{2\pi^2} |\delta_k|^2 = \frac{k^4}{a^4 H^4} \Delta_{\mathcal{R}}^2. \quad (3.13)$$

Accordingly, $\sigma$ is determined by

<figure><img src="image_6.png" /><figcaption>Fig. 2. $\beta(t_−)$ varying from the model parameters $H_−$ and $c_s$. The value of $\beta$ is illustrated by the color, decreasing from red to purple. The isoline of $\beta(t_−) = 0.1$ is shown. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

$$\sigma^2 = \int \frac{dk}{k} \Delta_{\delta}^2 W \left( \frac{kR_J}{a} \right) = \int_0^{k_J} \frac{dk}{k} \Delta_{\delta}^2, \quad (3.14)$$

where $W(\frac{kR_J}{a})$ is the window function. For simplicity, we take

$$W \left( \frac{kR_J}{a} \right) = \begin{cases} 1 & k \le k_J \\ 0 & k > k_J \end{cases},$$

and then obtain

$$\sigma = \frac{3}{4c_s^2} \sqrt{\frac{a^3 H^2}{48\pi^2 c_s M_P^2 a^3}}. \quad (3.15)$$

Therefore, the PBH mass fraction by (3.3) is

$$\beta(t) \approx \text{erfc} \left[ \left( 48c_s^{3/2} + 9.355c_s^{1/2} \right) \left( \frac{M_P}{|H_-|} \right) \left( \frac{a(t)}{a_-} \right)^{3/2} \right] \\ - \text{erfc} \left[ \frac{10.9427}{c_s^{1/2}} \left( \frac{M_P}{|H_-|} \right) \left( \frac{a(t)}{a_-} \right)^{3/2} \right], \quad (3.16)$$

for $c_s < 0.39$; and $\beta(t) = 0$ for $c_s \ge 0.39$. It can be seen that the value of $\beta$ increases as $a(t)$ becomes smaller along with the contracting and reaches its maximum value at the end of the matter contracting phase $t_−$.

From the above expression, one can read that the maximal value of $\beta$ is unity, which corresponds to the case that the universe is dominated by black holes. Accordingly, we expect that $\beta \ll 1$ so that the universe is dominated by the background dust-like fluid. By numerically solving Eq. (3.16), we therefore derive the parameter space of $|H_-|$ and $c_s$ for each fixed value of $\beta$ as shown in Fig. 2. From this figure, one can read that the non-vanishing mass fraction $\beta(t_−)$, when $c_s < 0.39$, is small in both the very low energy regime $|H_-| \ll M_P$ and very high one $|H_-| \gg M_P$. In the high energy regime, the density fluctuation is generally large $\sigma \gg 1$ according to Eq. (3.15), but $\beta(t_−)$ is small due to the constraint that PBHs ought to be within the cosmic apparent horizon. The constraint of $\beta$ from today's observation could be loose. For instance, if one assumes that PBHs constitute the totality of dark matter in the universe, then the value of $\beta$ could be of order

<footnote n="2">² In this case, the PBH is formed by an over-dense core with radius $R_J$, and a shell of the background fluid within the radius from $R_J$ to $|H^{-1}|$.</footnote>
O(0.1) at most which might have occurred at the moment of matter-radiation equality (see for example [68] for related analyses). Therefore, if one considers the case of $\beta \lesssim 0.1$ and $c_s = 10^{-4}$ as a specific example, then Eq. (3.16) yields $|H_-| \lesssim 10^{-1} M_p$ or $|H_-| \gtrsim 10^4 M_p$. Note that, $|H_-| \lesssim 10^{-1} M_p$ can be easily satisfied, but $|H_-| \gtrsim 10^4 M_p$ is disfavored from model construction.

## 4. Evolution of PBH in bouncing phase

It is important to point out that the subsequent evolutions of PBHs after their formation are not considered in the above section. Associated with the evolutions, including the accretion and the evaporation due to Hawking radiation, one can obtain constraints upon bouncing cosmology. This is the main subject in this section. For simplicity, we only consider the PBH evolutions during the bouncing phase.

### 4.1. The growth of PBH mass

At first, we calculate the growth of PBH by accreting the mass around in the contracting background. For simplicity, we still assume that the space outside the BH horizon is unperturbed Friedmann universe and the matter which flows into the horizon is due to the cosmic contracting. For a BH with mass $M$ and radius $R$, the mass increased equals to that flowing into the horizon, with the speed $|H|R$,

$$ \dot{M} = 4\pi R^2 \bar{\rho} |H|R. \quad (4.1) $$

We shall investigate the evolution in the contracting era of the bouncing phase $t_- < t < 0$. During this stage, the Hubble length diverges quickly (see Fig. 1), and all PBHs can be treated as Schwarzschild BHs and $R = 2GM$, according to the arguments in the last section. Therefore, Eq. (4.1) reduces to

$$ \dot{R} = -3R^3 H^3 = -3R^3 \gamma^3 t^3, \quad (4.2) $$

and the solution is

$$ R(t) = \sqrt{\frac{1}{\frac{1}{R^2(t_-)} + \frac{3}{2}\gamma^3(t^4 - H_-^4/\gamma^4)}}. \quad (4.3) $$

It is obvious that, the BH radius increases with $t$ approaching to the bouncing point. Eq. (4.3) also gives a constraint as follows,

$$ R^2(t_-) \le \frac{2\gamma}{3H_-^4}. \quad (4.4) $$

For the models violating the constraints above, the mass of PBH will grow to infinity before the bouncing point. According to (3.1), one can read that $R(t_-) \ge \sqrt{\frac{2}{3}\pi c_s |H_-^{-1}|}$, and hence, the above constraint yields explicitly that<fnref n="3" />

$$ \gamma \ge c_s^2 \pi^2 H_-^2, \quad (4.5) $$

for $c_s \le 0.39$. Eq. (4.5) implies that the bouncing phase is expected to proceed rapidly so that the mass of PBHs would not become divergent and the universe is safe from collapsing into the black hole completely.

### 4.2. Backreaction and theoretical constraints

The bouncing phase is driven by the background fluid which violates NEC in the frame of general relativity. From the Friedmann equations not considering PBHs, the density and pressure of the background fluid at the bouncing point are $\rho_{bg} = 0$ and $p_{bg} = -2M_p^2\gamma$ respectively. To include the contribution of backreaction from PBHs, the Friedmann equation at the bouncing point can be written as

$$ \frac{\ddot{a}}{a} = -\frac{1}{6M_p^2}(\rho_{BR} - 6M_p^2\gamma), \quad (4.6) $$

where $\rho_{BR} = \rho_{PBH} + \rho_\gamma + 3p_\gamma$ is the effective density of the backreaction, $\rho_{PBH}$ is the density of PBH, $\rho_\gamma$ and $p_\gamma$ are the density and pressure of Hawking radiation respectively. For simplicity, we assume all Hawking radiation particles are ultra-relativistic, which satisfy the relation $p_\gamma = \rho_\gamma/3$. If the back reactions neutralize the negative pressure of background $\rho_{BR} - 6M_p^2\gamma > 0$, the universe would not expand any more, as mentioned in Section 2. In the following, we will constrain the model through PBHs.

Recall that the energy density of PBHs is

$$ \rho_{PBH} = \frac{\langle M \rangle}{L^3} = \frac{4\pi M_p^2 \langle R \rangle}{L^3}, \quad (4.7) $$

which depends on the mean PBH mass $\langle M \rangle$ or the mean PBH radius $\langle R \rangle$, and the mean PBH separation $L$. At the beginning of the bouncing phase $t_-$, the PBH density is known

$$ \rho_{PBH}(t_-) = \beta(\rho_{PBH} + \rho_{bg}) = 3M_p^2 H_-^2 \beta(t_-), \quad (4.8) $$

according to Eq. (3.2), where the evolution of the earlier formed PBHs is not considered. If $\langle M \rangle$ at the moment $t_-$ is known, one can obtain the mean PBH separation $L(t_-)$ from Eqs. (4.7) and (4.8). To obtain $\langle M(t_-) \rangle$, one should first calculated the mass function of PBHs. Here we assume that all PBHs have the same mass, which is a simplification generally used [49], and all PBHs are of the Jeans scale $\langle R(t_-) \rangle \approx \sqrt{\frac{2}{3}\pi c_s |H_-^{-1}|}$. Therefore, the mean separation at $t_-$ is

$$ L(t_-) = \frac{R(t_-)}{(c_s^2 \beta)^{\frac{1}{3}}}. \quad (4.9) $$

At the present investigation, we have neglected the newly formed PBHs during the bouncing phase. Therefore, the mean separation follows $L(t) \propto a(t)$, and accordingly, $L(0) \simeq L(t_-)e^{-\frac{H_-^2}{2\gamma}}$. Note that the PBH evaporation due to Hawking radiation has been studied clearly in [43]. Following this, one can estimate the lifetime of PBH to be

$$ \tau \approx 13.7 \text{Gyr} \left( \frac{M}{5 \times 10^{14} \text{g}} \right)^3 \approx 5 \times 10^3 \frac{M^3}{M_p^4}. \quad (4.10) $$

We have adopted the approximation that for massive PBHs $\tau > |t_-|$, or $\gamma > H_-^4/(10^8 c_s^3 M_p^2)$, the Hawking radiation is subdominant. A closer analysis combing both the Hawking radiation and the growth of PBHs will be addressed in our follow-up study. As a result, the evolution of the PBHs follows Eq. (4.3), and their size approximately takes $R(0) \simeq R(t_-) \left( \frac{\gamma}{\gamma - c_s^2 \pi^2 H_-^2} \right)^{1/2}$ at the bouncing point. Otherwise, PBH would have evaporated completely before the universe reaches the bouncing point with $R(0) = 0$. As a result, we can estimate the energy density of PBHs at the bouncing point as

<footnote n="3">³ We note that, in fact, if $c_s$ is not very small, the radius of BH suffers an increasing after the bouncing begins, but this effect has been ignored in the present consideration for simplicity.</footnote>
$$
\begin{align}
\rho_{\text{PBH}}(0) &= \frac{4\pi M_{\text{p}}^2 R(0)}{L(0)^3} \nonumber \\
&= \begin{cases}
\frac{6M_{\text{s}}^2 H^2}{\pi} \beta(t_-) \sqrt{\frac{\Upsilon}{\Upsilon - c_s^2 \pi^2 H_-^2}} e^{\frac{3H^2}{2\Upsilon}}, & \Upsilon > \frac{H^4}{10^8 c_s^3 M_{\text{p}}^2} \\
0, & \Upsilon < \frac{H^4}{10^8 c_s^3 M_{\text{p}}^2}
\end{cases} \tag{4.11}
\end{align}
$$

Similarly, we analyze the backreaction due to the Hawking radia-
tion. For $\Upsilon > \frac{H^4}{10^8 c_s^3 M_p^2}$, the Hawking radiation can be ignored and
$\rho_\Upsilon(0) = 0$. In the limit that PBHs finish their evaporation at $t_-$,
one has $\rho_\Upsilon(t_-) = \rho_{\text{PBH}}(t_-)$ and then $\rho_\Upsilon \propto a^{-4}$. In the limit that
the effective evaporation happens at the bouncing point, $\rho_\Upsilon(0) =
\rho_{\text{PBH}}(0)$, and $\rho_{\text{PBH}} \propto a^{-3}$ before its evaporation. Therefore, one has

$$
\rho_{\gamma}(0) = \begin{cases} 0, & \Upsilon > \frac{H_{-}^{4}}{10^{8}c_{s}^{3}M_{p}^{2}}, \\ \rho_{\text{PBH}}(t_{-}) \left(e^{\frac{H_{-}^{2}}{2\Upsilon}}\right)^{n}, & \Upsilon < \frac{H_{-}^{4}}{10^{8}c_{s}^{3}M_{p}^{2}} \end{cases} \quad (4.12)
$$

where $\rho_{\text{PBH}}(t_{-})$ is given by Eq. (4.8) and $n$ is a parameter that
describes the moment of the complete evaporation of PBH. For
instance, $n=4$ corresponds to the limit that the evaporation hap-
pens at $t_{-}$; and, $n=3$ corresponds to this limit at the bouncing
point.

Inserting Eqs. (4.11) and (4.12) into the Friedmann equation (4.6), in order not to neutralize the negative pressure $\bar{p} = -2M_p^2\Upsilon$, the following constraint should be satisfied

$$
\frac{\rho_{\text{BR}}}{|\bar{p}|} = \begin{cases} \frac{\beta(t_{-})H^2}{\pi \Upsilon} \sqrt{\frac{\Upsilon}{\Upsilon - c_s^2 \pi^2 H_-^2}} e^{\frac{3H^2}{2\Upsilon}} < 1, & \Upsilon > \frac{H_-^4}{10^8 c_s^3 M_p^2} \\ \frac{\beta(t_{-})H^2}{\Upsilon} \left(e^{\frac{H_-^2}{2\Upsilon}}\right)^n < 1, & \Upsilon < \frac{H_-^4}{10^8 c_s^3 M_p^2} \end{cases} . \quad (4.13)
$$

The resulting constraints are numerically shown in Fig. 3 for the
model with $c_s = 10^{-3}$. It is seen that the lines of $\Upsilon = c_s^2 \pi^2 H_-$
and $\Upsilon = \frac{H^4}{10^8 c_s^3 M_p^2}$, which represent the limit case of (4.5) and the
boundary between the PBH growth and evaporation respectively,
play important roles in the constraints. The two lines shape the
forbidden parameter space in which PBH grows into infinity. Fur-
thermore, two lines intersect at the characteristic energy scale
$H^2 \simeq 10^9 c_s^5 M_p^2$. At the low energy $H^2 \ll 10^9 c_s^5 M_p^2$, $\Upsilon = c_s^2 \pi^2 H^2$ is
the asymptote for all isolines of $\rho_{\text{BR}}/|\bar{p}|$. Therefore, the constraint
$\rho_{\text{BR}}/|\bar{p}| < 1$ reduces to Eq. (4.5) at the low energy.

4.3. *Estimates of observational constraints*

After the bouncing phase, PBHs will evolve in the regular ex-
panding universe, and hence can be constrained by today's obser-
vations. In this subsection, we briefly estimate the constraints upon
the bouncing models though the observational limits of PBHs.

For massive PBHs with an initial mass $M > 10^{15}$ g, they could
survive until today and their density scales as $a^{-3}$ approximately
in the expanding universe [43]. As a result, the current PBH ($M >$
$10^{15}$ g) density parameter can be expressed as

$$
\Omega_{\text{PBH}} \approx \frac{\rho_{\text{PBH}}(0)}{3 M_{\text{p}}^2 H_0^2} \left(\frac{a_b}{a_0}\right)^3,
$$

where $a_0$ and $H_0$ denote the scale factor and Hubble parameter of
today, $a_b$ is the scale of bounce, and $\rho_{\text{PBH}}(0)$ is the PBH density at
the bouncing point. Since the PBHs may contribute to part of dark
matter, one immediately has the rough constraint $\Omega_{\text{PBH}} < \Omega_{\text{DM}} \simeq$
0.25 and thus can get a bound

<figure><img src="image_7.png" /><figcaption>**Fig. 3.** The back reaction $\rho_{\text{BR}}/|\bar{p}|$ varying from model parameters $H_{-}$ and $\Upsilon$, with $c_s = 10^{-3}$ and $n=4$ taken. The value of $\rho_{\text{BR}}/|\bar{p}|$ is illustrated by color, increasing from purple to red. The white thick curve is the isoline $\rho_{\text{BR}}/|\bar{p}|=1$, and only the parameter space over this curve is allowed. The light blue curve is the boundary between PBH growth and evaporation $\Upsilon = \frac{H^4}{10^8 c_s^3 M_p^2}$, the area over the boundary is the PBH growing region, and that under the boundary is the evaporation region. The gray region is the forbidden parameter space due to Eq. (4.5), and the upper boundary of this region is $\Upsilon = c_s^2 \pi^2 H^2$ illustrated by green curve. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

$$
\rho_{\text{PBH}}(0) < 3 M_{\text{p}}^2 H_0^2 \Omega_{\text{DM}} \left(\frac{a_b}{a_0}\right)^{-3} . \qquad (4.14)
$$

For the light PBHs with $M < 10^{15}$ g, they have been completely
evaporated today. Therefore, the density parameter of the Hawking
radiation $\Omega_{\text{HR}}$ is approximately estimated as

$$
\Omega_{\text{HR}} \approx \frac{\rho_{\text{PBH}}(0)}{3 M_{\text{p}}^2 H_0^2} \left(\frac{a_r}{a_0}\right)^4 \left(\frac{a_b}{a_r}\right)^3 > \frac{\rho_{\text{PBH}}(0)}{3 M_{\text{p}}^2 H_0^2} \left(\frac{a_b}{a_0}\right)^4, \quad (4.15)
$$

where $a_r$ is the scale factor when PBHs complete the evaporation,
and $a_b < a_r < a_0$. $\Omega_{\text{HR}}$ is not allowed to be larger than the total
density of today's radiation $\Omega_r \sim 10^{-4}$, measured by the CMB ex-
periments. As a consequence, one gets

$$
\rho_{\text{PBH}}(0) < 3 M_{\text{p}}^2 H_0^2 \Omega_{\text{r}} \left(\frac{a_b}{a_0}\right)^{-4} . \qquad (4.16)
$$

In Fig. 4, we provide the constraints on $\rho_{\text{PBH}}(0)$ along with $a_b/a_0$
from the bounds (4.14) and (4.16), respectively. The PBH mass $M$
and density $\rho_{\text{PBH}}(0)$ are functions of $c_s$, $H_-$ and $\Upsilon$, which are pro-
vided in the previous subsection. As a result, we find that the non-
singular bounce models can be constrained by observations due to
(4.14) and (4.16). We note that these constraints are quite loose.
More stringent constraints can be obtained by taking the observa-
tional bounds of $\Omega_{\text{PBH}}$ and $\Omega_{\text{HR}}$ from the CMB, gamma-ray burst
and BBN etc. [43].

5. Conclusion and outlook

In this paper, we have investigated the formation and evolu-
tion of PBHs in the matter bounce cosmology. Firstly, we described
<figure><img src="image_8.png" /><figcaption>**Fig. 4.** The constraints on $\rho_{PBH}(0)$ along with $a_p/a_0$. The blue shadowed region is given by (4.14) for PBHs with $M > 10^{15}$ g, and the red shadowed region is given by (4.16) for $M < 10^{15}$ g. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

the general matter bounce models by some parameters like $c_s$, $H_-$ and $\Upsilon$. The comoving curvature perturbation $R_k$ is also calculated during the matter contracting phase, which seeds the PBH formation. Then we had a discussion about the condition of PBH formation in the contracting background, which is different from that in the expanding universe. By taking a simple collapsing model, the threshold of the density fluctuation for forming a PBH is derived. Furthermore, in the comoving gauge, the density fluctuation and its threshold are all related to the curvature perturbation $R_k$. Therefore, we can calculate the mass fraction of PBHs $\beta$ in the contracting phase from the Press-Schechter theory, and constrain the bouncing models from PBHs. PBH formation depends on the model parameters $c_s$ and $H_-$.

For instance, PBHs can form in the contracting phase only if $c_s \le 0.39$, since the BHs are not allowed to be larger than the cosmic apparent horizon. When $c_s \le 0.39$, $\beta$ is small and the model is safe for $|H_-| \ll M_p$ or $|H_-| \gg M_p$, but only the low energy regime is favored from the model construction.

The subsequent evolution of PBHs in the bouncing phase is also investigated, with the PBH accretion and Hawking radiation considered. The growth behavior of PBH yields a constraint to the model $\Upsilon \ge c_s^2\pi^2H^2$, in case that the mass of PBHs grow to infinity before the bouncing point. Moreover, the back reaction of PBH and its Hawking radiation is calculated to constrain models, in order not to neutralize the negative pressure $p_{bg} = -2M_p^2\Upsilon$, since in such case the universe cannot expand again. The constraint reduces to $\Upsilon \ge c_s^2\pi^2H^2$ at the low energy scale $H^2 \ll 10^9 c_s^5 M_p^2$. Afterwards, a rough constraint of bouncing model though the PBH observations is given, and the constraint is stringent only when $\Upsilon$ is slightly larger than $c_s^2\pi^2H^2$. A more precise analysis taken into account both the PBH growth and the Hawking evaporation as well as the detailed constraints from cosmological observations will be addressed in our follow-up works.

We note that while our paper was being prepared, an independent work was being carried out by another group [59], which explores similar features of BH formation in bouncing cosmology.

## Acknowledgements

We are grateful to R. Brandenberger, J. Martin, P. Peter, T. Qiu, J. Quintin, D.-G. Wang, Z. Wang, X.-Y. Yang and Y. Zhang for valuable comments. We especially thank the anonymous Referee for very detailed and constructive suggestions. YFC and JWC are supported in part by the Chinese National Youth Thousand Talents Program (Grant No. KJ2030220006), by the USTC start-up funding (Grant No. KY2030000049) and by the National Natural Science

Foundation of China (Grant Nos. 11421303, 11653002). HLX is supported in part by the Fund for Fostering Talents in Basic Science of the National Natural Science Foundation of China (Grant No. J1310021) and in part by the Yan Jici Talent Students Program. HLX also thanks J. Martin and P. Peter for their warmest hospitality during visiting the Institut d'Astrophysique de Paris when this work was finalized. Part of numerical computations are operated on the computer cluster LINDA in the particle cosmology group at USTC.

## References

[1] R.H. Brandenberger, arXiv:1206.4196 [astro-ph.CO].

[2] D. Wands, Phys. Rev. D 60 (1999) 023507, arXiv:gr-qc/9809062.

[3] F. Finelli, R. Brandenberger, Phys. Rev. D 65 (2002) 103522, arXiv:hep-th/0112249.

[4] R.H. Brandenberger, Int. J. Mod. Phys. Conf. Ser. 01 (2011) 67, arXiv:0902.4731 [hep-th].

[5] M. Novello, S.E.P. Bergliaffa, Phys. Rep. 463 (2008) 127, arXiv:0802.1634 [astro-ph].

[6] J.L. Lehners, Phys. Rep. 465 (2008) 223, arXiv:0806.1245 [astro-ph].

[7] Y.F. Cai, Sci. China, Phys. Mech. Astron. 57 (2014) 1414, arXiv:1405.1369 [hep-th].

[8] D. Battefeld, P. Peter, Phys. Rep. 571 (2015) 1, arXiv:1406.2790 [astro-ph.CO].

[9] R. Brandenberger, P. Peter, arXiv:1603.05834 [hep-th].

[10] A.H. Guth, Phys. Rev. D 23 (1981) 347.

[11] A.A. Starobinsky, Phys. Lett. B 91 (1980) 99.

[12] K. Sato, Mon. Not. R. Astron. Soc. 195 (1981) 467.

[13] L.Z. Fang, Phys. Lett. B 95 (1980) 154.

[14] F.C. Adams, K. Freese, A.H. Guth, Phys. Rev. D 43 (1991) 965.

[15] E.J. Copeland, A.R. Liddle, D.H. Lyth, E.D. Stewart, D. Wands, Phys. Rev. D 49 (1994) 6410, arXiv:astro-ph/9401011.

[16] A. Borde, A. Vilenkin, Phys. Rev. Lett. 72 (1994) 3305, arXiv:gr-qc/9312022.

[17] J. Martin, R.H. Brandenberger, Phys. Rev. D 63 (2001) 123501, arXiv:hep-th/0005209.

[18] Y.F. Cai, T. Qiu, Y.S. Piao, M. Li, X. Zhang, JHEP 0710 (2007) 071, arXiv:0704.1090 [gr-qc].

[19] Y.F. Cai, T. Qiu, R. Brandenberger, Y.S. Piao, X. Zhang, JCAP 0803 (2008) 013, arXiv:0711.2187 [hep-th]; Y.F. Cai, X. Zhang, JCAP 0906 (2009) 003, arXiv:0808.2551 [astro-ph].

[20] Y.F. Cai, T.T. Qiu, R. Brandenberger, X.M. Zhang, Phys. Rev. D 80 (2009) 023511, arXiv:0810.4677 [hep-th].

[21] R. Brandenberger, Phys. Rev. D 80 (2009) 043516, arXiv:0904.2835 [hep-th].

[22] Y.F. Cai, E.N. Saridakis, JCAP 0910 (2009) 020, arXiv:0906.1789 [hep-th].

[23] X. Gao, Y. Wang, W. Xue, R. Brandenberger, JCAP 1002 (2010) 020, arXiv:0911.3196 [hep-th].

[24] Y.F. Cai, S.H. Chen, J.B. Dent, S. Dutta, E.N. Saridakis, Class. Quantum Gravity 28 (2011) 215011, arXiv:1104.4349 [astro-ph.CO].

[25] J. de Haro, J. Amoros, Phys. Rev. Lett. 110 (2013) 071104, arXiv:1211.5336 [gr-qc].

[26] Y.F. Cai, S. Capozziello, M. De Laurentis, E.N. Saridakis, Rep. Prog. Phys. 79 (2016) 106901, arXiv:1511.07586 [gr-qc].

[27] C. Lin, R.H. Brandenberger, L. Perreault Levasseur, JCAP 1104 (2011) 019, arXiv:1007.2654 [hep-th].

[28] T. Qiu, J. Evslin, Y.F. Cai, M. Li, X. Zhang, JCAP 1110 (2011) 036, arXiv:1108.0593 [hep-th].

[29] D.A. Easson, I. Sawicki, A. Vikman, JCAP 1111 (2011) 021, arXiv:1109.1047 [hep-th].

[30] Y.F. Cai, D.A. Easson, R. Brandenberger, JCAP 1208 (2012) 020, arXiv:1206.2382 [hep-th].

[31] Y.F. Cai, E. McDonough, F. Duplessis, R.H. Brandenberger, JCAP 1310 (2013) 024, arXiv:1305.5259 [hep-th].

[32] Y.F. Cai, E. Wilson-Ewing, JCAP 1403 (2014) 026, arXiv:1402.3009 [gr-qc].

[33] S. Alexander, C. Bambi, A. Marciano, L. Modesto, Phys. Rev. D 90 (2014) 123510, arXiv:1402.5880 [gr-qc].

[34] S. Alexander, Y.F. Cai, A. Marciano, Phys. Lett. B 745 (2015) 97, arXiv:1406.1456 [gr-qc].

[35] R.H. Brandenberger, AIP Conf. Proc. 1268 (2010) 3, arXiv:1003.1745 [hep-th].

[36] J. Quintin, Z. Sherkatghanan, Y.F. Cai, R.H. Brandenberger, Phys. Rev. D 92 (2015) 063532, arXiv:1508.04141 [hep-th].

[37] L. Battarra, M. Koehn, J.L. Lehners, B.A. Ovrut, JCAP 1407 (2014) 007, arXiv:1404.5067 [hep-th].

[38] Y.F. Cai, E. Wilson-Ewing, JCAP 1503 (2015) 006, arXiv:1412.2914 [gr-qc].

[39] Y.F. Cai, F. Duplessis, D.A. Easson, D.G. Wang, Phys. Rev. D 93 (2016) 043546, arXiv:1512.08979 [astro-ph.CO].

[40] Y.F. Cai, A. Marciano, D.G. Wang, E. Wilson-Ewing, arXiv:1610.00938 [astro-ph.CO].

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130
[41] Y.F. Cai, J. Quintin, E.N. Saridakis, E. Wilson-Ewing, JCAP 1407 (2014) 033, arXiv:1404.4364 [astro-ph.CO].

[42] B.J. Carr, Astrophys. J. 201 (1975) 1.

[43] B.J. Carr, K. Kohri, Y. Sendouda, J. Yokoyama, Phys. Rev. D 81 (2010) 104019, arXiv:0912.5297 [astro-ph.CO].

[44] B.J. Carr, S.W. Hawking, Mon. Not. R. Astron. Soc. 168 (1974) 399.

[45] D.H. Lyth, K.A. Malik, M. Sasaki, I. Zaballa, JCAP 0601 (2006) 011, arXiv:astro-ph/0510647.

[46] I. Zaballa, A.M. Green, K.A. Malik, M. Sasaki, JCAP 0703 (2007) 010, arXiv:astro-ph/0612379.

[47] M. Kopp, S. Hofmann, J. Weller, Phys. Rev. D 83 (2011) 124025, arXiv:1012.4369 [astro-ph.CO].

[48] T. Harada, C.M. Yoo, K. Kohri, Phys. Rev. D 88 (2013) 084051, Erratum: Phys. Rev. D 89 (2014) 029903, arXiv:1309.4201 [astro-ph.CO].

[49] A.S. Josan, A.M. Green, K.A. Malik, Phys. Rev. D 79 (2009) 103520, arXiv:0903.3184 [astro-ph.CO].

[50] B.J. Carr, A.A. Coley, Int. J. Mod. Phys. D 20 (2011) 2733, arXiv:1104.3796 [astro-ph.CO].

[51] B.J. Carr, Springer Proc. Phys. 170 (2016) 23, arXiv:1402.1437 [gr-qc].

[52] J. Martin, P. Peter, Phys. Rev. D 68 (2003) 103517, arXiv:hep-th/0307077.

[53] D.M. Solomons, P. Dunsby, G. Ellis, Class. Quantum Gravity 23 (2006) 6585, arXiv:gr-qc/0103087.

[54] J. Quintin, Y.F. Cai, R.H. Brandenberger, Phys. Rev. D 90 (2014) 063507, arXiv: 1406.6049 [gr-qc]. 66

[55] J. de Haro, Y.F. Cai, Gen. Relativ. Gravit. 47 (2015) 95, arXiv:1502.03230 [gr-qc]. 67

[56] P. Qian, Y.F. Cai, D.A. Easson, Z.K. Guo, arXiv:1607.06578 [gr-qc]. 68

[57] Y.F. Cai, W. Xue, R. Brandenberger, X.m. Zhang, JCAP 0906 (2009) 037, arXiv: 0903.4938 [hep-th]. 69

[58] C. Gao, X. Chen, Y.G. Shen, V. Faraoni, Phys. Rev. D 84 (2011) 104047, arXiv: 1110.6708 [gr-qc]. 71

[59] J. Quintin, R.H. Brandenberger, JCAP 1611 (11) (2016) 029, arXiv:1609.02556 [astro-ph.CO]. 72

[60] W.H. Press, P. Schechter, Astrophys. J. 187 (1974) 425. 74

[61] J.E. Gunn, J.R. Gott III, Astrophys. J. 176 (1972) 1. 75

[62] P.J.E. Peebles, Astrophys. J. 147 (1967) 859. 76

[63] R.C. Tolman, Proc. Natl. Acad. Sci. USA 20 (1934) 410. 77

[64] H. Bondi, Mon. Not. Astron. Soc. 107 (1947) 343. 78

[65] G. Lemaître, Ann. Soc. Sci. Brux. I A53 (1933) 51. 79

[66] D.H. Lyth, Phys. Rev. D 31 (1985) 1792. 79

[67] A.R. Liddle, D.H. Lyth, Cosmological Inflation and Large Scale Structure, Univ. Pr., Cambridge, UK, 2000, 400 pp. 80

[68] N. Bartolo, et al., arXiv:1610.06481 [astro-ph.CO]. 81

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130