Evolution of gravitational waves in the high-energy regime of brane-world cosmology

We discuss the cosmological evolution of gravitational waves (GWs) after inflation in a brane-world cosmology embedded in five-dimensional anti-de Sitter (AdS_5) bulk spacetime. In a brane-world scenario, the evolution of GWs is affected by the non-standard cosmological expansion and the excitation of the Kaluza-Klein modes (KK-modes), which are significant in the high-energy regime of the universe. We numerically solve the wave equation of GWs in the Poincare coordinates of the AdS_5 spacetime. Using a plausible initial condition from inflation, we find that, while the behavior of GWs in the bulk is sensitive to the transition time from inflation to the radiation dominated epoch, the amplitude of GWs on the brane is insensitive to this time if the transition occurs early enough before horizon re-entry. As a result, the amplitude of GWs is suppressed by the excitation of KK-modes which escape from the brane into the bulk, and the effect may compensate the enhancement of the GWs by the non-standard cosmological expansion. Based on this, the influence of the high-energy effects on the stochastic background of GWs is discussed.


Introduction
The stochastic background of gravitational waves (SBGW) generated during the accelerating phase of the early universe is one of the most fundamental predictions of the inflationary scenario and this can provide a direct way to probe the very early universe. In particular, information about extra-dimensions is expected to be imprinted on the SBGW. Motivated by recent developments in particle physics, the possibility that our universe is described by brane with three spatial dimensions embedded in a higher-dimensional space has been extensively discussed. According to this scenario, gravity propagates in the extra spatial dimensions, while the standard model particles are confined to the three dimensional brane. At low-energy scales, four-dimensional general relativity is successfully recovered and the extra-dimensional effects should be fairly small. On the other hand, at high-energy scales, the localization of gravity is not always guaranteed and a significant deviation from the standard four-dimensional theory is expected. If this is true, the SBGW can provide a direct probe of the extra-dimensional effects.
In this letter, we focus on the brane-world model proposed by Randall and Sundrum [1]. In this model, a three dimensional brane is embedded in the five dimensional anti-de Sitter (AdS 5 ) bulk spacetime (see Ref. [2] for a comprehensive review of this brane-world scenario). The important model parameter is the AdS 5 curvature radius ℓ that determines the scale on which Newton's force law is modified. Table-top experiments on Newton's force law impose the constraint on ℓ that ℓ < 0.1mm [3]. In the context of cosmology, the evolution of the universe is significantly modified by extra-dimensional effects when the Hubble horizon becomes shorter than the AdS 5 curvature scale, Hℓ > 1. Thus there is a critical frequency corresponding to the wavelength of the gravitational waves (GWs) that cross the horizon when Hℓ = 1. This frequency is given by f crit = 0.2mHz (0.1mm/ℓ) 1/2 [4]. Since the GWs with short wavelengths re-enter the horizon at high energies, it is expected that the SBGW above the critical frequency is crucially affected by the extra-dimensional effects. Theoretically, there are two main effects in the high-energy regime: (i) the cosmological expansion becomes slow compared with that in fourdimensional general relativity, which enhances the amplitude of the SBGW and (ii) the excitation of Kaluza-Klein modes (KK-modes) suppresses the amplitude of the GWs on the brane. An interesting point is that these two effects affect the amplitude of the SBGW in an opposite way. Thus, quantitative calculations of the effects of the KK-modes excitation are essential in predicting the spectrum of SBGW.
There have been many attempts to understand the effect of the KK-mode excitation. One way is to work with Gaussian normal coordinates [5,6]. Using these coordinates, we have performed a numerical study of the high-energy effects and showed that the excitation of KK-modes dominates over the effect of the non-standard cosmic expansion [7]. However, the Gaussian normal coordinates have a coordinate singularity and the analysis is limited to the relatively low-energy universe (Hℓ < 1). Recently, Ichiki and Nakamura have performed a numerical calculation in the high-energy regime (Hℓ > 1) using a null coordinate system [8,9]. While they reached the same conclusion as our previous study, they took the initial condition so that the perturbation is constant along a null hypersurface. At present, however, it remains unclear whether the initial condition imposed by them is appropriate for the one determined from inflation.
In this letter, we use the Poincaré coordinate system to perform numerical calculations. The main advantage of this coordinate system is that it covers the whole bulk spacetime. Using these coordinates, several analytic and semianalytic methods ahve been proposed [10,11,12], which work well in the lowenergy regime. Also, the quantum fluctuations in the inflationary epoch have been discussed in the Poincaré coordinate system [13,14]. Thus, we can set more plausible initial conditions from the inflationary universe than those in the Gaussian normal coordinates.
We set up the basic equations and prepare the numerical calculation in Section 2. The numerical results are presented in Section 3. We first discuss the validity of the initial conditions. Then we attempt to construct the spectrum of SBGW in the high energy regime of the universe. Finally, Section 4 is devoted to summary and discussion.

Background metric and evolution equation
In the Randall-Sundrum single-brane model [1], a three-dimensional brane is embedded in five dimensional anti-de Sitter spacetime (AdS 5 bulk). In this Letter, we specifically consider the AdS 5 bulk without a black hole mass and assume that the matter content on the brane is simply described by a homogeneous and isotropic perfect fluid whose equation of state satisfies p = wρ.
In our previous study [7], the Gaussian normal (GN) coordinate system was used to solve the wave equation of GWs. In this coordinate system, the brane is located at a fixed position y = 0. The lapse function n(t, y) is related to the warp factor a(t, y) through the relation, n(t, y) =ȧ(t, y)/ȧ(t, y = 0): a(t, y) = a 0 (t) cosh n(t, y) = e −y/ℓ + (2 + 3w) where ℓ is curvature scale of the AdS 5 bulk and λ > 0 is the tension of the brane. The quantity a 0 (t) ≡ a(t, y = 0) denotes the scale factor on brane.
Although the GN coordinates (1) are well-behaved near the brane and so it is convenient to impose the boundary condition on the brane, difficulties arise when investigating the behavior of GWs in the bulk due to the coordinate singularity at y = y h (t), where a(t, y h ) = 0 (Eq. (2)). This corresponds to the past null infinity of the AdS 5 spacetime. Furthermore, a space-like t =const. hypersurface approaches null near the coordinate singularity. Hence, a sophisticated treatment of boundary conditions near the singularity is required. As a result, the previous numerical investigation was restricted to the analysis at relatively low-energy scales.
In this paper, in order to extend our previous study to the analysis in the high-energy regime, we use the Poincaré coordinate system (τ, x, z). In the Poincaré coordinate system, the brane is moving in the static AdS 5 bulk [15]. The metric is given by where h ij is the tensor perturbation satisfying the transverse and the traceless conditions. In this metric, the trajectory of the brane is described as (τ b , z b ): where the variable t is the cosmic time on brane, which has the same meaning of time t as used in the GN metric (1). From the junction conditions, the Friedmann equation and the conservation law become [16]: where the dot denotes a derivative with respect to t and H is the Hubble parameter defined by H ≡ȧ 0 /a 0 . The function T (t) is given bẏ Hereafter, unless otherwise mentioned, we focus on the evolution of GWs during the radiation dominated epoch after inflation, i.e., w = 1/3. In this case, the solutions for the scale factor a 0 (t) and the normalized energy density defined by ǫ(t) ≡ ρ(t)/λ are expressed as (e.g., [17]): The variables a * and t * are numerical constants, whose meanings will be given in next section.
Next consider the tensor perturbation h ij in the metric (4). For convenience, we decompose the quantity h ij in spatial Fourier modes as h ij = h k (τ, z)e ik·xê ij , whereê ij represents a transverse-traceless polarization tensor. Then the evolution equation for the perturbation h k (τ, z) becomes Here we simply omit the subscript k. It is known that the above equation has the following general solution (e.g., [11,13]): 2 denote the Hankel functions of first and second kind respectively and the coefficientsh 1,2 (m) are arbitrary functions of m. The above expression implies that the GWs propagating in the bulk are generally described as a superposition of the zero mode (m = 0) and the KK-modes (m > 0). In the brane-world where the AdS 5 bulk is bounded by the brane, the evolution equation of GWs must be solved by imposing the boundary condition at the brane. The boundary condition at the brane is determined from the junction condition [2]. Imposing Z 2 symmetry on the brane, it is given by

Initial condition
The initial condition for the perturbed quantity h just after inflation is determined by the quantum fluctuations in the inflationary epoch. According to Ref. [16], the GN coordinates (1) are a useful spatial slicing in the inflationary epoch and the KK-modes defined in this slicing are shown to be highly suppressed during the inflation. Thus, the zero-mode solution in the GN coordinates gives a dominant contribution to the metric fluctuation which is given by where C is a normalization constant and η is the conformal time. On the other hand, the mode function given in Poincaré coordinates (10) can be rewritten in terms of the GN coordinate defined with respect to the inflationary brane as [13]: 2 (mη sinh(y/ℓ)) e i ω η cosh(y/ℓ) .
Comparing (12) with (13), we see that the zero-mode solution given in the inflationary epoch cannot be simply expressed in terms of the zero-mode solution in the Poincaré coordinates, indicating that there should be mixtures of KK-modes to express the zero-mode solution in the inflationary epoch. Nevertheless, in the long-wavelength limit k → 0, both the zero-mode solutions become constant over the time and the bulk space and they coincide with each other. Since we are specifically concerned with the evolution of long-wavelength GWs after inflation, the constant mode, i.e., h = const. and dh/dτ = 0, seems a natural and a physically plausible initial condition for our numerical calculation in the Poincaré coordinate.
However, a subtle point arises when we consider the evolution of GWs in the radiation dominated epoch. In this case, the mode decomposition becomes generally impossible in GN coordinates due to the non-separable form of the background metric (1)-(3), and the constant mode is not necessarily a solution.
From the viewpoint of the Poincaré coordinate in the AdS 5 bulk, the mixture of KK-modes could be significant in the high-energy regime of the universe and this is even true in the long-wavelength GWs. Indeed, even in the lowenergy regime (ρ/λ ≪ 1), the mixture of KK-modes has been shown to be essential for the recovery of the standard four-dimensional result [11].
Hence, the constancy of the GW amplitudes after inflation cannot be always guaranteed even on super-horizon scales. Depending on the choice of the bulk coordinate, the mode h(t 0 , z) = const. may not be a good approximation to the initial condition for numerical calculation if one tries to impose the initial condition in the radiation dominated epoch. In order to clarify these subtleties, the validity of the initial condition h = const., must be checked. This point will be carefully discussed in section 3.1.

Numerical simulation
To solve the wave equation (9) numerically, one problem is that the computational domain should be finite. We must introduce an artificial cutoff (regulator) boundary in the bulk at z = z reg and impose the boundary condition. Here, we impose the Neumann condition at the regulator boundary, i.e., (∂h/∂z) z=zreg = 0. The location of the boundary is set to z reg = 30-100ℓ, which is far enough away from the physical brane to avoid artificial suppression of light KK-modes. We checked that the amplitude of GWs on the brane is fairly insensitive to the location of regulator boundaries. Further, we stop the numerical calculations before the influence of the boundary condition at z = z reg can reach the physical brane z b . With these treatments, all the results presented in Section 3 are free from the effect of regulator boundary.
The numerical calculations of wave equation (9) were carried out by employing the Pseudo spectral method [18]. To be precise, we adopt a Tchebychev collocation method with Gauss-Lobatto collocation points. To implement this, instead of using the Poincaré coordinates (τ, z) directly, we use the following new coordinates (t, ξ): so that the locations of both the physical and the regulator branes are kept fixed. Adopting this coordinate system, the perturbed quantity h(t, ξ) is first transformed into the Tchebychev space through the relation, h(t, ξ) = N n=0 h n (t)T n (ξ) in a finite and a compact domain, −1 ≤ ξ ≤ 1. Here, the functions T n (ξ) denote the Tchebychev polynomials, defined by T n (ξ) ≡ cos(n cos −1 (ξ)). We then discretize the ξ-axis to the N + 1 points (collocation points) using the inhomogeneous grid ξ n = cos(nπ/N). With this grid, fast Fourier transformation can be applied to perform the transformation between the amplitude h(t, ξ) and the coefficients h n (t). In this letter, we specifically set the collocation point as N = 512 or 1024. The partial differential equation (9) is now reduced to a set of ordinary differential equations for the coefficients h n (t). Hence, one can obtain the time evolution of h n (t) by simply adopting the Predictor-Corrector method based on the Adams-Bashforth-Moulton scheme.

Results
Given the initial condition h(t 0 , ξ) = const., the remaining free parameters in our numerical simulation are the wave number k and the initial time t 0 . For convenience, we set the wave number k to k = a * H * , that is, the GW just crosses the Hubble horizon at the time t = t * (see Eq. (8)). Also, we introduce the quantity s 0 given by s 0 ≡ a(t 0 )H(t 0 )/k, which represents the physical scale of the long-wavelength GW normalized by the Hubble horizon at an initial time t 0 . Thus, the free parameters may be represented by the dimensionless energy at the horizon-crossing time, ǫ * ≡ ǫ(t * ) and the normalized wavelength at initial time, s 0 . In the following, the numerical results are presented for various choices of the parameters (s 0 , ǫ * ).

Sensitivity to the initial condition
Let us first consider the initial condition after inflation and check the validity of the assumption h(t 0 , ξ) = const. for long-wavelength GWs. For this purpose, we plot the time evolution of GWs in the Poincaré coordinate system in Figure  1. The upper panel shows the case of the de Sitter brane by setting the equation of state on the brane to w = −1. The lower panel shows the case of the radiation-dominated Friedmann brane (w = 1/3). In both panels, we set the comoving wave number to k = √ 3/ℓ or ǫ * = 1.0 with s 0 = 100.
In the upper panel of Figure 1, the universe on the brane experiences accelerated cosmic expansion and the wavelength of GWs becomes longer than the Hubble horizon. The resultant GW amplitude remains constant not only on the brane but also in the bulk. On the other hand, in the case of the Friedmann brane (lower), the wavelength of GW becomes shorter than the Hubble horizon at τ ≈ −1. In the bulk, a complicated oscillatory behavior of perturbations was found in the region that is causally connected to the physical brane. This indicates that the excitation of KK-modes occurs even if the wavelength of GWs is still outside the Hubble horizon.
These results are somewhat surprising from the viewpoint of the AdS 5 bulk, because the different behaviors simply arise from the difference in the motion of the brane. While the trajectory of the brane is described by a straight line (d 2 z b /dτ 2 = 0) in the case of the de Sitter brane, the trajectory of the Friedmann brane describes an arc with a non-zero curvature d 2 z b /dτ 2 < 0. The situation might be very similar to the moving mirror problem in an electromagnetic field (e.g., Chap. 4.4 of Ref. [19]), where the acceleration or deceleration of the mirror yields the creation of photons due to vacuum polarization. In our present case, the excitation of the KK-modes arises due to the deceleration of the moving brane.
The results in Figure 1 suggest that the initial condition h(t 0 , ξ) = const. may be validated if we set the initial condition just after the end of inflation. However, the constancy of the long-wavelength mode would not be guaranteed in the case of the radiation dominated epoch, as discussed in section 2.2. This implies that the choice of the initial time t 0 is crucial when setting the initial condition at the radiation dominant epoch. Thus, for quantitative investigation of the GWs generated during inflation, the sensitivity to the choice of the initial time should be examined.
In Figure 2, the dependence of the evolution of GWs on the initial time is plotted by varying the parameter s 0 . Left panels show the snapshots of the amplitude h(τ, z) in the bulk when the wavelength of GWs just becomes five times larger than the Hubble horizon, i.e., a 0 H/k = 5. On the other hand, right panels show the time evolution of GWs projected on the brane. Clearly, in the bulk, the amplitude of GWs is very sensitive to the choice of the parameter s 0 , or equivalently, the initial time t 0 . The resultant wave form away from the physical brane does not show any convergence even in the low-energy case (ǫ * = 0.1). By contrast, on the brane, the GW amplitudes tend to converge if we set the initial time t 0 early enough (or set s 0 large enough).
Although we do not fully understand the reason for this convergence, as far as the GWs on the brane are concerned, the evolution of GW amplitudes becomes insensitive to the choice of the initial time when setting the parameter s 0 large enough, s 0 50, for instance.

Influence of high-energy effects on spectrum of gravitational wave background
Having validated the setup of the initial conditions, we now attempt to clarify the high-energy effects of the GWs and evaluate the spectrum of the SBGW on the brane. To quantify these, we wish to discriminate the influence of KKmode excitation in the bulk from the non-standard cosmological expansion caused by the ρ 2 -term in the Friedmann equation (6). For this purpose, we introduce the reference wave h ref , which is a solution of the four-dimensional wave equation just replacing the scale factor and the Hubble parameter defined in the standard Friedmann equation with those defined in the modified Friedmann equation (Eqs. (6) and (8)): Comparing the numerical simulations with the solution of the wave equation (15), the effect of the excitation of KK-modes can be quantified. In both panels, the horizon re-entry time is set to a 0 = 1. As we increase the energy scale at horizon crossing time, the GW amplitude h 5D becomes significantly reduced compared to the reference wave, h ref . Since the late-time evolution of GWs simply scales as h ∝ 1/a 0 in both h 5D and h ref , the results may be interpreted as the excitation of KK-modes during horizon re-entry, which are caused by an escape of five-dimensional graviton from the brane to the bulk. Note that the normalized energy density at horizon re-entry time, ǫ * is related to the observed proper frequency 2π f = k/a 0 (t today ) as where the critical frequency f crit is defined by ℓH * = 1 or ǫ * = √ 2 − 1, which typically yields f crit = 0.2mHz (0.1mm/ℓ) 1/2 [4]. Thus, one expects from Figure 3 that the deviation from the standard four-dimensional prediction for the spectrum of SBGW becomes more prominent above the critical frequency, f > f crit .
In order to estimate the influence of KK-mode excitation on the spectrum of SBGW, the frequency dependence of the GW amplitudes is examined based on the results in Figure 3. Figure 4 shows the ratio of the amplitudes |h 5D /h ref | as a result of the KK-mode excitation. The ratio is evaluated at the low-energy regime long after the horizon crossing time and is plotted as a function of the frequency f /f crit . As a reference, we also show the convergence properties of the ratio |h 5D /h ref | by varying the initial time; s 0 = 25, 50, 100, 150 and 200. Clearly, the ratio |h 5D /h ref | monotonically decreases with the frequency and the suppression of amplitude h 5D becomes significant above the critical frequency f crit (vertical solid line). Using the data points in the asymptotic region ǫ * ≥ 5 (or H * ℓ ≥ √ 35), we try to fit the ratio of amplitudes with s 0 = 200 to a power-law function. The result based on the least-square method becomes (17) with α = 0.76 ± 0.01 and β = 0.67 ± 0.01, which is shown as a dashed line in Figure 4 1 .
The power-law fit (17) can be immediately translated to the spectrum of SBGW as follows. To evaluate the spectrum, we conventionally use the density parameter of SBGW Ω GW (f ) defined by [20] : where ρ GW is the energy density of SBGW and ρ c is the critical density. Above the critical frequency f > f crit , the cosmological expansion at horizon re-entry time is dominated by the ρ 2 -term in the Friedmann equation (6) (17), if we combine the effects of the KK-mode excitation, the SBGW spectrum becomes nearly flat, i.e., above the critical frequency. This result seems almost indistinguishable from the standard four-dimensional prediction and slightly contradicts with the result by Ichiki and Nakamura [9], who predict a relatively steeper spetrum, Ω GW ∝ f −0.46 from the numerical simulation based on a null coordinate system. At present, the reason for this discrepancy is unclear, however, it might be ascribed to the differences of initial conditions arising from the different choices of the bulk coordinates. Although it is premature to discuss the detectability of the SBGW, the significance of KK-mode excitations should play an important clue to probe the signature of extra-dimensions.

Summary and Discussion
We numerically studied the behavior of GWs in the high-energy regime of the universe in the Randall-Sundrum brane-world model. In contrast to the previous study in the GN coordinates, we used the Poincaré coordinate system in which the brane is moving in the AdS 5 background. We numerically investigated the evolution of GWs in the radiation dominated universe.
We have confirmed that the constant mode, h(z, τ ) =const., in the Poincaré coordinates just after inflation is a plausible initial condition from inflation. In the radiation dominated epoch, however, the constant mode does not hold in the bulk. Thus, in principle, we have to start our numerical calculations just after inflation. Fortunately, however, the resultant amplitude of the GWs on the brane becomes insensitive to the initial time t 0 of our numerical simulation if we set t 0 sufficiently early. In other words, our numerical result is insensitive to the transition time from inflation to the radiation dominated epoch if the perturbation re-enters the horizon long enough after the transition.
Based on these discussions, we turn to focus on the frequency dependence of the high-energy effects in order to predict the spectrum of SBGW. We found that the excitation of KK-modes compensates the effect of non-standard cosmological expansion. As a result, the density parameter of SBGW scales as Ω GW ∝ f 0 above the critical frequency f crit = 0.2mHz (0.1mm/ℓ) 1/2 . The result seems indistinguishable from the prediction in the standard four-dimensional theory and slightly contradicts with the result by Ichiki and Nakamura [9]. The most striking result of our work is that the behavior of the GWs in the bulk sensitively depends on the transition time. This feature implies that the constant initial condition in the Poincaré coordinate does not agree with the constant initial condition on the null hypersurface even on super-horizon scales. In this sense, the result obtained by Ichiki and Nakamura [9] may not necessarily agree with our results. In addition, the present numerical analysis is restricted to the frequency range 0.3 f /f crit 15. It is thus premature to discuss the detectability of the SBGW by extrapolating our numerical results to the high-frequency bands observed via future GW detectors. Nevertheless, the significance of KK-mode excitations above the critical frequency holds the clue to probe the signature of extra-dimensions and/or to constrain on braneworld models. In this respect, a more quantitative and precise prediction for the spectrum of SBGW should deserve a further investigation.
Recently, Kobayashi and Tanaka [21] analytically estimated the effect of KKmode excitation at low-energy scales. According to their result, the amplitude of the GWs on the brane depends on the transition time as well as the energy scale ǫ * . This would not be a contradiction with the present numerical result because we set the initial condition in the high-energy regime Hℓ > 1, where the low-energy approximation used by the authors does not work. Thus, in order to understand our result, we need an analytical study of the KK-mode excitation in the high-energy regime. This is a challenging future work.   (17), where fitting was examined using the data with s 0 = 200 at the asymptotic region ǫ * ≥ 5. The dotted line shows the result by Ichiki and Nakamura [9].