Effective energy density of relic gravitons

The energy density, the pressure and the anisotropic stress of the relic gravitons do not have a unique gauge-invariant and frame-invariant expression since the equivalence principle ultimately forbids the localization of the energy-momentum of the gravitational field. All the strategies proposed so far suggest compatible descriptions inside the effective horizon but they lead to sharply different answers when the wavelengths of the gravitons exceed the Hubble radius by remaining shorter than the typical extension of causally connected patches. We present a coherent discussion of the energy-momentum pseudo-tensors of the relic gravitons in the different kinematical regimes with the aim of scrutinizing the mutual consistency of the competing suggestions. The various proposals are systematically compared by deducing and analyzing the explicit form of the observables in realistic physical situations. General lessons are drawn on the most plausible parametrization of the energy-momentum pseudo-tensor of the relic gravitons.


Introduction
The energy and the momentum of the gravitational field cannot be localized [1,2]. In fact, assuming the equivalence principle in its stronger formulation, the laws of physics are those of special relativity in a freely falling (non-rotating) laboratory that occupies a small portion of spacetime. As long as the coordinates can be transformed to a freely falling frame, it is possible to eliminate locally the gravitational field. It is then conceptually difficult to propose a unique and local definition of the gravitational energy density. As a consequence there are no reasons why the energy-momentum tensor of the gravitational field itself should be either unique or covariantly conserved. For the same reason the energy-momentum tensor of the gravitational waves of cosmological origin is not unique. A variety of pseudo-tensors can be concocted and they should be ultimately equivalent at least in the short wavelength limit, i.e. when the typical frequencies are much larger than the rate of variation of the corresponding geometry. This equivalence is not conclusive since there are physical situations where the frequencies of the waves are smaller than the rate of variation of the corresponding background. For instance, if we consider cosmological backgrounds during an accelerated stage of expansion the particle horizon diverges while the event horizon is proportional to the Hubble rate. The wavelengths of the gravitons become larger than the Hubble radius but are still shorter than the typical size of causally connected regions: by a mode being beyond the horizon we only mean that the physical wavenumber k/a is much less than the expansion rate H and this does not have anything to do with causality [3].
Some of the most notable strategies developed through the years can be summarized, in short, as follows. The Landau and Lifshitz [4] proposal is rooted in the second-order corrections of the Einstein tensor supplemented by the observation that Bianchi identities must be valid to all orders in the perturbative expansion. The Brill-Hartle strategy [5] is instead based on the properties of a specific covariant averaging scheme aimed at separating the terms that evolve faster than the rate of variation of the corresponding background. The Brill-Hartle scheme has been used to derive the Isaacson effective pseudo-tensor providing a sound description of gravitational radiation in the high-frequency limit [6] (see also [7]). The suggestion of Ford and Parker [8,9] follows instead from the effective action of the relic gravitons derived by perturbing the gravitational action to second order in the tensor amplitude. Other apparently different strategies are related to the ones mentioned above. For instance the approach of Refs. [10,11] is the Landau-Lifshitz approach appropriately discussed in the case of a cosmological background.
Through the years the various suggestions have been tested in different frameworks either for the solution of the concrete problems of backreaction [11,12] or for the analysis of the implications of the different proposals [13,14,15,16] not necessarily in connection with the cosmological problems. Babak and Grishchuk [17] came up with a possible definition of a true energy-momentum tensor of the gravitational field. By treating gravity as a nonlinear tensor field in flat space-time [18], Ref. [17] claimed a result with all the necessary properties of a true energy-momentum tensor of the gravitational field itself (i.e. symmetry, uniqueness, gauge-invariance and covariant conservation). By taking the results of Ref. [17] at face value the problem of localizing the energy and momentum of the gravitational field would be completely solved. The perspective of Ref. [17] has been subsequently questioned by Butcher, Hobson and Lasenby [19,20,21] who suggested that the proposal of Ref. [17] does not have a definite physical significance. In spite of the reasonable concerns of Refs. [19,20,21], what matters for the present ends is that the geometrical object most closely related to the Babak-Grishchuk suggestion is (again) the Landau-Lifshitz pseudo-tensor [4] as explicitly recognized by the authors of Ref. [17].
In this paper we intend to clarify the analogies and the important differences characterizing the various approaches developed so far. After scrutinizing the limitations and the ambiguities of the diverse proposals some general lessons will be drawn. In the light of a number of reasonable criteria (i.e. gauge-invariance, frame-invariance, positivity of the energy density) a rather plausible strategy to assign the energy-momentum tensor of the relic gravitons is rooted in the original Ford-Parker suggestion [8,9] where the background metric and the corresponding perturbations are treated as independent fields; the effective energy-momentum pseudo-tensor follows by functional derivation of the effective action with respect to the background metric. In view of the general discussion it is practical to separate three complementary aspects of the problem: i) the strategy for the derivation of the energymomentum pseudo-tensor; ii) the averaging scheme; iii) the connection to the observables. This will be the overall logic followed in the present investigation.
The pseudo-tensors explored through the years must implicitly assume an averaging scheme which is often difficult to formulate in general terms. As long as the relic gravitons potentially populating the present universe did start their evolution either from a quantum mechanical initial state the expectation values of their energy density and of their pressure can be computed without imposing any extrinsic averaging scheme. Indeed during an inflationary stage of expansion the classical fluctuations are diluted away [22,23,24,25] while the quantum fluctuations reappear continuously so that the relic gravitons are parametrically amplified thanks to the pumping action of the gravitational field itself, a perspective invoked in Refs. [26,27] (see also [28]) even before the formulation of the conventional inflationary paradigm. By following the tenets of the quantum theory of parametric amplification (originally developed in the case of optical photons [29]) a fair estimate of the mean energy density and pressure of the relic gravitons is obtained by averaging the various expressions over the same initial state (e.g. the vacuum). Within each of the various parametrizations of the energy-momentum pseudo-tensor the quantum averages will then be used to compare the competing proposals. The other schemes (like the Brill-Hartle average [5] and its descendants [13,14,16]) reproduce the results of the quantum averaging when the wavelengths are shorter than the Hubble radius but are not defined in the opposite limit.
The present paper is organized as follows. In section 2 the various proposals are presented in a common perspective and with the purpose of easing their mutual comparison. By focussing the attention on the case of cosmological backgrounds the explicit expressions for the energy density, for the pressure and for the anisotropic stress are obtained within the different physical suggestions. The quantum averaging and its basic properties is discussed in section 3. Two physically relevant examples are presented in section 4 with the aim of illustrating the basic features and the patent ambiguities of the different parametrizations. The connection of the various proposals with the observables customarily employed in the analysis of relic graviton backgrounds is discussed in section 5. In the same framework we also evaluate the spectral density in the Brill-Hartle-Isaacson approach by assuming that the tensor amplitude is just an isotropic random field varying in time; we connect the related spectral density to the power spectrum. At the end of section 5 we show that the effective energy-momentum pseudo-tensor can be derived in a different (conformally related) frame even if the obtained results are ultimately frame-invariant. Finally section 6 contains the concluding remarks.

Different answers for a similar question
Let us consider a conformally flat background geometry g µν = a 2 (τ )η µν , where a(τ ) is the scale factor and τ denotes the conformal time coordinate; η µν = diag(1, −1, −1, −1) is the Minkowski metric. The metric fluctuations are introduced as g µν ( x, τ ) = g µν + δ (1) t g µν , where δ (1) t g µν denotes the (first-order) tensor fluctuation. The same notation will be used for any tensor in four-dimensions (i.e. the Ricci or Einstein tensor) so that δ (1) t A µν and δ (2) t A µν will denote the first-and second-order tensor fluctuations of the generic tensor A µν . The firstand second-order tensor fluctuations of the metric and of the square root of its determinant are given by: where h i j is a rank-two tensor in three dimensions which is divergenceless and traceless, i.e. ∂ i h i j = 0 = h i i . The prime will denote a derivation with respect to the conformal time coordinate. With this notation H = a /a where H = a H and H =ȧ/a is the Hubble rate in the cosmic time parametrization [note that a(τ ) dτ = d t]. To avoid the proliferation of superscripts we shall sometimes make explicit the derivation with respect to τ and write h k ∂ τ h k instead of h k h k . It is relevant to mention that the tensor amplitude h ij defined in Eqs. (2.1)-(2.2) is invariant under infinitesimal coordinate transformations; if the tensor amplitude is defined as in Eqs. (2.1)-(2.2) its quadratic combinations will be automatically gauge-invariant.

The effective action of the relic gravitons
The effective energy-momentum pseudo-tensor of the relic gravitons follows from the observation [8,9] that the tensor fluctuations and the background metric can be regarded as independent variables. Neglecting, for simplicity, the presence of the sources the action for the relic gravitons is essentially the Einstein-Hilbert action perturbed to second-order i.e.
where P = √ 8πG and the quantity appearing inside the curly bracket is the Einstein-Hilbert action where the total derivatives have been excluded. The second-order fluctuation implicitly indicated in Eq. (2.3) can also be expressed as where Z αβ = Γ ρ αβ Γ σ ρσ −Γ σ βρ Γ ρ σα and Z αβ denotes the corresponding background quantity. To zeroth-order we have that Z 00 = 0 and Z ij = 2H 2 δ ij . To first-order in the tensor amplitude we have instead δ (1) Z 00 = 0 while δ (1) Z ij = 2H 2 h ij . Finally the explicit second-order contributions are: Inserting Eqs. (2.5) and (2.6) into Eq. (2.4) the effective action of the relic gravitons is: The possible presence of background sources does not change the result of Eq. (2.7). In fact δ (2) t S must always be evaluated by imposing the validity of the background evolution and the tensor modes decouple from the matter fields at least if the anisotropic stress of the sources vanishes. Since the effective action of the relic gravitons in a conformally flat metric is given by Eq. (2.7), their energy-momentum pseudo-tensor can be introduced from the functional derivative of S t with respect to g µν by considering h ij and g µν as independent variables: From Eq. (2.7) the explicit form of Eq. (2.8) becomes: where we used the notation F µν = T (gw) µν to distinguish Eq. (2.9) from the other proposals examined below. The indices of F µν are raised and lowered with the help of the background metric (i.e. F ν µ = g αν F αν ); the energy density and the pressure are: The associated anisotropic stress is traceless (i.e. Π (F ) i i = 0) and it is: In terms of Eqs. (2.10), (2.11) and (2.12) the components of F ν µ are denotes the energy flux. The energy density, the pressure and the energy flux combine in the following identity (2.14) The first term at the right hand side of Eq. (2.14) vanishes because of the evolution of the tensor amplitude following from the extremization of the action (2.7) with respect to h ij . The second term at the right hand side of Eq. (2.14) is not vanishing, in general; but if we regard the energy flux as an operator constructed from the corresponding quantum fields, its expectation value over the vacuum is generally vanishing (see section 3).

The second-order variation of the Einstein tensor
The Landau-Lifshitz strategy that is based on the analysis of the nonlinear corrections to the Einstein tensor consisting, to lowest order, of quadratic combinations of the tensor amplitude h ij . While the derivation of Eq. (2.9) does not require the systematic use of a the evolution of the tensor amplitude the opposite is true in the Landau-Lifshitz framework where the energy-momentum pseudo-tensor L ν µ can be expressed as: Furthermore, since the Bianchi identity ∇ µ G µ ν = 0 must be valid to all orders, we must also demand that δ (2) t (∇ µ G µ ν ) = 0 ultimately leading to a relation analog to Eq. (2.14). From the second-order fluctuations of the Einstein tensor the energy density, the pressure and the anisotropic stress are: t G 00 (and to L 00 ): where D R is the total derivative term. According to the logic of this approach the term h k must be replaced by −2Hh k +∇ 2 h k that follows from the evolution of the tensor amplitude.
As a result of this lengthy but straightforward procedure when Eqs. (2.20) and (2.21) are combined in L 00 a further total derivative term emerges so that the final result for the energy density is: All in all Eq. (2.22) shows, as anticipated, that Eq. (2.16) is determined up to the total derivative term (i.e. D 00 ) and the same happens in the case of the pressure terms whose associated total derivatives are qualitatively similar to D R and D 00 but will not be explicitly reported. Finally the explicit form of the condition δ (2) t (∇ µ G µ ν ) = 0 following from the validity of the Bianchi identity is: In Eq. (2.24) the shifted pressure P (L) gw does not coincide with p (L) gw (see Eq. (2.17)); the same comment holds for Q (L) which differs from S (L) introduced in Eq. (2.19). In explicit terms we have that the shifted pressure and the shifted vector are: In comparison with p (L) gw the value of P (L) gw is shifted by the second-order fluctuations of the Christoffel connection The shifted pressure entering Eq. (2.24) should be regarded as the true physical pressure as it will emerge from the explicit examples of sections 4 and 5. The first term at the right hand side of Eq. (2.24) vanishes because of the evolution of the first-order amplitude. The second term at the right-hand side of Eq. (2.24) vanishes when averaged over the quantum state of the relic gravitons (see the discussion in section 3).

The covariant approach
The covariant approach, in its original formulation, assumes the Brill-Hartle scheme [5] that implicitly selects the frequencies exceeding the Hubble rate. In the covariant approach the metric is decomposed as where ∇ µ denotes the covariant derivative with respect to the background metric g µν ; the indices of h µν are raised and lowered with the help of g µν . Within the approach of Eq. (2.28) the cosmological fluctuations correspond to u µ h µν = 0 where u µ is the fluid four-velocity. In the case of a conformally flat background geometry g µν = a 2 (τ )η µν the conditions u µ h µν = 0 and ∇ µ h µν = 0 imply h µ µ = 0; if a(τ ) is constant the three conditions must all be separately imposed 2 . The tensor amplitude h µν of Eq. (2.28) is related to the tensor amplitude h ij of Eqs. Within the covariant approach the energy-momentum tensor following from the Brill-Hartle average is given by: To compare the covariant approach with the other proposals we insert Eq. (2.29) inside Eq. (2.30) and the result for the various components of B ν µ is: where the energy density, the pressure and the anisotropic stress are: 2 By projecting the condition ∇ µ h µν = 0 along u ν we obtain, for a cosmological background with flat spatial sections, that (∇ µ h µν )u ν = H(g αβ − u α u β ) h αβ where H is the Hubble rate. If we then impose, according to Eq. (2.28) that u µ h µν = 0, the condition (∇ µ h µν )u ν = 0 also demands h µ µ = 0 provided H = 0 [i.e. a(τ ) must not be constant].

Equation (2.30)
is the result of an averaging procedure which excludes by construction the long wavelengths and can therefore be applied only inside the Hubble radius. If Eq. (2.30) is blindly applied beyond the Hubble radius various ambiguities arise and they will be discussed in sections 4 and 5. The covariant approach can be however extended for typical wavelengths larger than the Hubble radius by employing a different averaging scheme. In this case the results applicable to cosmological background geometries will coincide exactly with Eq. (2.9) (see the discussion at the end of section 3).

Mutual relations between the different prescriptions
The expressions of the energy density obtained in the cases examined above do not coincide in general terms. To appreciate the differences it is useful to write their mutual relations: (2.37) The expressions of ρ (F ) gw and ρ (L) gw are very similar but they differ by a crucial term containing H. A similar remark holds in the case of the relation between ρ (B) gw and ρ (F ) gw since they both differ by terms that are negligible beyond the Hubble radius i.e. when the frequency of the graviton is smaller than the rate of variation of the background. The corresponding pressures obey qualitatively similar relations that can be easily deduced from the results given above.

The quantum averaging
The classical and quantum fluctuations of cosmological backgrounds obey the same evolution equations, but while classical fluctuations are given once forever (on a given space-like hypersurface) quantum fluctuations keep on reappearing all the time. If the kinematical and dynamical problems of a decelerated cosmology are fixed by means of a phase of accelerated expansion lasting (at least) 65 efolds, the classical fluctuation are exponentially suppressed during inflation [22,23,24,25] (see also [30,31,32]). At a purely classical level it is then plausible to conclude that any finite portion of the event horizon gradually loses the memory of an initially imposed anisotropy or inhomogeneity so that the metric attains the observed regularity regardless of the initial boundary conditions. Since in this situation the power spectra of the scalar and tensor modes of the geometry follow from the quantum mechanical expectation values of two field operators evaluated at the same time (but at different spatial locations), it is also very reasonable to apply the quantum averaging for the estimates of the expectation values of the different components of the pseudo-tensors derived above. In the quantum description the classical fields and their derivatives are promoted to the status of quantum mechanical operators i.e. h ij →ĥ ij and h ij →ĥ ij : ij (k) (with λ = ⊕, ⊗) accounts for the two tensor polarizations 3 and the mode functions (for each separate polarization) obey The sum over the polarizations is given by λ e λ i j e (λ) m n (k) = 4S i j m n (k) and S i j m n is defined as is the traceless projector. The field operators (3.1) and (3.2) consist of a positive and of a negative frequency part, i.e.ĥ ij (x) =ĥ ij (x). If |vac is the state that minimizes the tensor Hamiltonian when all the modes are inside the effective horizon (for instance at the onset of inflation) the operator h ij (x) = 0]. The two-point functions associated withĥ ij andĥ ij are therefore given by: and j 0 (kr) are spherical Bessel functions of zeroth order [33,34]; P T (k, τ ) is the standard tensor power spectrum while Q T (k, τ ) is usually not discussed but its presence is essential in the present context: In Eqs. (3.5)-(3.6) the expectation values have been computed over the vacuum state. The averages of the field operators can also be obtained directly in Fourier space from the corresponding Fourier transforms; by representing the quantum operators in Fourier space asĥ the explicit expressions ofĥ ij ( k, τ ) and ofĥ ij ( k, τ ) follow from Eqs. (3.1) and (3.2) so that the corresponding expectation values are The expressions of Eqs. (3.9) and (3.10) hold for quantum mechanical operators but can be easily viewed as classical expectation values of isotropic random fields, as we shall discuss in section 5. The expectation values of the energy density and of the pressures in the different parametrizations examined above will be computed in the remaining part of this section by using the following notations 4 : where X = F, L, B. From Eq. (3.11) it is also practical to introduce the spectral energy density and the spectral pressure defined as The quantum averaging implies a correct ordering of the operators: for instance the quantum version of the classical expression 2H∂ τ h k h k reads, as usual, H(∂ τĥk ĥ k +ĥ k ∂ τĥ k ). The choice of |vac is not mandatory: if the vacuum state is replaced by some other initial state the present considerations apply in the same way provided the same initial state is used for all the expectation values in the various descriptions.

The effective energy momentum pseudo-tensor
In spite of the specific parametrization of the energy-momentum pseudo-tensor, it is a general property of the quantum mechanical expectation values that the averages of the anisotropic stresses and of the total derivatives are all vanishing: Similarly the expectation values of the three-divergences of the energy fluxes vanish, i.e.
While the results of Eqs. (3.13) and (3.14)-(3.15) hold for all the cases examined above, the spectral energy densities and the spectral pressures of Eq. (3.12) are all different. In the case X = F (see Eq. (3.12)) the spectral energy density and the spectral pressure become: Recalling the explicit form of the power spectra given in Eq. (3.7), ρ (F ) gw (k, τ ) and p (F ) gw (k, τ ) can also be expressed as: When the typical frequencies of the gravitons are much larger than the rate of variation of the geometry, Eq.
According to Eq. (3.20), when the modes are inside the Hubble radius the barotropic index of the relic gravitons coincides approximately with 1/3 to leading order in H 2 /k 2 1 (i.e. kτ 1). In the opposite limit (i.e. k 2 /H 2 1) the frequency of the waves is much smaller than the rate of variation of the geometry and Eq. (3.3) can be solved by iteration: where τ ex is the conformal time corresponding to the exit of a given wavelength from the Hubble radius (i.e. kτ ex = O(1)). From Eqs. (3.21) and (3.22) we have that in the limit k H If the background expands, the terms proportional to |G k (τ ex )| 2 quickly becomes negligible and the effective barotropic index goes to −1/3; conversely if the background contracts the first term becomes subleading and the effective barotropic index tends asymptotically towards 1: All in all we can then say that when the typical frequency of the gravitons exceeds the rate of variation of the geometry (i.e. k H) the high-frequency gravitons behave as a perfect relativistic fluid and the barotropic index is 1/3. In the opposite limit (i.e. k H) the effective barotropic index becomes −1/3 if the background expands while it becomes 1 if the background contracts. Equation (2.14) can also be averaged term by term and the result will be an evolution equation for the expectation values of the energy density and of the pressure, i.e.
where the contribution of the energy flux originally present in Eq. (2.14) disappears because of Eq. (3.15).

The Landau-Lifshitz pseudo-tensor
The quantum averaging of the Landau-Lifshitz pseudo-tensor leads to the same results of the effective energy-momentum tensor in the high-frequency limit but the results are sharply different when the frequency is smaller than the Hubble rate. Using the same notations of Eq. (3.11) we have that, in this case, the spectral distributions are: Equations (3.29) and (3.30) give the pressure and the shifted pressure respectively; note that P (L) gw enters the conservation equation obeyed by the mean values:

The Brill-Hartle scheme and the quantum averaging
The quantum average of the Brill-Hartle-Isaacson results of Eqs. (2.32) and (2.33) leads to the following spectral energy and pressure: (3.33) As usual in the limit k H we have that p (B) gw (k, τ )/ρ (B) gw (k, τ ) → 1/3 since |G k (τ )| 2 = (k 2 + H 2 )|F k (τ )| 2 when the corresponding wavelengths are shorter than the Hubble radius. In the opposite limit, however, Eq. (3.33) leads to a bizarre result: the energy density is asymptotically vanishing (i.e. ρ (B) gw (k, τ ) → 0) and the spectral pressure becomes much larger than ρ (B) gw (k, τ ) (i.e. p B (k, τ )/ρ B (k, τ ) 1) and it is formally divergent. For contracting backgrounds the opposite is true: the spectral pressure gets progressively more negligible (i.e. p (B) gw (k, τ ) → 0) so that p B (k, τ )/ρ B (k, τ ) 1. These apparent inconsistencies (further explored in the concrete examples of sections 4 and 5) can be expected as long as Brill-Hartle average automatically selects all the modes that are inside the Hubble radius and it is therefore not surprising that they lead to quantitative ambiguities when the wavelengths exceed the Hubble radius. It is possible to obtain a covariant expression that also applies beyond the Hubble radius. In this case, however, the energy density and pressure do not follow from the Brill-Hartle scheme. To prove this statement let us neglect, for simplicity, the potential sources and write the full second-order action in the covariant case: where the conditions of Eq. (2.28) have been consistently imposed. By extremizing the action (3.34) with respect to the variation of h µν we obtain the evolution equation for the covariant tensor amplitude If we now consider the background metric and the fluctuating amplitude as independent variables the energy-momentum tensor following from Eq. (3.34) reads: If we now apply the tenets of the Brill-Hartle procedure [5] the covariant gradients average out to zero. Therefore we can flip the covariant derivative from one amplitude to the other. If we do this with the terms inside the squared bracket of Eq. (3.36) we can obtain terms like h αβ ∇ ρ ∇ ρ h αβ . Using Eq. (3.35) all these terms will produce various Riemann tensors that will be neglected so that, at the very end, the only term surviving the average will be the first contribution of Eq. (3.36), i.e. ) and express it in the conformally flat case (i.e. g µν = a 2 η µν and h ij = −a 2 h ij ) we obtain, after a lengthy but straightforward calculation, the same expression of the effective energy-momentum tensor F µν given in Eqs. (2.10)-(2.11).

The spectral energy density
The salient properties of the different pseudo-tensors in the case of expanding backgrounds are summarized, for the sake of conciseness, in Tab. 1. While the basic features of F µν have been deduced without assuming any specific evolution of the background, the physical properties of L µν in the long wavelength limit demand a more concrete analysis of some specific examples.  Table 1: Summary of the salient properties of the different pseudo-tensors in the case of expanding backgrounds; w gw denotes the ratio of the spectral pressure and of the spectral energy density in the different cases.
backgrounds 5 but in what follows the ambiguities of Tab. 1 will be addressed mainly in the case of expanding backgrounds. When the wavelengths are all inside the Hubble radius, the frequency range of the spectrum roughly ranges between the aHz and 100 MHz (i.e. 10 −18 Hz and 10 8 Hz). In backreaction problems (see e.g. [11,12]) the averaged energy density and pressure beyond the Hubble radius are determined by integrating the spectral energy density and the spectral pressure over d ln k between the fixed extrema k ex and k re corresponding to the wavelengths that exit and reenter the Hubble radius.

The expanding branch of the de Sitter space-time
If de Sitter space is exact (i.e. in the absence of slow-roll corrections) the scale factor is given by a i (τ ) = (−τ 1 /τ ) with H = −1/τ ; the scalar modes are absent but the propagating tensors are characterized by the following mode function: where the boundary conditions follow from Eq. (3.3). The spectral energy density is in general a function of k and τ but if we introduce x = |kτ | the spectral energy and pressure are both functions of the single dimensionless variable x: where H 1 a 1 ≡ H 1 = 1/τ 1 [recall that a 1 = a(−τ 1 ) = 1]. According to Eq. (4.2) the spectral energy is always positive semi-definite and that the effective barotropic index interpolates between −1/3 (when x 1) and 1/3 (when x 1). Both results have been anticipated in Tab. 1 on the basis of the general considerations of section 3. Using then Eq. (4.2) we have: so that Eq.  .17)) the spectral energy density and the spectral pressure are: According to Eq. (4.4) spectral energy density does not have a definite sign since it is positive inside the Hubble radius but negative outside: in agreement with previous results [10,11,12]. Since in the de Sitter case H 2 = H we also have that the pressure and the shifted pressure coincide, i.e. p (L) gw (x) = P (L) gw (x); the effectice barotropic index is then given by which is formally the same result of Eq. (4.3) with the difference that the signs are inverted: the averaged energy density is negative while the corresponding pressure is positive. Taken at face value the result of Eq. (4.6) violates the weak energy conditions but it is difficult to attribute a profound physical meaning to this occurrence as long as there exist other pseudo-tensors (like F µν ) not violating the weak energy condition. In the Brill-Hartle case Eqs. (2.32)-(2.33) and (3.33) imply that the corresponding spectral energy density and pressure are According to Eq. (4.7) the energy density is positive semidefinite but the effective barotropic index diverges in the limit x → 0:

The expanding de Sitter background matched to radiation
If the mode function normalized during the de Sitter phase but it evolves through radiation the modes not only exit the Hubble radius but they can also reenter. The spectral energy density and pressure can then be expressed in terms of two dimensionless variables 6 i.e. x = kτ and y = kτ 1 . The scale factor for τ ≥ −τ 1 is linear (as in the case of a radiation-dominated regime) a r (τ ) = (τ + 2τ 1 )/τ 1 and it is continuously matched to a i (τ ) = (−τ 1 /τ ). The scale factor and the its rate of variation are both continuous in −τ 1 , i.e. a i (−τ 1 ) = a r (−τ 1 ) and H i (−τ 1 ) = H r (−τ 1 ). With these precisions we have that for τ ≥ −τ 1 the mode functions are given by where H = 1/(τ + 2τ 1 ) and G k = F k . The complex coefficients c ± (k, τ 1 ) appearing in Eq. (4.9) obey |c + (k, τ 1 )| 2 − |c − (k, τ 1 )| 2 = 1 and are given by: c + (y) = e 2 i y (2y 2 + 2 i y − 1) 2y 2 , c − (y) = 1 2y 2 , (4.10) where we introduced the dimensionless variable y = kτ 1 . The spectral energy density and the pressure have an exact expression which is however not so revealing. For instance in the case of F µν we have: gw (x, y) = H 4 1 y 3 8π 2 (x + 2y) 6 2y 4 + 1 2x 2 + 8xy + 8y 2 + 1 + 4xy 2 − 2x + 8y 3 − 2y sin 2(x + y) − 4xy + 6y 2 + 1 cos 2(x + y) , gw (x, y) = H 4 1 y 3 24π 2 (x + 2y) 6 2y 4 + 1 2x 2 + 8xy + 8y 2 + 3 − 2 4x 2 y + 10xy 2 + 3x + 4y 3 + 3y sin 2(x + y) − x 2 8y 2 − 4 + 4x 8y 2 − 1 y + 32y 4 + 2y 2 + 3 cos 2(x + y) , (4.11) The spectral energy density and pressure appearing in Eq. (4.11) depend on the two dimensionless variables x = kτ and y = kτ 1 ; these expressions can be usefully compared with the results of Eq. (4.2) holding in the pure de Sitter case. The frequencies amplified in the transition from de Sitter space-time to the radiation epoch always obey the condition y 1. When x < 1 the amplified frequencies are still smaller than the rate of variation of the geometry: this means that to make sure that the wavelengths are larger than the Hubble radius during the radiation stage, Eq. (4.11) should be expanded for y 1 and for x 1 with the condition y < x. The leading order result of this double expansion is then given by: As expected Eqs. (4.12) and (4.13) imply that the effective barotropic index is −1/3 while the energy density is positive semi-definite. The wavelengths that exited the Hubble radius during the de Sitter phase and reentered during the radiation epoch correspond to the limit y 1 and x 1; as expected, the barotropic index goes in this case to 1/3 even if the approach is not monotonic but oscillating as it can also be argued from Eq. (4.11). Equations (4.12) and (4.13) confirm, once more, the summary of Tab. 1 in the case of F µν .
The results of Eqs. (4.11), (4.12) and (4.13) will now be compared with the analog expressions derived from Eqs. (2.16) and (2.17) in the Landau-Lifshitz approach; the exact results in this case are: gw (x, y) = H 4 1 y 3 8π 2 (x + 2y) 6 2y 4 + 1 2x 2 + 8xy + 8y 2 − 7 + 12xy + 10y 2 + 7 cos 2(x + y) − 2 6xy 2 − 3x + 12y 3 + y sin 2(x + y) , gw (x, y) = H 4 1 y 3 24π 2 (x + 2y) 6 2y 4 + 1 2x 2 + 8xy + 8y 2 − 5 + 12x 2 2y 2 − 1 + 4x 24y 2 − 7 y + 96y 4 − 18y 2 + 5 cos 2(x + y) + 2 12x 2 y + 38xy 2 + 5x + 28y 3 + 5y sin 2(x + y) . (4.14) According to Eq. (4.14), for x 1 and y 1 the effective barotropic index is always 1/3; this conclusion is compatible with Eq. (4.11) in the same physical limit. However in the limit y 1, x 1 and y < x the results are: According to Eqs. (4.15) and (4.16) the spectral energy density is negative while the barotropic index is given by −7/5. If we consider the shifted pressure P (L) gw (k, τ ) (see Eqs. (2.27) and (3.30)) the result is different (4.17) and it leads, as expected, to the more standard (i.e. −1/3) barotropic index. However, while in the case of Eqs. (4.12) and (4.13) the energy density is positive and the pressure is negative (as it is common when the spatial gradients dominate) in the Landau-Lifshitz case the situation is reversed since the energy density is negative and the pressure is positive. Finally, in the case of the Brill-Hartle proposal the energy density is positive semi-definite and the effective barotropic index is 1/3 when x 1 and y 1. However, in the limits y 1, x 1 and y < x we have instead showing that the spectral pressure is much larger than the energy density, as already discussed in the pure de Sitter case and as it follows from the general arguments illustrated after Eq. (3.33). All in all the effective energy momentum tensor obtained from the second-order variation of the action leads to an energy density that is always gauge-invariant and positive semidefinite exactly as argued in Tab. 1. In the Landau-Lifshitz parametrization the weak energy condition is violated while the ambiguities of the Brill-Hartle approach (when applied for frequencies smaller than the rate of variation of the geometry) demand a completion of the energy-momentum pseudo-tensor in the low-frequency limit.

Observables in the concordance scenario
The expectation values of the energy density [i.e. ρ (X) gw with X = F, L, B] lead to the corresponding spectral energy densities in critical units ; Ω (X) gw (k, τ ) together with the power spectra P T (k, τ ) and Q T (k, τ ) are the pivotal observables customarily employed in the concordance scenario to assess the energy density of the relic gravitons. The (less conventional) spectral pressure in critical units can be instead defied as: Ω (X) gw (k, τ ) and Σ (X) gw (k, τ ) will now be computed in the different parametrizations explored so far and in the realistic situation where the evolution begins with a quasi-de Sitter phase, continues through a radiation-dominated epoch and finally arrives at a matter-dominated stage of expansion.

The spectral energy density in critical units
To leading order in H/k < 1 (and even without an explicit form of the mode functions) we have that ; thus the expressions of Ω (F ) gw (k, τ ) and Σ (F ) gw (k, τ ) inside the Hubble radius become: The expressions of Ω (F ) gw (k, τ ) and Σ (F ) gw (k, τ ) for typical wavelengths larger than the Hubble radius are equally immediate since, in this limit, Eqs. (3.21)- (3.22) imply  For wavelengths shorter that the Hubble radius the general results of Eqs. (5.5) and (5.6) are, in practice, the same for all the various prescriptions and the only differences appear from the first correction to the leading-order result; for illustration the next-to-leading order correction of the spectral energy density is reported in two relevant cases: ) All in all, as long as we are inside the Hubble radius, the spectral energy density and pressure are unambiguous and do not crucially change from one pseudo-tensor to the others. The same is not true when the corresponding wavelengths are larger than the Hubble radius.

Explicit results in the concordance scenario
To examine more closely the implications of the different proposals we consider a realistic evolution where a quasi-de Sitter stage of expansion ends at a time −τ r and it is replaced by the radiation-dominated stage: where β = (1 − ) −1 is a numerical factor required for the continuity of the scale factors in the quasi-de Sitter stage and = −Ḣ/H 2 is the conventional slow-roll parameter. The inflationary phase ends for τ = −τ r and the scale factor is normalized as a r (−τ r ) = 1. The evolution dictated by Eq. (5.13) lasts until τ m when the matter dominated stage begins: where a m (τ m ) = a r (τ m ) and a m (τ m ) = a r (τ m ). From Eqs. (5.13) and (5.14) the power spectra can be derived before (i.e. τ < τ m ) and after (i.e. τ > τ m ) the dominance of matter: T (k, τ, τ r , τ m ) = P T (k, τ r ), (5.17) where P T (k, τ r ) denotes the (constant) inflationary power spectrum: )  . Similarly we can also deduce the spectral energy density during the matter stage: 1440 |kτ | 2 P T (k, τ r ), |k τ | 1.

Stationary random process inside the Hubble radius
If the tensor amplitude is not a quantum field operator but it describes a stationary random process the spatial variation can be approximately neglected by focussing on the conformal time dependence: In this case the spectral energy density is only determined by the temporal variation of the tensor amplitude and can be deduced within the Brill-Hartle scheme with the caveat that the result will only apply inside the Hubble radius. If the random process is stationary, by definition the autocorrelation function will only depend on the time difference at which the two amplitudes are evaluated, i.e. h λ (τ )h λ (τ ) = δ λ λ Γ(τ − τ ). The temporal autocorrelation implies that in Fourier space Since the expression of Ω gw (k, τ ) inside the Hubble radius is unambiguous, Eq. (5.27) implies then a specific relation between the power spectrum P T , the spectral amplitude S h and the spectral energy density in critical units: If we pass from the angular frequencies ω to the frequencies ν (and recall that in the natural units adopted here ω = k = 2πν) Eq. (5.28) can also be phrased as: The result of Eq. (5.29) demonstrates that the tensor amplitudes can be considered as isotropic random fields characterized by stationary autocorrelation functions. In this case Eqs. (3.9), (3.10) and (5.26) must be viewed as averages of classical stochastic processes not necessarily related to quantum field operators.

Frame-invariance of the effective action
The effective energy-momentum pseudo-tensor shall now be evaluated in a generalized Jordan frame where the scalar-tensor action reads: where A(ϕ) and B(ϕ) are dimensionless and depend on the scalar field ϕ. The second-order variation of Eq. (5.30) can be easily obtained by repeating the same steps leading to Eqs.
(2.4) and (2.5)-(2.6) and the result is: Equation (5.31) contains comparatively more terms than the analog results valid in the case A → 1 (see e.g. (2.4)) where various contributions disappear and are replaced by a pair of total derivatives that do not affect the final result. After some lengthy but straightforward algebra the explicit form of the second-order action reads: ij ; a J is the scale factor appearing in the J-frame, i.e. G αβ = a 2 J η αβ . The expression inside the squared bracket of Eq. (5.32) vanishes identically since it corresponds to the (ij) component of the background equations derived from the extremization of the action (5.31) with respect to the variation of the metric. By considering the tensor amplitude h (J) ij and the background metric as independent variables the effective energy-momentum tensor in the J-frame follows from Eq. (5.32) and it is: The conformal rescaling A G αβ = g αβ brings the action (5.32) from the J-frame the Einstein frame: where the first equality follows from the conformal rescaling of the background (i.e. AG αβ = g αβ ) while the second equality is implied by the relation between the first-order tensor fluctuations in the two frames (i.e. Aδ (2.7). As a consequence, the energy densities in the two frames are related as: where ρ (E) gw coincides with Eq. (2.10). Since the energy density of a radiation plasma also scales as ρ (J) r = A 2 ρ (E) r , Eq. (5.35) implies that ρ (J) gw /ρ (J) r = ρ (E) gw /ρ (E) r . This observation ultimately implies that the spectral energy density in critical units is the same in the two conformally related frames (i.e. Ω (J) gw = Ω (E) gw ). Let us remark, as we close, that the class of scalar-tensor theories of Eq. (5.30) is purely illustrative and the effective action of the relic gravitons may also inherit further parity-violating contribution [35,36]; in this case a more general form of the effective action has been proposed in [37] and it is relevant for the description of the polarized backgrounds of relic gravitons. This development is however not central to the present considerations.

Concluding remarks
The energy density of the relic gravitons is not univocally and unambiguously defined. The various suggestions proposed so far coincide when the rate of variation of the background geometry is smaller than the frequency of the corresponding gravitons. However, in cosmological backgrounds the rate of variation of the space-time curvature can also exceed the typical frequencies of the gravitons. The energy-momentum pseudo-tensor of the relic gravitons should fulfil four plausible criteria: it should be frame-invariant and gauge-invariant, it should not violate the weak energy condition and it should be derived in general terms, i.e. without explicitly demanding that the rate of variation of the background geometry is either faster or slower than the frequencies of the corresponding gravitons. An energy-momentum pseudo-tensor with these features follows from the effective action of the relic gravitons by considering the tensor fluctuations and the background metric as independent variables. In its simplest realization the effective action coincides with the result of Ford-Parker and it is defined in all the relevant physical regimes. The spectral energy density in critical units derived within this approach is gauge-invariant and frame invariant since its value in two conformally related frame does not change.
If we assume, a priori, that the typical frequencies of the gravitons must exceed the rate of variation of the geometry we are implicitly following the logic of the Brill-Hartle-Isaacson pseudo-tensor whose results are applicable when the wavelengths of the corresponding gravitons are shorter than the Hubble radius. This proposal can be extended to encompass wavelengths larger than the Hubble radius; if this is done the Brill-Hartle-Isaacson result coincides with the effective energy-momentum tensor derived from the second-order variation of the action. Finally the Landau-Lifshitz pseudo-tensor does not assume that the frequencies must exceed the rate of variation of the geometry but its expression depends explicitly on the expansion rates. Hence the actual results for the energy density and of the pressure easily follow from the specific evolution of the mode functions but they are difficult to assess in general terms. In various realistic and semi-realistic situations the energy density computed in the Landau-Lifshitz approach always becomes negative when the typical wavelengths are larger than the Hubble radius. It seems difficult to attribute a profound physical significance to this occurrence: since we explicitly demonstrated that there exist effective pseudo-tensors not leading to a negative energy density, there are no reasons to conclude that relic gravitons must inevitably violate the weak energy condition as they evolve beyond the Hubble radius.
All in all the effective action of the relic gravitons discussed here leads to a computable energy-momentum pseudo-tensor that can be assessed in the asymptotic physical regimes even without a detailed knowledge of the background evolution. In this context the energy density is positive semi-definite and the whole description can be easily extended to a conformally related frame. The other strategies examined in this investigation give reasonable results only when the relevant wavelengths are shorter than the Hubble radius. Even if the present conclusions have been reached in the framework of a quantum mechanical averaging scheme rooted in the properties of the relic gravitons, we argued that the same conclusions can be obtained by considering the tensor amplitudes as isotropic random fields characterized by stationary autocorrelation functions.