Interpretating Pulsar Timing Array data of Gravitational Waves with Ekpyrosis-Bouncing Cosmology

The generic action is taken to be

$$S=\int d^{4}x\sqrt{-g}\left[M_{p}^{2}\frac{R}{2}+K(\phi,X)+\mathcal{L}_{\textnormal{EFT}}+\mathcal{L}_{\sigma}\right]~{},$$

(2.1)

where $R/2$ is the Einstein-Hilbert action. The action (2.1) is composed by three part, which we shall explain below.

The k-essence action $K(\phi,X)$ is responsible for the background dynamics. We assume that the universe starts with an Ekpyrotic contraction epoch at $\phi\ll-1$, and finally evolves into an radiation-dominated expansion epoch at $\phi\gg 1$. This suggests the following asymptotic behavior of $K$:

$$\lim_{\phi\to-\infty}K(\phi,X)=X+2V_{0}e^{\sqrt{\frac{2}{q}}\phi}~{};~{}\lim_{\phi\to\infty}K(\phi,X)=X^{2}~{}.$$

(2.2)

Moreover, during the bouncing epoch, NEC is violated. We will use a ghost-condensate type action [84], i.e. $-X+X^{2}$, to violate NEC. Therefore, the k-essence action shall take the form

$$K(\phi,X)=[1-g(\phi)]M_{p}^{2}X+\beta_{2}X^{2}-M_{p}^{4}V(\phi)~{},$$

(2.3)

with $g$ being the smooth top-hat function

	$\displaystyle g(\phi)$	$\displaystyle=\beta_{1}\left[\frac{1+\tanh\lambda_{1}(\phi-\phi_{-})}{2}\right]\left[\frac{1+\tanh\lambda_{2}(\phi_{+}-\phi)}{2}\right]$
		$\displaystyle+\frac{1+\tanh\lambda_{4}(\phi-\phi_{+})}{2}~{},~{}\phi_{-}<\phi_{+}~{},$		(2.4)

and a potential

$$V(\phi)=-2V_{0}e^{\sqrt{\frac{2}{q}}\phi}\frac{1-\tanh\lambda_{3}(\phi-\phi_{-})}{2}~{},$$

(2.5)

with the restriction $\lambda_{3}>1/\sqrt{2q}$. One may check that in the ansatz, (2.3) along with (2.1) and (2.5), the asymptotic behavior (2.2) is satisfied (the potential $V(\phi)$ is exponentially suppressed so that the action is dominated by the $X^{2}$ term when $\phi>\phi_{+}$). Besides, when $\phi_{-}<\phi<\phi_{+}$, the action (2.3) becomes

$$K(\phi,X)\simeq(1-\beta_{1})M_{p}^{2}X+\beta_{2}X^{2}-M_{p}^{4}V(\phi)~{}.$$

(2.6)

which permits an NECV as long as $\beta_{1}>1$.

There will be ghost or gradient instability problems in the generic bouncing models within the framework of Horndeski theory [85, 86, 87]. The simplest way to evade the instabilities is to introduce an EFT operator in the context of non-singular cosmology [88, 89, 90, 91] (alternative approaches can be found in [92, 93, 94, 95, 79, 96, 97, 98, 99, 100, 101, 102, 103, 80]). Since the EFT operator changes the scalar sound speed $c_{s}^{2}$ only, and has no impact on the background dynamics and tensor perturbations, we may simply use it to resolve the instability problem for bouncing cosmology, and don’t need to worry about its contribution to both the background evolution and linear perturbations.

However, the Ekpyrotic contraction period will generate a blue-tilted scalar power spectrum [104, 105, 106, 107, 108, 109]. Therefore, it is useful to introduce a curvaton field $\sigma$ coupled to the scalar field $\phi$, which can lead to the desired scale-invariant scalar power spectrum, while leaving the bouncing background unchanged [83, 110, 111, 112, 113, 114, 70]. ³³3Curvaton has been discussed a lot in recent years with many new models proposed, see [115, 116, 117, 118, 119, 120].As shown explicitly in [83], we can adopt the follwing action for curvaton:

$$\mathcal{L}_{\sigma}=-\frac{M_{p}^{2}}{2}f^{2}e^{\sqrt{2/q}\phi}(\partial\sigma)^{2}~{},$$

(2.7)

acquiring a scale-invariant scalar power spectrum which is constant in time. Here, $f$ is a dimensionless parameter. Unfortunately, in the simple set (2.7), the curvaton field will be dominant in the expansion stage where $\phi\gg 1$. In view of that, a realistic choice of curvaton action can be

$$\mathcal{L}_{\sigma}=-M_{p}^{2}f^{2}\frac{\sqrt{\frac{2}{q}}+\phi e^{\sqrt{\frac{q}{2}}\phi}}{e^{-\sqrt{\frac{q}{2}}\phi}+\phi^{2}e^{\sqrt{\frac{q}{2}}\phi}}(\partial\sigma)^{2}~{},$$

(2.8)

which exactly has the asymptotic form (2.7) when $\phi\ll 1$. At the expansion stage, the action (2.8) approximates to $\mathcal{L}_{\sigma}\propto\phi^{-1}(\partial\sigma)^{2}$, which we shall explain in details in Sec. 5.3. Besides, one may check that the denominator and numerator in (2.8) is always positive.

In the following, we will start with the action (2.1), along with (2.3) and (2.8), to pursue its background and perturbation evolutions and to see how the primordial GWs can meet with the PTA data.

According to (2.1), (2.3) and (2.8), the Friedmann’s equations are

$$3H^{2}M_{p}^{2}=\rho=\frac{1}{2}(1-g)M_{p}^{2}\dot{\phi}^{2}+\frac{3}{4}\beta_{2}\dot{\phi}^{4}+M_{p}^{4}V(\phi)+\frac{M_{p}^{2}f^{2}}{2\cosh(\sqrt{2/q}\phi)}\dot{\sigma}^{2}~{},$$

(2.9)

$$-2\dot{H}M_{p}^{2}=\rho+p=(1-g)M_{p}^{2}\dot{\phi}^{2}+\beta_{2}\dot{\phi}^{4}+\frac{M_{p}^{2}f^{2}}{\cosh(\sqrt{2/q}\phi)}\dot{\sigma}^{2}~{}.$$

(2.10)

Since we now have two fields, we need one additional dynamical equation. We obtain the dynamics of $\sigma$ by varying the action (2.1) with respect to $\sigma$:

$$\ddot{\sigma}+\dot{\sigma}\frac{d}{dt}\left[\ln\left(a^{3}M_{p}^{2}f^{2}\cosh^{-1}(\sqrt{2/q}\phi)\right)\right]=0~{}.$$

(2.11)

Apart from the constant solution $\dot{\sigma}=0$, there’s another branch of solution

$$\dot{\sigma}a^{3}M_{p}^{2}f^{2}\cosh^{-1}(\sqrt{2/q}\phi)=\textnormal{const}~{}.$$

(2.12)

Before proceeding, we shall examine whether the curvaton field has negligible contributions to the background dynamics. In the contraction epoch, we have $\cosh^{-1}(\sqrt{2/q}\phi)\simeq e^{\sqrt{2/q}\phi}$, so the curvaton has no significant contribution to the background, as proved in [83]. In the radiation dominated epoch, we have $\dot{\phi}\propto t^{-1/2}$ such that $\phi$ grows as $t^{1/2}$. Thus $\dot{\sigma}\propto t^{-1}$ according to (2.8) and the energy density of curvaton field evolves as $\rho_{\sigma}\propto\phi^{-1}\dot{\sigma}^{2}=t^{-5/2}$. So $\rho_{\sigma}$ decays more rapidly than the background energy density and there is no backreaction problem in the radiation dominated epoch. In both cases, the curvaton field has no change to introduce a large back-reaction to the background.

Finally, it would be useful to write down the dynamical equation of $\phi$ in the numerical process. Since we’ve omitted the contributions of the curvaton field on the background level, the equation simplifies to

$$\left[(1-g)M_{p}^{2}+3\beta_{2}\dot{\phi}^{2}\right]\ddot{\phi}+3H\dot{\phi}\left((1-g)M_{p}^{2}+\beta_{2}\dot{\phi}^{2}\right)+M_{p}^{4}V^{\prime}(\phi)-\frac{M_{p}^{2}}{2}g^{\prime}(\phi)\dot{\phi}^{2}=0~{}.$$

(2.13)

In the far past $\phi\ll-1$ and $\dot{\phi}\ll M_{p}$ (so that the $X^{2}$ term is subdominate to the $X$ term), with the help of (2.2), the Friedmann equations (2.9), (2.10) become

$$3H^{2}M_{p}^{2}=\frac{1}{2}M_{p}^{2}\dot{\phi}^{2}-2V_{0}M_{p}^{4}e^{\sqrt{\frac{2}{q}}\phi}~{},~{}-2\dot{H}=\dot{\phi}^{2}~{},$$

(2.14)

which gives the Ekpyrotic attractor solution

$$\phi\simeq-\sqrt{\frac{q}{2}}\ln\left[\frac{2V_{0}M_{p}^{2}t^{2}}{q(1-3q)}\right]~{},~{}\dot{\phi}=\frac{\sqrt{2q}}{-t}~{};~{}H=\frac{q}{t}<0~{}.$$

(2.15)

The effective equation-of-state (EoS) parameter in the Ekpyrotic contraction epoch is solely determined by the parameter $q$: $w_{c}=-1+2/(3q)$. The spacetime geometry evolves accordingly:

$$a(\tau)=a_{-}\left(\frac{\tau_{e}-\tau}{\tau_{e}-\tau_{-}}\right)^{\frac{q}{1-q}}~{},~{}\mathcal{H}=\frac{-q}{(1-q)(\tau_{e}-\tau)}<0~{},$$

(2.16)

where $\tau\equiv\int dt/a$ is the conformal time and $\mathcal{H}=aH$ is the conformal Hubble parameter. $\tau_{-}$ is the end time of the contraction epoch and $a_{-}\equiv a(\tau_{-})$ is the corresponding scale factor. Finally, $\tau_{e}>\tau_{-}$ is an integration constant. It is convenient to have the relation between $t$ and $\tau$:

$$t=\int a(\tau)d\tau=-a_{-}(1-q)\left[\frac{(\tau_{e}-\tau)}{(\tau_{e}-\tau_{-})^{q}}\right]^{\frac{1}{1-q}}<0~{},$$

(2.17)

up to an integration constant.

The evolution of fluctuations in bouncing epoch is generically involved. Fortunately, the duration of the bouncing epoch is usually taken to be short. For example, in [121] it is pointed out that the sourced anisotropy in the bouncing epoch grows exponentially. Thus a short bounce is favored for the model to be free from anisotropic stress problem. In this case, the bouncing epoch can be parametrized in the following

$$\mathcal{H}\simeq\alpha^{2}(\tau-\tau_{B})~{},~{}a=a_{B}e^{\frac{1}{2}\alpha^{2}(\tau-\tau_{B})^{2}}~{},~{}\tau_{-}<\tau<\tau_{+}~{},$$

(2.18)

where $\tau_{B}$ is an integration constant and $\tau_{+}$ is the end of bouncing epoch (or equivalently, the beginning of expanding epoch). The validation of the linear approximation requires $|\alpha(\tau-\tau_{B})|<1$.

While the parameterization (2.18) can be understood as a taylor expansion around the bouncing point, the dynamics of scalar field $\phi$ in bouncing phase is more involved. In [72] the authors present an approximated formulae

$$\dot{\phi}\simeq M_{p}\sqrt{\frac{2(\beta_{1}-1)}{3\beta_{2}}}e^{-t^{2}/T^{2}}~{},$$

(2.19)

with $T$ a free parameter, which roughly equals to a quarter of the duration of bouncing phase. The validation of (2.19) is examined by the authors in [72] by numerics. We will use (2.19) to describe the dynamics of $\phi$ in bouncing phase in the rest of the paper.

After the bouncing epoch, the universe enters into an expansion phase. The Lagrangian of $\phi$ is $K(\phi,X)=X^{2}$, so that $\rho=3X^{2}$ and $p=X^{2}$ and the Equation-of-state parameter for $\phi$ becomes $w_{\phi}=1/3$. The universe will then be radiation-dominated. The corresponding background dynamics becomes

$$a=a_{+}\frac{\tau-\tilde{\tau}}{\tau_{+}-\tilde{\tau}}~{},~{}\mathcal{H}=\frac{1}{\tau-\tilde{\tau}}~{}.$$

(2.20)

In order to maintain the continuity of $a$ and ${\cal H}$, the background dynamics in each epoch can be furtherly matched by the junction condition at the transition surface $\tau=\tau_{-}$ and $\tau=\tau_{+}$. Firstly, we set the bouncing point $H=0$ to be the zero point of $\tau$, (i.e., $\tau_{B}=0$), so that $\tau_{-}<0$, $\tau_{+}>0$ and

$$\mathcal{H}=\alpha^{2}\tau~{},~{}a=a_{B}e^{\frac{1}{2}\alpha^{2}\tau^{2}}~{},~{}\tau_{-}<\tau<\tau_{+}~{},$$

(2.21)

the continuity at $\tau=\tau_{-}$ gives

$$\tau_{e}=\tau_{-}-\frac{q}{(1-q)\alpha^{2}\tau_{-}}~{},$$

(2.22)

and the continuity at $\tau=\tau_{+}$ gives $\tilde{\tau}=\tau_{+}-(\alpha^{2}\tau_{+})^{-1}$. For the sake of simplicity, let’s work in a symmetric bounce, i.e., $|\tau_{-}|=\tau_{+}$. Then we have $a_{-}=a_{+}=a_{B}e^{\alpha^{2}\tau_{+}^{2}/2}$. Implemented with those conditions, we summarize the evolution of scale factor

$$a(\tau)=\left\{\begin{array}[]{cl}&a_{B}e^{\frac{1}{2}\alpha^{2}\tau_{+}^{2}}\left[1+\frac{1-q}{q}\alpha^{2}|\tau_{-}|(\tau_{-}-\tau)\right]^{\frac{q}{1-q}}~{},~{}\tau<\tau_{-}~{},\\ \\ &a_{B}e^{\frac{1}{2}\alpha^{2}\tau^{2}}~{},~{}\tau_{-}<\tau<\tau_{+}~{},\\ \\ &a_{B}e^{\frac{1}{2}\alpha^{2}\tau_{+}^{2}}\left[1+\mathcal{H}_{+}(\tau-\tau_{+})\right]~{},~{}\tau>\tau_{+}~{}.\end{array}\right.$$

(2.23)

and of conformal Hubble parameter

$$\mathcal{H}(\tau)=\left\{\begin{array}[]{cl}&\alpha^{2}\tau_{-}\left[1+\frac{1-q}{q}\alpha^{2}|\tau_{-}|(\tau_{-}-\tau)\right]^{-1}~{},~{}\tau<\tau_{-}~{},\\ \\ &\alpha^{2}\tau~{},~{}\tau_{-}<\tau<\tau_{+}~{},\\ \\ &\alpha^{2}\tau_{+}\left[1+\alpha^{2}\tau_{+}(\tau-\tau_{+})\right]^{-1}~{},~{}\tau>\tau_{+}~{}.\end{array}\right.$$

(2.24)

The background value is thus fixed by four parameters: $\alpha$, $\tau_{+}$, $a_{B}$ and $q$.

We numerically evaluate the background dynamics using the following parameter set $\beta_{1}=1.1$, $\beta_{2}=1$, $\phi_{-}=-0.1$, $\phi_{+}=3.7$, $V_{0}=10^{-1}$, $q=0.005$, and present the evolution of Hubble parameters in Fig. 1. In the figure, the Hubble parameter goes from below zero to above zero, indicating that bounce indeed occurs. While before the bounce the EoS is very large as in the Ekpyrosis phase $w=-1+2/(3q)$, after the bounce the universe enters a radiation-dominated era where the EoS is nearly $1/3$. Note that all model parameters are normalized to be dimensionless.

Now we turn to the tensor perturbations generated in our model. The quadratic action of tensor perturbation is

$$S_{2,T}=\int d\tau d^{3}x\frac{a^{2}}{8}M_{p}^{2}\left[\gamma_{ij}^{\prime 2}-c_{T}^{2}(\partial\gamma_{ij})^{2}\right]~{},$$

(3.1)

where $\gamma_{ij}$ represents the tensor perturbation which is dimensionless. Here for simplicity, we don’t distinguish the polarization modes for the tensor sector. Since the matter sector is minimally coupled to gravity, the propogation speed of GWs is $c_{T}^{2}=1$. The dynamical equation is accordingly

$$\nu_{k}^{\prime\prime}+\left(k^{2}-\frac{a^{\prime\prime}}{a}\right)\nu_{k}=0~{},$$

(3.2)

where $\nu_{k}\equiv a\gamma_{k}/(2\sqrt{M_{p}})$ is the mode function of tensor perturbation.

We will need the initial condition to solve (3.2). The universe starts in an Ekpytoric contraction configuration described by (2.15) and (2.16), during which the dynamical equation is:

$$\nu_{k}^{\prime\prime}+\left[k^{2}-\frac{q(2q-1)}{(1-q)^{2}(\tau_{e}-\tau)^{2}}\right]\nu_{k}=0~{},$$

(3.3)

whose general solution is the Hankel function $H_{\nu}^{(1)}[k(\tau_{e}-\tau)]$ and $H_{\nu}^{(2)}[k(\tau_{e}-\tau)]$. Imposing the vacuum initial condition $\nu_{k}\sim e^{-ik\tau}/\sqrt{2k}$, the mode function evolves as

$$\nu_{k}(\tau)\simeq\frac{\sqrt{\pi(\tau_{e}-\tau)}}{2}H_{\nu}^{(1)}[k(\tau_{e}-\tau)]~{},~{}\nu\equiv\frac{1-3q}{2(1-q)}~{}.$$

(3.4)

Note that $\nu_{k}$ has a dimension of $[M]^{-1/2}$, as expected from the definition. We will set Equation (3.4) as the initial condition. The dimensionless primordial tensor power spectrum is defined as

$$P_{T}=2\cdot\frac{k^{3}}{2\pi^{2}M_{p}^{3}}|\gamma_{k}|^{2}=\frac{4k^{3}}{\pi^{2}M_{p}^{2}}\left|\frac{\nu_{k}}{a}\right|^{2}~{},$$

(3.5)

where the factor $2$ comes from the two tensorial polarizations.

Our task is to evaluate the tensor power spectrum using (3.2) and (3.4) for $k<\mathcal{H}_{+}$ modes. As we argued above, the modes of interest are super-horizon at the beginning of radiation dominated epoch, $\tau=\tau_{+}$. These modes are conserved before the horizon re-entry, so it suffices to evaluate the tensor spectrum at $\tau=\tau_{+}$. The dynamics of tensor perturbation during the contraction phase is dictated by (3.4). On super-horizon scale, the tensor perturbation and its derivative at $\tau=\tau_{-}$ is

$$\nu_{k}(\tau_{-})=-\frac{i\sqrt{\tau_{e}-\tau_{-}}}{2\sqrt{\pi}}\Gamma(\nu)\left[\frac{2}{k(\tau_{e}-\tau_{-})}\right]^{\nu}=-\frac{2^{\nu-1}i}{\sqrt{\pi}}\Gamma(\nu)k^{-\nu}\left[\frac{\alpha^{2}(1-q)\tau_{+}}{q}\right]^{\nu-\frac{1}{2}}~{},$$

(3.6)

$$\nu_{k}^{\prime}(\tau_{-})=\frac{i\Gamma(\nu)2^{\nu}}{(1-2\nu)\sqrt{\pi}}k^{-\nu}\left[\frac{\alpha^{2}(1-q)\tau_{+}}{q}\right]^{\nu+\frac{1}{2}}~{},$$

(3.7)

where we use the fact $|\tau_{-}|=-\tau_{-}=\tau_{+}$.

With (2.21), the dynamical equation for tensor perturbation in the bouncing phase approximates to

$$\nu_{k}^{\prime\prime}+\left(k^{2}-\alpha^{2}\right)\nu_{k}=0~{},$$

(3.8)

whose general solution is

$$\nu_{k}=b_{T,1}e^{\sqrt{\alpha^{2}-k^{2}}\tau}+b_{T,2}e^{-\sqrt{\alpha^{2}-k^{2}}\tau}~{},~{}\tau_{-}<\tau<\tau_{+}~{}.$$

(3.9)

We can interpret (3.8) in the following way. Modes with $k^{2}<\alpha^{2}$ will experience a tachyonic instability [70, 122], where the corresponding fluctuations are amplified. On the other hand, modes with $k^{2}>\alpha^{2}$ will simply exhibit an oscillatory behavior.

Finally, in the radiation dominated epoch, we have $a^{\prime\prime}/a=0$, thus all modes of our interest oscillate in this epoch, and their amplitude remain fixed. The power spectrum in radiation dominated epoch is invariant, and it would be sufficient for us to evaluate $P_{T}$ at $\tau=\tau_{+}$.

We depict the evolution of Hubble horizon, $|H|^{-1}$, in Fig. 2. It is easy to see that only modes with $k\leq|\mathcal{H}_{-}|=\mathcal{H}_{+}$ has chances to be super-horizon in the bouncing scenario. Those modes will cross the horizon firstly in the contraction phase, and re-enter then re-exit the horizon during the bouncing phase, which finally re-enter the horizon in the radiation dominated phase. Modes with $k>\mathcal{H}_{+}$ will be sub-horizon in the whole cosmic evolution.

In fact, the cross of the Hubble horizon of the primordial fluctuations is equivalent to its classicalization, which is crucial for them to be observable. Modes which never cross the Hubble horizon has a spurious divergence thus need specific treatment. As we are interested in the intepretation of recent PTA signals from PGWs originated from the contraction phase, we shall terminate our power spectrum at the scale $k_{\ast}=\mathcal{H}_{+}$, where $k_{\ast}$ refers to the upper limit of the PTA detectability. Moreover, for range of smaller $k$ (larger scales), the CMB constraint must be satisfied as well.

The remaining task is to evaluate the power spectra for the modes with $k\leq\mathcal{H}_{+}=\alpha^{2}\tau_{+}$. Here, $\alpha\tau_{+}$ should be at most of $\mathcal{O}(1)$ order, otherwise the scale factor shall change intensively in the bouncing epoch, in constrast with our short bounce assumption. Therefore, we have $\mathcal{H}_{+}\leq\alpha$, thus modes of our interest shall satisfy the condition $k<\alpha$. As a result, these modes will experience a tachyonic growth in the bouncing epoch as suggested by (3.9). In this case, the exponential growing mode will dominate over the exponential dacaying mode, so it suffices to evaluate $b_{T,1}$.

We calculate $b_{T,1}$ by applying the junction condition for tensor fluctuations, i.e., the mode function and its first derivative is continuous on the transition surface $\tau=\tau_{-}$:

$$-\frac{2^{\nu-1}i}{\sqrt{\pi}}\Gamma(\nu)k^{-\nu}\left[\frac{\alpha^{2}(1-q)\tau_{+}}{q}\right]^{\nu-\frac{1}{2}}=b_{T,1}e^{\sqrt{\alpha^{2}-k^{2}}\tau_{-}}+b_{T,2}e^{-\sqrt{\alpha^{2}-k^{2}}\tau_{-}}~{},$$

(3.10)

$$\frac{i\Gamma(\nu)2^{\nu}}{(1-2\nu)\sqrt{\pi}}k^{-\nu}\left[\frac{\alpha^{2}(1-q)\tau_{+}}{q}\right]^{\nu+\frac{1}{2}}=\sqrt{\alpha^{2}-k^{2}}\left(b_{T,1}e^{\sqrt{\alpha^{2}-k^{2}}\tau_{-}}-b_{T,2}e^{-\sqrt{\alpha^{2}-k^{2}}\tau_{-}}\right)~{}.$$

(3.11)

After some calculation one finds that it leads to

$$b_{T,1}=e^{-\sqrt{\alpha^{2}-k^{2}}\tau_{-}}\frac{2^{\nu-2}i}{\sqrt{\pi}}\Gamma(\nu)k^{-\nu}\left[\frac{q}{\alpha^{2}(1-q)\tau_{+}}\right]^{\frac{q}{1-q}}\left[\frac{\alpha^{2}(1-q)^{2}\tau_{+}}{\sqrt{\alpha^{2}-k^{2}}q^{2}}-1\right]~{},$$

(3.12)

where we apply $1-2\nu=2q/(1-q)$.

The tensor power spectrum at $\tau=\tau_{+}$ is thus

$$P_{T}=2^{-\frac{1+q}{1-q}}\frac{\Gamma^{2}(\nu)H_{+}^{2}}{\pi^{3}M_{p}^{2}a(\tau_{+})^{2}}k^{\frac{2}{1-q}}e^{2\sqrt{\alpha^{2}-k^{2}}\tau_{+}}\left[\frac{q}{\alpha^{2}(1-q)\tau_{+}}\right]^{\frac{2q}{1-q}}\left[\frac{\alpha^{2}(1-q)^{2}\tau_{+}}{\sqrt{\alpha^{2}-k^{2}}q^{2}}-1\right]^{2}~{}.$$

(3.13)

In terms of Hubble parameters:

$$P_{T}=2^{-\frac{1+q}{1-q}}\frac{\Gamma^{2}(\nu)}{\pi^{3}M_{p}^{2}}e^{\frac{2\sqrt{\alpha^{2}-k^{2}}\mathcal{H}_{+}}{\alpha^{2}}}\left(\frac{q}{1-q}\right)^{\frac{2q}{1-q}}\left(\frac{k/a_{B}}{H_{+}}\right)^{\frac{2}{1-q}}\left[\frac{(1-q)^{2}\mathcal{H}_{+}}{\sqrt{\alpha^{2}-k^{2}}q^{2}}-1\right]^{2}~{}.$$

(3.14)

We have $H_{+}=\alpha^{2}t_{+}$, since $\alpha t_{+}$ can be at most $\mathcal{O}(1)$, $H_{+}$ cannot exceed $\alpha$ too much. As a result, $\mathcal{H}_{+}\simeq a_{B}H_{+}$ must be much smaller than $\alpha$, since $a_{B}\ll a_{\rm today}$. Similarly, we have $k\leq k_{\ast}=\mathcal{H}_{+}\ll\alpha$. Thus, unless we take a vanishing $q$, the terms in bracket will approximates to unity. The expression simplifies to

$$P_{T}\simeq 2^{-\frac{1+q}{1-q}}\frac{\Gamma^{2}(\nu)}{\pi^{3}}\left(\frac{H_{+}}{M_{p}}\right)^{2}\left(\frac{q}{1-q}\right)^{\frac{2q}{1-q}}\left(\frac{k/a_{B}}{H_{+}}\right)^{\frac{2}{1-q}}~{}.$$

(3.15)

It’s easy to see that the tensor power spectrum is governed by three parameters: the Hubble parameter at the end of bouncing phase $H_{+}$, the model parameter $q$ which essentially determine the scale dependence of the primordial fluctuations, and the scale factor $a_{B}$ which is relevant to the post-bouncing phases.

From the above results we can see that, for $q\lesssim 1$, $2/(1-q)>0$, and we can get a blue-tilted tensor power spectrum. For current PTA data, it requires $n_{T}=1.8\pm 0.3$ [65], which suggests $q<0.0476$.

Moreover, in the minimal setup, the spectral energy density parameter $\Omega_{GW}(k)$ is related to the primordial tensor spectrum, defined as the energy density of the GWs per unit logarithmic frequency, by

$$\Omega_{GW}(k)\equiv\frac{1}{3H^{2}}\frac{d\rho_{\rm GW}}{d\ln k}\simeq 10^{-6}P_{T}(k)~{}.$$

(3.16)

The dynamics for the fluctuation of the Ekpyrotic field $\phi$ is [103]

$$\mu_{k}^{\prime\prime}+\left(k^{2}-\frac{z_{s}^{\prime\prime}}{z_{s}}\right)\mu_{k}=0~{},$$

(5.1)

where $\mu_{k}\equiv z_{s}\delta\phi$ is the mode function for $\phi$ and

$$\frac{z_{s}^{2}}{a^{2}}\equiv\frac{\dot{\phi}^{2}K_{X}+\dot{\phi}^{4}K_{XX}}{H^{2}M_{p}^{2}}=\frac{3\beta_{2}\dot{\phi}^{4}+2(1-\beta_{1})M_{p}^{2}\dot{\phi}^{2}}{4M_{p}^{2}H^{2}}~{}.$$

(5.2)

As we discussed in Sec. 3.3, only the modes that becomes super-horizon in the contraction phase are of observational interest to us. For these modes, the power spectrum from $\phi$ at $\tau=\tau_{-}$ is evaluated as [103]

$$P_{\phi}(k,\tau_{-})\simeq\frac{\Gamma^{2}(\nu)(1-2\nu)^{1-2\nu}}{3\pi^{3}2^{8-4\nu}}\left(\frac{k}{|\mathcal{H}_{-}|}\right)^{\frac{2}{1-q}}\frac{H_{-}^{2}}{M_{p}^{2}}~{},~{}\nu\equiv\frac{1-3q}{2(1-q)}~{},$$

(5.3)

which is strongly blue. Correspondingly, we define the tensor-to-scalar ratio of the scalar field $\phi$, $r_{\phi}\equiv P_{T}/P_{\phi}$, therefore, at $\tau=\tau_{-}$ it is just a number $r_{\phi}=96$.

The dynamics in the bouncing phase is more involved. At first glance, the parameter $z_{s}$ diverges at the bouncing point $H=0$ according to (5.2). However, if we use (2.19), the expression for $z_{s}$ becomes

$$\frac{z_{s}^{2}}{a^{2}}=\frac{M_{p}^{2}}{4H^{2}}\frac{(1-\beta_{1})^{2}}{3\beta_{2}}\left(e^{-\frac{t^{2}}{T^{2}}}-e^{-\frac{2t^{2}}{T^{2}}}\right)~{},$$

(5.4)

which converges at the limit $t\to 0$ since $H=\alpha^{2}t$ ⁵⁵5The regularity of $z_{s}^{2}$ explains the prefactor $\sqrt{2(\beta_{1}-1)/3\beta_{2}}$ in (2.19). . Moreover, at leading order we have $z_{s}/a=\rm{const}$ and thus $z_{s}^{\prime\prime}/z_{s}=a^{\prime\prime}/a$. In this case, the scalar fluctuation shares the same dynamics as the tensor fluctuation and as a result, the tensor-to-scalar ratio $r_{\phi}$ is invariant in the bouncing phase. We can thus immediately write down the power spectrum for $\phi$ at the end of bouncing phase:

$$P_{\phi}(k,\tau_{+})=P_{T}(k,\tau_{+})/r_{\phi}=\frac{1}{96}P_{T}(k,\tau_{+})~{}.$$

(5.5)

Before closing this section, we comment on the difference between our result and that in [72, 103]. The non-singular bouncing models in [72, 103] are constructed with a cubic Galileon action $\gamma G(X)\Box\phi$, in which the denominator in the expression of $z_{s}^{2}$ (5.2) become $(HM_{p}-\gamma\dot{\phi}G^{\prime}(X)/2)^{2}$. In this case, the dominant term in the denominator of (5.2) during bouncing phase is $(\gamma\dot{\phi}G^{\prime})^{2}$ instead of $H^{2}M_{p}^{2}$ in our scenario. Without the inclusion of cubic Galileon term, we simply get a constant tensor-to-scalar ratio $r_{\phi}=96$ in our model.

As we show above, the single-field Ekpyrotic bouncing scenario predicts a blue-tilted scalar spectra, so we introduce the curvaton field $\sigma$ to get the scale-invariant power spectrum on CMB scale. The dynamical equation for the fluctuation of curvaton is

$$u_{k}^{\prime\prime}+\left(k^{2}-\frac{z_{c}^{\prime\prime}}{z_{c}}\right)u_{k}=0~{},$$

(5.6)

where $u_{k}\equiv z_{c}\delta\sigma$ is the mode function for the curvaton field and

$$z_{c}\equiv a\sqrt{\frac{M_{p}^{2}f^{2}}{\cosh(\sqrt{2/q}\phi)}}\simeq\frac{aM_{p}f}{\sqrt{e^{-\sqrt{2/q}\phi}}}=aM_{p}fe^{\frac{\phi}{\sqrt{2q}}}~{},~{}\phi\ll-1~{}.$$

(5.7)

With the help of (2.15) and (2.17), we have

$$e^{\frac{\phi}{\sqrt{2q}}}=\frac{\sqrt{q(1-3q)}}{\sqrt{2V_{0}}M_{p}t}=-\frac{\sqrt{q(1-3q)}}{\sqrt{2V_{0}}M_{p}a_{-}(1-q)}\left[\frac{(\tau_{e}-\tau)}{(\tau_{e}-\tau_{-})^{q}}\right]^{-\frac{1}{1-q}}~{}.$$

(5.8)

The background geometry (2.16) then implies

$$z_{c}=\frac{-f\sqrt{q(1-3q)}}{\sqrt{2V_{0}}(1-q)}(\tau_{e}-\tau)^{-1}~{},$$

(5.9)

so we get the following dynamical equation

$$u_{k}^{\prime\prime}+\left(k^{2}-\frac{2}{(\tau_{e}-\tau)^{2}}\right)u_{k}=0~{}.$$

(5.10)

Imposing the vacuum initial condition, the solution to (5.10) is

$$u_{k}(\tau)=\frac{\sqrt{\pi(\tau_{e}-\tau)}}{2}H_{3/2}^{(1)}[k(\tau_{e}-\tau)]=\frac{e^{ik(\tau_{e}-\tau)}}{\sqrt{2k}}\left[1+\frac{i}{k(\tau_{e}-\tau)}\right]~{},$$

(5.11)

where $H_{3/2}^{(1)}$ is the Hankel function of the first kind. For modes that become super-horizon in the contraction epoch, the corresponding power spectrum at the end of contraction phase is

$$P_{\sigma}(k,\tau_{-})=\frac{k^{3}}{2\pi^{2}}\frac{|u_{k}|^{2}}{z_{c}^{2}}\simeq\frac{V_{0}(1-q)^{2}}{2\pi^{2}f^{2}q(1-3q)a_{-}^{2}}~{}.$$

(5.12)