 Review
 Open Access
Stochastic backgrounds of relic gravitons: a theoretical appraisal
 Massimo Giovannini^{1, 2}Email author
 Received: 8 July 2009
 Accepted: 19 February 2010
 Published: 19 February 2010
Abstract
Stochastic backgrounds or relic gravitons, if ever detected, will constitute a prima facie evidence of physical processes taking place during the earliest stages of the evolution of the plasma. The essentials of the stochastic backgrounds of relic gravitons are hereby introduced and reviewed. The pivotal observables customarily employed to infer the properties of the relic gravitons are discussed both in the framework of the ΛCDM paradigm as well as in neighboring contexts. The complementarity between experiments measuring the polarization of the Cosmic Microwave Background (such as, for instance, WMAP, Capmap, Quad, Cbi, just to mention a few) and wide band interferometers (e.g. Virgo, Ligo, Geo, Tama) is emphasized. While the analysis of the microwave sky strongly constrains the lowfrequency tail of the relic graviton spectrum, wideband detectors are sensitive to much higher frequencies where the spectral energy density depends chiefly upon the (poorly known) rate of postinflationary expansion.
PACS codes: 04.30.w, 14.70.Kv, 04.80.Nn, 98.80.k
Keywords
 Stokes Parameter
 Associate Legendre Function
 Tensor Mode
 Spectral Energy Density
 Curvature Perturbation
1 The spectrum of the relic gravitons
1.1 The frequencies and wavelengths of relic gravitons
Terrestrial and satellite observations, scrutinizing the properties of the electromagnetic spectrum, are unable to test directly the evolution of the background geometry prior to photon decoupling. The redshift probed by Cosmic Microwave Background (CMB in what follows) observations is of the order of z_{dec} ≃ 1087 and it roughly corresponds to the peak of the visibility function, i.e. when most of the CMB photons last scattered free electrons (and protons). After decoupling the ionization fraction drops; the photons follow null geodesics whose slight inhomogeneities can be directly connected with the fluctuations of the spatial curvature present before matterradiation decoupling. It is appropriate to mention that the visibility function has also a second (smaller) peak which arises because the Universe is reionized at late times. The reionization peak affects the overall amplitude of the CMB anisotropies and polarization. It also affects the peak structure of the linear polarization. The present data suggest that the typical redshift for reionization is z_{reion} ~11 and the corresponding optical depth is 0.087. The optical depth at reionization actually constitutes one of the parameters of the concordance model (see below Eq. (2.7)).
The temperature of CMB photons is, today, of the order of 2.725 K. The same temperature at photon decoupling must have been of the order of about 2962 K, i.e. 0.25 eV (natural units ħ = c = k_{B} = 1 will be adopted; in this system, K = 8.617 × 10^{5}eV). The CMB temperature increases linearly with the redshift: this fact may be tested empirically by observing at high redshifts clouds of chemical compounds (like CN) whose excited levels may be populated thanks to the higher value of the CMB temperature [1, 2].
The initial conditions for the processes leading to formation of CMB anisotropies are set well before matter radiation equality and right after neutrino decoupling (taking place for temperatures of the order of the MeV) whose associated redshift is around 10^{10}. The present knowledge of particle interactions up to energy scales of the order of 200 GeV certainly provides important (but still indirect) clues on the composition of the plasma.
where g_{ ρ }denotes the weighted number of relativistic degrees of freedom at the onset of the radiation dominated evolution and 106.75 corresponds to the value of the standard model of particle interactions. Fermions and bosons contribute with different factors to g_{ ρ }. By assuming that all the species of the standard model are in local thermodynamic equilibrium (for instance for temperature higher than the top quark mass), g_{ ρ }will be given by g_{ ρ }= 28 + (7/8)90 = 106.75 where 28 and 90 count, respectively, the bosonic and the fermionic contributions.
In Eq. (1.1) it has been also assumed that H ≃10^{5} M_{P} as implied by the CMB observations in the conventional case of singlefield inflationary models. In the latter case, absent other paradigms for the generation of (adiabatic) curvature perturbations, the condition H ≃ 10^{5} M_{P} is required for reproducing correctly the amplitude of the temperature and polarization anisotropies of the CMB. For more accurate estimates the quaside Sitter nature of the inflationary expansion must be taken into account. In the ΛCDM paradigm, the basic mechanism responsible for the production of relic gravitons is the parametric excitation of the (tensor) modes of the background geometry and it is controlled by the rate of variation of spacetime curvature.
In the present article the ΛCDM paradigm will always be assumed as a starting point for any supplementary considerations. The reasons for this choice are also practical since the experimental results must always be stated and presented in terms of a given reference model. Having said this, most of the considerations presented here can also be translated (with the appropriate computational effort) to different models.
Given a specific scenario for the evolution of the Universe (like the ΛCDM model), the relic graviton spectra can be computed. The amplitude of the relic graviton spectrum over different frequencies depends upon the specific evolution of the Hubble rate. The theoretical error on the amplitude increases with the frequency: it is more uncertain (even within a specified scenario) at high frequencies rather than at small frequencies.
The experimental data, at the moment, do not allow either to rule in or to rule out the presence of a primordial spectrum of relic gravitons compatible with the ΛCDM scenario. The typical frequency probed by CMB experiments is of the order of ν_{p} = k_{p}/(2π) ≃ 10^{18} Hz = 1 aHz where k_{p} is the pivot frequency at which the tensor power spectra are assigned. We are here enforcing the usual terminology stemming from the prefixes of the International System of units: aHz (for atto Hz i.e. 10^{18} Hz), fHz (for femto Hz, i.e. 10^{15} Hz) and so on. CMB experiments will presumably set stronger bounds on the putative presence of a tensor background for frequencies (aHz). This bound will be significant also for higher frequencies only if the whole postinflationary thermal history is assumed to be known and specified.
The typical frequency window of wideband interferometers (such as Ligo and Virgo) is located between few Hz and 10 kHz, i.e. roughly speaking, 20 orders of magnitude larger than the frequency probed by CMB experiments. The typical frequency window of LISA (Laser Interferometric Space Antenna) extends down to 10100 μHz. While Virgo and Ligo are now operating, the schedule of LISA is still under discussion. The frequency range of wideband interferometers will be conventionally denoted by ν_{LV}. To compute the relic graviton spectrum over the latter range of frequencies, the evolution of our Universe should be known over a broad range of redshifts. We do have some plausible guesses on the evolution of the plasma from the epoch of neutrino decoupling down to the epoch of photon decoupling. The latter range of redshifts corresponds to an interval of comoving frequencies going from ν_{p} ≃ 10^{18} Hz up to ν_{bbn} ≃ 0.01 nHz ~10^{11}Hz (at most).
• the range ν_{p} <ν <ν_{bbn} will be generically referred to as the lowfrequency domain; in this range the spectrum of relic gravitons basically follows from the minimal ΛCDM paradigm;
• the range ν_{bbn} <ν <ν_{max} will be generically referred to as the highfrequency domain; in this range the spectrum of relic gravitons is more uncertain.
The highfrequency branch of the relic graviton spectrum, overlapping with the frequency window of wideband detectors (see shaded box in the left plot of Fig. 1), is rather sensitive to the thermodynamic history of the plasma after inflation as well as, for instance, to the specific features of the underlying gravity theory at small scales. This is why we said that the theoretical error in the calculation of the relevant observables increases, so to speak, with the frequency.
In Fig. 1 (plot at the right) the electromagnetic spectrum is reported in its salient features. It seems instructive to draw a simple minded parallel between the electromagnetic spectrum and the spectrum of relic gravitons. Consider first the spectrum of relic gravitons (see Fig. 1, plot at the left): between 10^{18} Hz (corresponding to ν_{p}) and 10 kHz (corresponding to ν_{LV}) there are, roughly, 22 decades in frequency. A similar frequency gap (see Fig. 1, plot at the right), if applied to the well known electromagnetic spectrum, would drive us from lowfrequency radio waves up to xrays or γrays. As the physics explored by radio waves is very different from the physics probed by γ rays, it can be argued that the informations carried by low and high frequency gravitons originate from two very different physical regimes of the theory. Low frequency gravitons are sensitive to the large scale features of the given cosmological model and of the underlying theory of gravity. High frequency gravitons are sensitive to the small scale features of a given cosmological model and of the underlying theory of gravity.
The interplay between long wavelength gravitons and CMB experiments will be specifically discussed in the subsection 1.2. The main message will be that, according to current CMB experiments, long wavelength gravitons have not been observed yet. The latter occurrence imposes a very important constraint on the low frequency branch of the relic graviton spectrum of the ΛCDM scenario whose salient predictions will be introduced in subsection 1.3. According to the minimal ΛCDM paradigm a very peculiar conclusion seems to pop up: the CMB constraints on the lowfrequency tail of the graviton spectrum jeopardize the possibility of any detectable signal for frequencies comparable with the window explored by wide band interferometers (see subsection 1.4). The natural question arising at this point is rather simple: is it possible to have a quasiflat lowfrequency branch of the relic graviton spectrum and a sharply increasing spectral energy density at highfrequencies? This kind of signal is typical of a class of completions of the ΛCDM paradigm which have been recently dubbed TΛCDM (for tensorΛCDM). The main predictions of these models will be introduced in subsection 1.5. We shall conclude this introductory section with a discussion of two relevant constraints which should be applied to relic graviton backgrounds in general, i.e. the millisecond pulsar and the bigbang nucleosynthesis constraint (see subsection 1.6).
1.2 Long wavelength gravitons and CMB experiments
where k_{p} is the pivot wavenumber and is the amplitude of the power spectrum of curvature perturbations computed at k_{p}; n_{s} and n_{T} are, respectively, the scalar and the tensor spectral indices. The value of k_{p} is conventional and it corresponds to an effective harmonic ℓ_{eff} ≃ 30. The perturbations of the spatial curvature, conventionally denoted by ℛ are customarily employed to characterize the scalar fluctuations of the geometry since ℛ is approximately constant (in time) across the radiationmatter transition. As it is clear from Eqs. (1.2) and (1.3) there is a difference in the way the scalar and the tensor spectral indices are assigned: while the scaleinvariant limit corresponds to n_{s} → 1 for the curvature perturbations, the scale invariant limit for the long wavelength gravitons corresponds to n_{T} → 0.
so, as anticipated, ν_{p} is of the order of the aHz. The amplitude at the pivot scale is controlled exactly by r_{T}. In the first release of the WMAP data the scalar and tensor pivot scales were chosen to be different and, in particular, k_{p} = 0.05 Mpc^{1} for the scalar modes. In the subsequent releases of data the two pivot scales have been taken to coincide.
The change in determination of the parameters of the tensor background for three different choices of cosmological data sets.
Data  r _{T}  n _{s}  α _{S}  Ω_{Λ}  Ω_{M0}  k_{eq}Mpc 

WMAP5 alone  < 0.58 
 0.050 ± 0.034 

 0.00993 
WMAP5 + Acbar  < 0.54 
 0.048 ± 0.027 

 0.01004 
WMAP5 + LSS + SN  < 0.54 



 0.01008 
WMAP5 + CMB data  < 0.64 



 0.01013 
Same as in Tab. 1 but assuming no running in the (scalar) spectral index (i.e. α_{S} = 0).
Data  r _{T}  n _{s}  α _{S}  Ω_{Λ}  Ω_{M0}  k_{eq}Mpc 

WMAP5 alone  < 0.43  0.986 ± 0.22  0 

 0.00936 
WMAP5 + Acbar  < 0.40 
 0  0.767 ± 0.032  0.233 ± 0.032  0.00944 
WMAP5+ LSS + SN  < 0.20  0.968 ± 0.015  0  0.725 ± 0.015  0.275 ± 0.015  0.00999 
WMAP5+ CMB data  < 0.36  0.979 ± 0.020  0  0.775 ± 0.032  0.225 ± 0.032  0.00922 
The inferred values of the scalar spectral index (i.e. n_{s}), of the dark energy and dark matter fractions (i.e., respectively, Ω_{Λ} and Ω_{M0}), and of the typical wavenumber of equality k_{eq} are reported in the remaining columns. While different analyses can be performed, it is clear, by looking at Tabs. 1 and 2, that the typical upper bounds on r_{T}(k_{p}) range between, say, 0.2 and 0.4. More stringent limits can be obtained by adding supplementary assumptions.
where ϵ measures the rate of decrease of the Hubble parameter during the inflationary epoch. The overdot will denote throughout the paper a derivation with respect to the cosmic time coordinate t while the prime will denote a derivation with respect to the conformal time coordinate τ.
where α_{T} now measures the running of the tensor spectral index.
As already mentioned, among the CMB experiments a central role is played by WMAP [3–7] (see also [8–10] for first year data release and [11, 12] for the third year data release. In connection with [3–7], the WMAP 5year data have been also combined with observations of the Acbar satellite [20–23] (the Arcminute Cosmology Bolometer Array Receiver (ACBAR) operates in three frequencies, i.e. 150, 219 and 274 GHz). The TT, TE and, partially EE angular power spectra have been measured by the WMAP experiment. Other (i.e. non spaceborne) experiments are now measuring polarization observables, in particular there are
as well as various other experiments at different stages of development. Other planned experiments have, as specific target, the polarization of the CMB. In particular it is worth quoting here the recent projects Clover [36], Brain [37], Quiet [38], Spider [39] and EBEX [40] just to mention a few. In the near future the Planck explorer satellite [41] might be able to set more direct limits on r_{T} by measuring (hopefully) the BB angular power spectra.
Following the custom the TT correlations will simply denote the angular power spectra of the temperature autocorrelations. The TE and the EE power spectra denote, respectively, the cross power spectrum between temperature and polarization and the polarization autocorrelations.
1.3 The relic graviton spectrum in the ΛCDM model
Having defined the frequency range of the spectrum of relic gravitons, it is now appropriate to illustrate the possible signal which is expected within the ΛCDM scenario.
where is the critical energy density. In the present review the ln will denote the natural logarithm while the log will always denote the common logarithm.
Since ρ_{crit} depends upon (i.e. the present value of the Hubble rate), it is practical to plot directly (ν, τ_{0}) at the present (conformal) time τ_{0}. The proper definition of Ω_{GW}(ν, τ_{0}) in terms of the energymomentum pseudotensor in curved spacetime is postponed to section 5. The salient features of the relic graviton spectra arising in the context of the ΛCDM scenario can be appreciated by looking carefully at Fig. 2.
Both in Eqs. (1.8) and (1.10) Ω_{M0} and Ω_{R0} denote, respectively, the present critical fraction of matter and radiation with typical values drawn from the best fit to the WMAP 5yr data alone and within the ΛCDM paradigm. In Eq. (1.10) g_{ ρ }denotes the effective number of relativistic degrees of freedom entering the total energy density of the plasma. While ν_{eq} is still close to the aHz, ν_{bbn} is rather in the nHz range. In Fig. 2 (plot at the left) the spectral index n_{T} is frequency independent; in the plot at the right, always in Fig. 2, the spectral index does depend on the wavenumber. These two possibilities correspond, respectively, to α_{T} = 0 and α_{T} ≠ 0 in Eq. (1.6). In the regime ν <ν_{eq} a numerical calculation of the transfer function is mandatory for a correct evaluation of the spectral slope. In the approximation of a sudden transition between the radiation and matterdominated regimes the spectral energy density goes, approximately, as . The spectra illustrated have been computed within the approach developed in [46, 47] and include also other two effects which can suppress the amplitude of the quasiflat plateau, i.e., respectively, the late dominance of the cosmological constant and the progressive reduction in the number of relativistic species. The latter two effects can be estimated analytically (see the final part of section 6) and they are, however, numerically less relevant than neutrino free streaming.
Apart from the modification induced by the neutrino freestreaming the slope of the spectral energy density for ν > ν_{eq} is quasi flat and it is determined by the wavelengths which reentered the Hubble radius during the radiationdominated stage of expansion. The suggestion that relic gravitons can be produced in isotropic FriedmannRobertsonWalker models is due to Ref. [48] (see also [49]) and was formulated before the inflationary paradigm. After the formulation of the inflationary scenario the focus has been to compute reliably the low frequency branch of the relic graviton spectrum. In [50–52] the lowfrequency branch of the spectrum has been computed with slightly different analytic approaches but always assuming an exact de Sitter stage of expansion prior to the radiationdominated phase. The analytical calculation (whose details will be described in section 6) shows that in the range ν_{p} <ν <ν_{eq}, the spectral energy density of the relic gravitons (see Eq. (1.7)) should approximately go as Ω_{GW}(ν, τ_{0}) ≃ ν^{2}. Within the same approximation, for ν > ν_{eq} the spectral energy density is exactly flat (i.e. Ω_{GW}(ν, τ_{0}) ≃ ν^{0}). This result, obtainable by means of analytic calculations (see also [53–56]), is a bit crude in the light of more recent developments. To assess the accurately spectral energy density it is necessary to take into account that the infrared branch is gradually passing from a quasiflat slope (for ν > ν_{eq}) to the slope ν^{2} which is the one computed within the sudden approximation [53–56]. It is useful to quote some of the previous reviews which covered, in a more dedicated perspective, the subject of the stochastic backgrounds of relic gravitons. The review article by Thorne [57] does not deal solely with relic graviton backgrounds while the reviews of Refs. [58–60] are more topical.
where, as in Eqs. (1.2) and (1.3), denotes the amplitude of the power spectrum of curvature perturbations evaluated at the pivot wavenumber ν_{p}. It is worth noticing that between ν_{bbn} and ν_{max} there are approximately 20 orders of magnitude in frequency. In the ΛCDM scenario the spectrum has, in this range, always the same slope (i.e. n_{T} is frequencyindependent in Eq. (1.2)).
Some details of the calculations leading to the spectral energy densities illustrated in Fig. 2 can be found in sections 5 and 6. Without dwelling on the details it is however clear, as anticipated, that the constraints on the long wavelength gravitons make it difficult (if not impossible) to have a detectable spectral energy density at the scale of wideband interferometers. The latter statement, valid in the minimal ΛCDM scenario, will be sharpened in the following subsection.
1.4 Short wavelength gravitons and wideband interferometers
In the ΛCDM scenario the spectral energy density of the relic gravitons has its larger amplitude in the lowfrequency branch. As the frequency increases the spectral energy density diminishes so that it is plausible to expect a rather small amplitude over the frequencies corresponding to wideband interferometers (see, for instance, Fig. 2 for ν ≃ ν_{LV} = 100 Hz).
Wideband interferometers operate in a window ranging from few Hz up to 10 kHz (see also Fig. 1). The available interferometers are Ligo [61], Virgo [62], Tama [63] and Geo [64]. In loose terms these instruments are Michelson interferometers with two important differences: the mirrors are suspended and FabryPérot cavities are used to increase the optical path of the photons. It would be too pretentious to describe in detail, in the present script, also the experimental apparatus and we therefore suggest Ref. [65] where the basics of wideband interferometers are introduced in a selfcontined perspective.
The sensitivity of a given pair of wideband detectors to a stochastic background of relic gravitons depends upon the relative orientation of the instruments. The wideness of the band (important for the correlation among different instruments) is not as large as 10 kHz but typically narrower and, in an optimistic perspective, it could range up to 100 Hz. The putative frequency of wideband detectors will therefore be indicated as ν_{LV}, i.e. in loose terms, the Ligo/Virgo frequency. There are daring projects of wideband detectors in space like the Lisa [66], the Bbo [67] and the Decigo [68] projects. The common feature of these three projects is that they are all spaceborne missions and that they are all sensitive to frequencies smaller than the mHz (1 mHz = 10^{3} Hz). While wideband interferometers are now operating and might even reach their advanced sensitivities during the incoming decade, the wished sensitivities of spaceborne interferometers are still on the edge of the achievable technologies.
In Fig. 3, the common logarithm of the spectral energy density is illustrated as a function of the common logarithm of r_{T}.
The variable β is used in Eq. (1.12) just because this is the notation endorsed by the Ligo collaboration and there is no reason to change it. At the same time, in the present review, β will be used also with different meanings. In section 6, β quantifies the theoretical error on the maximal frequency of the relic graviton spectrum(see e.g. Eq. (6.48) and discussion therein). In section 7, β parametrizes a portion of the azimuthal structure of the Stokes parameters. Since none of these variables appear in the same context, potential clashes of conventions are avoided.
The parametrization of Eq. (1.12) fits very well with Fig. 3 where the pivot frequency ν_{LV} = 100 Hz coincides with the pivot frequency appearing in the parametrization (1.12). For the scaleinvariant case (i.e. β = 3 in eq. (1.12)) the Ligo collaboration sets a 90% upper limit of 1.20 × 10^{4} on the amplitude appearing in Eq. (1.12), i.e. Ω_{GW,3}. Using different sets of data (see [69, 71]) the Ligo collaboration manages to improve the bound even by a factor 2 getting down to 6.5 × 10^{5}. Thus Fig. 3 together with the upper limit of Eq. (1.12) shows that the current Ligo sensitivity is still too small to detect the relic graviton background arising within the ΛCDM paradigm.
1.5 Beyond the ΛCDM paradigm and highfrequency gravitons
where T denotes the observation time and SNR is the signal to noise ratio. Equation (1.13) is in close agreement with the sensitivity of the advanced Ligo apparatus [61] to an exactly scaleinvariant spectral energy density [77–81]. Equation (1.13) together with the plots of Fig. 3 suggest that the relic graviton background predicted by the ΛCDM paradigm is not directly observable by wideband interferometers in their advanced version.
CMB observations probe the aHz region of the spectral energy density of Fig. 2. Wideband interferometers probe a frequency range between few Hz and 10 kHz. In both ranges, the signal of the ΛCDM scenario might be too small to be directly detectable.
The comparison of Fig. 2 and 4 suggests, in short, the following subjects of reflection:

the theoretical error in the estimate of the spectral energy density increases with the frequency;

departures from the standard postinflationary thermal history can be directly imprinted in the primordial spectrum of the relic gravitons;

in the incoming decade the observations of wideband interferometers could be analyzed in conjunction with more standard data sets (i.e. CMB data supplemented by largescale structure data and by the observations of type Ia supernovae) to constrain the spectral energy density of the relic gravitons both at small and at high frequencies.
The presence of postinflationary phases stiffer than radiation is, after all, rather natural and this was the original spirit of [82]. We do not know which was the rate of the postinflationary expansion and since guesses cannot substitute experiments it would be productive to use the TΛCDM paradigm as reference model for a unified analysis of the lowfrequency data stemming from CMB and of the highfrequency data provided by wideband interferometers. Already in [82] (see also [83, 84]) a rather special candidate for a postinflationary phase stiffer than radiation was the case when the sound speed equals the speed of light, i.e. the case when the energy density of the sources driving the geometry is dominated by the kinetic term of a (minimally coupled) scalar field. This particular case was also prompted by various classes of quintessence models. A specific example of this dynamics was provided in [85].
A more detailed account of the techniques leading to Fig. 4 will be swiftly presented in section 6 and can be found in [46, 47]. Without going through the details it is however important to stress that the calculations should be accurate enough not only in the highfrequency region but also in the lowfrequency part of the spectrum. Indeed, as stressed above, one of the purposes of the TΛCDM scenario is to convey the idea that lowfrequency and highfrequency measurements of the relic graviton background can be analyzed in a single theoretical framework.
1.6 The millisecond pulsar bound and the nucleosynthesis constraint
where ν_{pulsar} roughly corresponds to the inverse of the observation time during which the pulsars timing has been monitored. The spectral energy densities illustrated in Figs. 2 and 4 satisfy the pulsar timing bound.
The most constraining bound for the highfrequency branch of the relic graviton spectrum is represented by bigbang nucleosynthesis. Gravitons, being relativistic, can potentially increase the expansion rate at the BBN epoch. The increase in the expansion rate will affect, in particular, the synthesis of ^{4}He. To avoid the overproduction of ^{4}He the expansion rate the number of relativistic species must be bounded from above.
where ν_{bbn} and ν_{max} are given, respectively, by Eqs. (6.61) and (8.4). Thus the constraint of Eq. (1.18) arises from the simple consideration that new massless particles could eventually increase the expansion rate at the epoch of BBN. The extrarelativistic species do not have to be, however, fermionic [89] and therefore the bounds on ΔN_{ ν }can be translated into bounds on the energy density of the relic gravitons.
The spectral energy densities illustrated in Figs. 2 and 4 are both compatible with the bigbang nucleosynthesis bound. Thus the bigbang nucleosyntheis constraint does not forbid a potentially detectable signal in the highfrequency branch of the relic graviton spectrum. Potential deviations of the thermal history of the plasma must anyway occur before bigbang nucleosynthesis.
2 The polarization of relic gravitons and of relic photons
2.1 Basic notations
where Ω_{b0}, Ω_{c0}, Ω_{Λ} denote, respectively, the (present) critical fractions of baryons, CDM particles and dark energy; h_{0} fixes the present value of the Hubble rate; n_{s}, as already mentioned in section 1, is the spectral index of curvature perturbations and ϵ is the reionization optical depth.
At the beginning of the previous section we started by stressing analogies and differences between relic gravitons and relic photons. The most important one is that both gravitons and photons carry two polarizations. This observation is important for a quantitative understanding of the present endevours aimed at measuring the Emode and the Bmode polarization of the CMB. In the present section the description of the polarization of the gravitons will be developed by stressing, when possible, the analogy with polarization observables of the electromagnetic field.
2.2 Linear and circular tensor polarizations
The transformation properties of the circular polarization under a rotation in the plane orthogonal to the direction of propagation are closely analog to the transformation properties, under the same rotation, of the polarization of the electromagnetic field. This analogy will now be exploited to introduce the Emode and Bmode polarization.
Before proceeding with the discussion it is appropriate to recall a very basic aspect of rotations which can have, however, some confusing impact of the polarization analysis especially in the case of the tensor modes. Consider, for simplicity, a coordinate system characterized by two basis vectors, i.e. cos ϑ and sin ϑ. If we now perform a clockwise (i.e. righthanded) rotation of the axes and , the rotated basis will be given as in Eqs. (2.22) and (2.23) by replacing and . Some authors, for different reasons, instead of rotating the coordinate system prefer to rotate the polarization vector. If angles are in the righthanded sense for the rotation of the axes, they are in the lefthanded sense for the rotation of the vectors.
2.3 Polarization of the CMB radiation field
Equations (2.24)(2.25) and (2.37) express the fact that the polarization of the graviton and of the radiation field do change for a rotation on the plane orthogonal to the direction of propagation of the radiation (either gravitational or electromagnetic). It is possible to construct polarization observables which are invariant for rotations on the plane orthogonal to the direction of propagation of the radiation: because of their properties under parity transformations they are called Eand Bmodes.
2.4 E and Bmodes
Overall, the existence of linear polarization allows for 6 different power spectra.
In the minimal version of the ΛCDM paradigm the adiabatic fluctuations of the scalar curvature lead to a polarization which is characterized exactly by the condition a_{2, ℓm}= a_{2, ℓm}, i.e. = 0. This observation implies that, in the ΛCDM scenario, the nonvanishing angular power spectra are given by the TT, EE and TE correlations. In the TΛCDM scenario the TT, EE and TE angular power spectra are supplemented by a specific prediction for the Bmode autocorrelation (see section 7).
2.5 Spin2 spherical harmonics
Spherical harmonics of higher spin appear in matrix elements calculations in nuclear physics (see e.g. the classic treatise of Blatt and Weisskopf [92], and, in a similar perspective the book of Edmonds [93]). The comprehensive treatments of Biedenharn and Louk [94] and of Varshalovich et al. [95] can also be usefully consulted.
The spins harmonics have been introduced, in their present form, by Newman and Penrose [96] and their group theoretical interpretation has been discussed in [97]. The spins spherical harmonics have been applied to the discussion of CMB polarization induced by relic gravitons in a number of papers [98–100]. They are rather crucial in the formulation of the socalled total angular momentum approach. Discussions of the spinweighted spherical harmonics in a cosmological context can also be found in [101, 102]. The spin weighted spherical harmonics will now be introduced by following the spirit of Ref. [97] which has been also used, with different conventions, in [98]. In subsection 2.6 the (equivalent) approach of [99, 100] will be more specifically outlined.
where P_{ℓ}(μ) are the Legendre polynomials and (μ) the associated Legendre functions. It is appropriate to mention here that the factor (1)^{ m }(i.e. CondonShortley phase) can either be included in the normalization factor or (as it has been done) in the definition of the associated Legendre functions appearing in Eq. (2.55). When using the recurrence relations of the associated Legendre functions the CondonShortley phase introduces a sign difference every time m is odd. The conventions expressed by Eqs. (2.54) and (2.55) will be followed throughout the present discussion and, in particular, in section 7 where the correlation functions of the Emodes and of the Bmodes will be specifically computed with different techniques.
This time, in , s = 1 since is a quantity of spin weight 1.
where, as already mentioned, . In Eqs. (2.63) and (2.64) there appear only ordinary (i.e. spinweight 0) spherical harmonics. This occurrence suggests a complementary approach to the problem: instead of expanding Δ_{±} ( , τ) in terms of spin2 spherical harmonics, fluctuations of spinweight 0 can be directly constructed (in real space) from Δ_{±} ( , τ) by repeated application of the ladder operators defined in Eqs. (2.52) and (2.53).
As discussed at the end of subsection 2.1 the sign of φ can be flipped. This possibility is not related to a parity transformation and it has to do with the way twodimensional rotations are introduced. This aspect will also be relevant in section 7 for explicit derivations.
The transformation (2.71) implies that the two basis vectors defined in Eq. (2.70) transform as and , i.e. while does not change flips its sign under space inversion. It follows that spaceinversion does not flip the sign of Δ_{E}( ) but it does flip the sign of Δ_{B}( ), i.e. under the transformation (2.71), Δ_{E}( ) → Δ_{E}( ) while Δ_{B}( ) → Δ_{B}( ).
where L_{±} and L_{z} obey the well known commutation relations [L_{±}, L_{ z }] = ∓L_{±} and [L_{+}, L_{}] = 2L_{ z }.
Looking at Eq. (2.79) it is tempting draw a parallel between the (orbital) ladder operators and the ladder operators raising (or lowering) the spin weight of a given function (see Eqs. (2.52) and (2.53). This problem has been discussed and solved in [97]. It is possible to formulate the parallel in terms of a putative O(4) group. Half of the generators will be connected with the orbital angular momentum operators, while the other half will allow to increase (or decrease) the spin weight of a given function. The two sets of generators commute. The operators are not directly, though, the ladder operators stemming from the second set of generators. This has to do with the fact that in Eq. (2.51) the third Euler angle (i.e. γ) has been fixed to zero. The are ladder operators defined within a putative O(4) group in the case γ ≠ 0. When γ → 0 the dependence upon γ drops and we are left with Eqs. (2.52) and (2.53).
2.6 Polarization on the 2sphere
satisfying P_{ ab }= P_{ ba }, and g^{ ab }P_{ ab }= 0, where is a unit vector in the direction (ϑ, φ). The sign of the offdiagonal entries in Eq. (2.81) is opposite with respect to the one obtained in Eq. (2.34). This is just because we want to match with the conventions adopted, for instance, in [100–102]. To avoid possible confusions, furthermore, the Latin indices a, b, c, d, .... run over the twodimensional space.
is the LeviCivita symbol on the 2sphere. Notice that N_{ℓ} differs from defined in Eqs. (2.46) (see also Eqs. (2.63) and (2.64)) by a factor . This difference will be ultimately relevant to relate and .
The two approaches to the spin weighted spherical harmonics described in the present section are equivalent and can be used interchangeably depending upon the specific problem.
3 The action of the relic gravitons
3.1 Secondorder fluctuations of the EinsteinHilbert action
3.2 Lagrangian densities
All the three Lagrangian densities of Eqs. (3.14), (3.16) and (3.17) lead to the same EulerLagrange equations.
3.3 Hamiltonian densities
as it can be explicitly verified by using Eqs. (3.19), (3.20) and (3.25) into Eq. (3.26).
3.4 Evolution equations in different regimes
This form of the evolution equation for the tensor modes is the one required to compute the effects related to the finite value of the anisotropic stress.
4 Quantization of the tensor modes
There are analogies between the quantum state of relic gravitons and the quantum treatment of visible light. Quantum effects are not crucial to treat firstorder interference of the radiation field (i.e. Young interferometry) [104]. Firstorder interference in quantum optics correspond to the calculation of the twopoint function of the relic gravitons. Quantum effects arise, in optics, from secondorder interference, i.e. when computing (and measuring) the interference between the intensities of the radiation field. Secondorder interference effects are associated with the possibilities of counting photons and have been pioneered by HanburyBrown and Twiss in the early fifties [104, 105]. HanburyBrownTwiss interferometry is based on photon counting statistics.
Having said that we are not even close (experimentally) to study graviton counting statistics (as we do it with the photons), second order interference effects would allow, in principle, to assess the coherence properties of relic graviton backgrounds. The quantum state of the relic gravitons can be described in terms of a generalized coherent state usually called squeezed state. Squeezed states can be described in terms of quadrature operators where one of the modes of the radiation field is always broadened by the time evolution, while the other one is squeezed.
4.1 Heisenberg description
Since, by construction, the Hamiltonians of Eqs. (3.19) and (3.20) are related by canonical transformations, the mode functions of Eqs. (4.11) and (4.12) will have both to obey Eq. (4.13). In different terms, the commutation relations between field operators should be preserved by the time evolution and this is equivalent to the Wronskian normalization condition of Eq. (4.13).
4.2 Generalized coherent states of relic gravitons
There is a slight difference in the normalizations adopted between Eqs. (4.9)(4.10) and Eqs. (4.25)(4.26). This difference is due to the fact that, in Eqs. (4.9)(4.10) the mode functions f_{ k }are normalized, asymptotically, in such a way that f_{ k }→ 1/ . In Eqs. (4.15)(4.16) the factors and have been included in the definition creation and annihilation operators.
which has exactly the same physical content of Eqs. (4.32) and (4.33). When the Universe expands, g_{ k }(τ) decreases and that the solution associated with A_{2}(k) becomes progressively subleading. However, this observation does not imply that g_{ k }(τ) disappears since the evolution must be unitary. This feature of squeezed quantum state suggests the possibility of associating an effective entropy to the process of graviton production[106–110].
where the creation and destruction operators are the ones computed in τ_{0}, i.e. by definition of Schrödinger description. The state 0⟩ is annihilated both by and by . These twomodes appear simultaneously since gravitons are produced from the vacuum whose total momentum vanish. The x and have been dubbed, in the literature, as superfluctuant operators (see, e. g., [106–108]).
the colons denote normal ordering and denotes the operator corresponding to the intensity of the radiation field. The normal ordering is related to the fact that, in the optical domain, most measurements of the electromagnetic field are based on the absorption of photons via the photoelectric effect. Needless to say that there is no analog of photoelectric detection for (single) relic gravitons. In this sense the following considerations should be regarded as a conditional predictions based on the analogy between squeezed states of photons and squeezed states of gravitons.
where ⟨ ⟩ = sinh^{2} r_{ k }is the multiplicity. The coherent state leads to a radiation field with Poissonian statistics. Thermal states (as well as squeezed states) have a statistics which is, according to the quantum optical terminology, superpoissonian. The latter statement is often dubbed by saying that if g^{(2)}(0) > 1 photons are bunched while, in the opposite case (i.e. g^{(2)}(0) < 1) the photons are said to be antibunched. The quantum optical language is much more effective for a mathematical description of the semiclassical limit than the usual considerations related to the limit ħ → 0. Squeezed states are genuine quantum states with many particles. They are, in some sense, like coherent states with the crucial difference that their statistics is superPoissonian. The possibility of scrutinizing the statistical properties of manygravitons systems would rely on our ability of resolving single gravitons which is not even close to the present technological capabilities.
5 Relic graviton backgrounds: observables
In the literature relic graviton backgrounds are characterized in terms of different quantities and, in particular, the most common ones are:

the power spectrum (k, τ);

the spectral energy density of the relic gravitons Ω_{GW}(k, τ);

the spectral amplitude (ν, τ).
It is understood that all the mentioned quantities can be expressed either in terms of the wavenumber or in terms of the frequency since k = 2πν.
The three listed variables can be related in different regimes. For instance the power spectrum has a simple relation to the spectral energy density when the relevant wavelengths are inside the Hubble radius. In section 6 it will be argued that, for numerical applications, the transfer function for the spectral energy density is more practical to compute than the transfer function for the power spectrum or for the spectral amplitude itself. The power spectrum is actually a strongly oscillating function of the conformal time coordinate for wavelengths shorter than the Hubble radius (i.e. kτ > 1); in the same limit the spectral energy density is asymptotically constant.
5.1 The tensor power spectrum and the spectral amplitude
where (k) denotes the tensor power spectrum and where the factor 2 in front of the averages arises as a consequence of the appearing in Eq. (5.11). In Eqs. (5.11), (5.12) and (5.13) the conformal time coordinate is absent. In Eq. (5.10) the conformal time appears explicitly. Indeed, Eqs. (5.11), (5.12) and (5.13) tacitly assume that h_{⊕}( , τ) = e_{⊕}( )T_{e}(k, τ) and that h_{⊗}( , τ) = e_{⊗}( )T_{e}(k, τ). This factorization is related to the concept of transfer function for the amplitude which will be discussed in section 6. The decomposition of Eqs. (5.11), (5.12) and (5.13) is useful when all the polarization have to be treated simultaneously typically in problems involving long wavelength gravitons (see Eqs. (7.85)(7.86) and discussion therein). Furthermore the decomposition of Eqs. (5.11), (5.12) and (5.13) allows to factorize the dependence upon the initial spectrum which is useful for numerical applications.
where the second equivalence defines the spectral amplitude (ν) by recalling, once more, that the comoving wavenumber is related to the comoving frequency as k = 2πν.
5.2 Energymomentum tensors for the relic gravitons
The superscript in the energy density and pressure (i.e. and ) is convenient since different prescriptions for assigning the energymomentum pseudotensor will be compared in a moment.
By comparing Eqs. (5.38)(5.39) with Eqs. (5.22)(5.23) we can remark that the first term appearing in Eq. (5.38) is absent from Eq. (5.22). Moreover, also and seems to be superfficially different. As it will be shown in a moment the equivalence of the two approaches is clear as soon as the relevant wavelengths are larger than the Hubble radius at a given time.
5.3 The energy density of the relic gravitons
In the limit kτ > 1 we will have ℋ^{2} ≪ k^{2}. Thus Eqs. (5.50) and (5.51) coincide (up to corrections (ℋ^{2}/k^{2})). In this limit it is also possible to express Ω_{GW}(k, τ) solely in terms of the power spectrum.
As long as the relevant wavelengths are shorter than the Hubble radius at a given time, different prescriptions for assigning the energymomentum pseudotensor lead to the same result (see also the discussion in section 6). In the opposite limit different choices may exhibit quantitative differences. The limit of short wavelengths in comparison with the Hubble radius is the relevant one when discussing wideband interferometers. Conversely, the initial conditions for the CMB anisotropies are set when the relevant wavelengths are larger than the Hubble radius before equality.
6 Relic gravitons from the ΛCDM scenario
6.1 Inflationary power spectra
Within the present notations, as already established in Eq. (3.12), .
Equation (6.8) implies, recalling Eq. (6.5), that r_{T} = 8n_{T}. Since there is a direct relation of the tensor spectral index to r_{T}, the number of the parameters can be reduced from two to one. In Tabs. 1 and 2 the values of r_{T} have been reported as they can be estimated in few different analyses of the cosmological data sets.
6.2 Transfer functions for inflationary power spectra
In Eq. (6.10), (τ) denotes the approximate form of the mode function (holding during the matterdominated phase); F_{ k }(τ) denotes, instead, the solution obtained by fully numerical methods. The averages appearing in Eq. (6.10) refer to the average over the oscillations: as the wavelengths are inside the Hubble radius, the solutions are all oscillating. The numerical average over the phases introduces some arbitrariness which can be cured by computing directly the transfer function for the spectral energy density.
where the typical value selected for is given by the sum of the photon component (i.e. = 2.47 × 10^{5}) and of the neutrino component (i.e. = 1.68 × 10^{5}): the neutrinos, consistently with the ΛCDM paradigm, are taken to be massless and their (present) kinetic temperature is just a factor (4/11)^{1/3} smaller than the (present) photon temperature.
Within the standard approach, Eq. (6.16) is customarily connected to the spectral energy density of the relic gravitons. In [46, 47] it has been observed that it is simpler and more accurate to compute directly the transfer function for the spectral energy density. In the following subsection this procedure will be illustrated in two different cases.
6.3 Transfer function for the spectral energy
To avoid unnecessary complications, the initial condition of the integrations illustrated in Figs. 5 and 6 have been set as (x_{i}) = x_{i}, i.e. the initial spectrum has been rescaled. The transfer function, by definition, must always depend only on the dynamics of the transition and not upon the features (e.g slope, amplitude) of the initial power spectrum.
Equation (6.28) permits the accurate evaluation of the spectral energy density of relic gravitons, for instance, in the minimal version of the ΛCDM paradigm.
Yet another relevant physical situation for the present considerations is the one where the background geometry, after inflation, transits from a stiff epoch to the ordinary radiationdominated epoch. In the primeval plasma, stiff phases can arise for various reasons. Zeldovich [129] (see also [130]) suggested this possibility in connection with the entropy problem. In [82–84, 74] it has been suggested that the stiff phase could take place after the inflationary phase with the main purpose of identifying a potential source of highfrequency gravitons. This possibility was also prompted by a possible postinflationary dominance of a quintessence field.
to a radiationdominated phase where c_{st} = 1/41. Note that, according to Eqs. (6.29) and (6.30), = w_{t} iff the (total) barotropic index is constant in time. In the limiting case w_{t} = 1 = and the speed of sound coincides with the speed of light. As argued in [132], barotropic indices w_{t} >1 would not be compatible with causality (see, however, [133]). The presence of a suitable stiff phase has been also discussed recently as an effective way of suppressing entropic fluctuations [134] which are observationally constrained by the WMAP 5yr data.
where k_{s} = . The value of k_{s} can be computed in an explicit model but it can also be left as a free parameter. Taking into account that the energy density of the inflaton will be exactly , the value of k_{s} (as well as the duration of the stiff phase) will be determined, grossly speaking, by H_{i}/ . In the context of quintessential inflation [85] (see also [83, 84]) ρ_{Ri} ≃ [135].
In Fig. (7) (plot at the left) the full line superimposed to the numerical points (illustrated by boxes) is the fit of Eq. (6.32).
6.4 Analytic results for the mixing coefficients
The analytic results for the mixing coefficients are rather useful to obtain the final expression of the various transfer functions. Indeed, defining as k_{*} the typical wavenumber of the transition (e.g. k_{*} = k_{eq} in the case of the radiationmattter transition), the slope of the transfer function of the spectral energy density can analytically obtained in the limit κ ≫ 1. This observation helps when we have, for instance, to fit the numerical data points with an analytical expression which will however reproduce the data not only for κ > 1 (as Figs. 5 and 7 clearly show).
The logarithms arising in Eqs. (6.43) and (6.44) explain why, in Eq. (6.32), the transfer function of the spectral energy density contains logarithms. In spite of the fact that semianalytical estimates can pin down the slope of the transfer functions in different intervals, they are insufficient for a faithful account of more realistic situations where the slowroll corrections are relevant and when other dissipative effects (such as neutrino fee streaming) are taken into account.
6.5 Exponential damping of the mixing coefficients
Equations (6.45) and (6.46) are derived by assuming that, right after inflation, the radiationdominated phase takes over. Furthermore, recalling the slowroll dynamics, and V ∝ . In Eqs. (6.45) and (6.46) denotes, as already established, the amplitude of the curvature power spectrum evaluated at the pivot scale.
where the typical values of the slowroll parameter have been derived by taking into account that, in the absence of running of the tensor spectral index, r_{T} = 16ϵ; since, according to the WMAP 5yr data alone, r_{T} < 0.43, ϵ ≤ 0.01.
For τ → ∞ (i.e. τ ≪ τ_{i}), a(τ) ≃ a_{i}/τ and the quasi deSitter dynamics is recovered. In the opposite limit (i. e. τ ≫ +τ_{i}), a(τ) ≃ a_{i} τ and the radiation dominance is recovered. In Fig. 6 (plot at the left) the exponential damping of the mixing coefficients is numerically illustrated. The curve at the top (full line) illustrates the case κ = 1. The cases κ = 2 and κ = 3 are barely distinguishable at the bottom of the plot. Notice, always in the right plot, the rather narrow range of times which are reported in a linear scale. In the plot at the right the asymptotic values of the mixing coefficients are reported for different values of κ = k/k_{max}. By fitting the numerical data with with an equation of the form given in Eq. (6.48), the value of β = 6.33. Different examples can be presented on the same line of the one discussed in Fig. 6. While it is pretty clear that the decay is indeed exponential, the value of β may well vary. This can be summarized, for instance, in a rescaling of k_{max}, i.e. by positing, for instance that k_{max} → /β. Thus, the dynamics of the transition can slightly shift the numerical value of the upper cutoff by a numerical factor which depends upon the width of the transition regime.
6.6 Nearly scaleinvariant spectra
where, we recall, H_{0} = 3.24078 × 10^{18} h_{0} Hz.
By comparing Eqs. (6.50)(6.51) to Eqs. (6.56)(6.57), the amplitude for ν ≫ ν_{eq} differs, roughly, by a factor 2. This coincidence is not surprising since Eqs. (6.50)(6.51) have been obtained by averaging over the oscillations (i.e. by replacing cosine squared with 1/2) and by imposing that g_{ k } = kf_{ k }. These manipulations are certainly less accurate than the procedure used to derive the transfer function for the spectral energy density.
At the present time g_{ρ0 }= 3.36 and g_{s0} = 3.90. In general terms the effect parametrized by Eq. (6.67) will cause a frequencydependent suppression, i.e. a further modulation of the spectral energy density Ω_{GW}(ν, τ_{0}). The maximal suppression one can expect can be obtained by inserting into Eq. (6.67) the largest g_{s} and g_{ ρ }. So, in the case of the minimal standard model this would imply that the suppression (on Ω_{GW}(ν, τ_{0})) will be of the order of 0.38. In popular supersymmetric extensions of the minimal standard models g_{ ρ }and g_{ s }can be as high as, approximately, 230. This will bring down the figure given above to 0.29.
All the three effects estimated in the last part of the present section (i.e. free streaming, dark energy, evolution of relativistic degrees of freedom) have common features. Both in the case of the neutrinos and in the case of the evolution of the relativistic degrees of freedom the potential impact of the effect could be more pronounced. For instance, suppose that, in the early Universe, the particle model has many more degrees of freedom and many more particles which can free stream, at some epoch. At the same time we can say that all the aforementioned effects decrease rather than increasing the spectral energy density. Taken singularly, each of the effects will decrease Ω_{GW} by less than one order of magnitude. The net result of the combined effects will then be, roughly, a suppression of Ω_{GW}(ν, τ_{0}) which is of the order of 3 × 10^{2} (for 10^{16} Hz <ν < 10^{11} Hz) and of the order of 4 × 10^{2} for ν > 10^{11} Hz. These figures are comparable with the possible inaccuracies stemming from the calculation of the transfer function and, therefore, this is a further motivation, to use the transfer function of the spectral energy density. Finally the late time effects reduce a quantity which is already pretty small, i.e., as computed, (ν, τ_{0}) ≃ 10^{15} for ν ≫ ν_{eq}.
7 Bmodes induced by long wavelength gravitons
In the minimal realization of the ΛCDM scenario the scalar fluctuations of the geometry induce an Emode polarization which has been observed and which is now subjected to closer scrutiny [3–7] The tensor modes of the geometry not only induce an Emode polarization but also a Bmode polarization. The detected angular power spectra due to the presence of a putative (adiabatic) curvature perturbation are the temperature autocorrelation (TT angular power spectrum) the Emode autocorrelation (EE angular power spectrum) and their cross correlation (i.e. the TE angular power spectrum). The various angular power spectra of the temperature and polarization observables have been already defined in section 2 (see, in particular, Eqs. (2.48)(2.49) and discussions therein). Long wavelength gravitons contribute not only to the TT, EE and TE angular power spectra but also to the Bmode autocorrelations, i.e. the BB angular power spectra. The effect of long wavelength gravitons on the temperature and polarization observables can be studied by deriving the evolution equations of the brightness perturbations which are related, in loose terms, to the fluctuations of the Stokes parameters. The tensor nature of the fluctuation defined in Eq. (2.8) plays, in this respect, a decisive role. In particular the following two points should be borne in mind:

in the case of the scalar modes of the geometry the heat transfer equations have an azimuthal symmetry;

in the case of the tensor modes the fluctuations of the brightness do depend, both, upon μ = cos ϑ as well as upon φ ; this is ultimately, the rationale for the existence of a Bmode polarization.
where is the differential optical depth. In the differential optical depth enters not only the cross section but also the electron concentration and the ionization fraction x_{e}. The notation for the differential optical depth varies: some authors prefer κ' some other . Given the notations used for the conformal time coordinate we will stick to the choice made in Eq. (7.1).
The right left side of Eq. (7.1) constitutes the collisionless term while the right hand side is the collisional contribution. At the righthand side of Eq. (7.3) ℳ(Ω, Ω') is, in general, a matrix whose dimensionality depends upon the specific problem. As it will be shown ℳ(Ω, Ω') can be easily computed from Eq. (7.3). In similar terms f(Ω) should be understood as a column matrix whose components are the various Stokes parameters.
7.1 Collisionless Boltzmann equation for the tensor modes
Equations (7.17) and (7.18) are not symmetric for φ → φ: while Eq. (7.17) is left unchanged, Eq. (7.18) acquires a minus sign. Conversely, the evolution equations of the scalar modes of the geometry are symmetric for φ → φ.
7.2 Azimuthal structure of the collisional contribution
Equations (7.37) and (7.38) is just a useful warmup in view of the realistic situation where:

the components of f(ϑ, φ) are not 2 but 3, i.e. ℐ_{ x }(ϑ, φ) and ℐ_{ y }(ϑ, φ) are supplemented by U (ϑ, φ);

all the 3 components of f(ϑ, φ) do depend upon φ ; in the analog of Eqs. (7.37) and (7.38) on top of the integration over μ' an integration over φ' will appear.
As in Eqs. (7.33)(7.36), the second equality in each of Eqs. (7.45)(7.49) follows immediately from Eqs. (7.27) and (7.28). A final remark on the symmetry properties of the various entries of ℳ(Ω, Ω') is in order:

ℳ_{11}, ℳ_{12}, ℳ_{21}, ℳ_{22} and ℳ_{33} are all symmetric under the simultaneous transformation φ → φ and φ' → φ';

for the same transformation, the remaining entries flip their respective sign.
7.3 Different parametrizations of the full equation
Now the essential steps of the derivation are the following:

Eqs. (7.61)(7.63) must be inserted, respectively, at left hand side of Eqs. (7.54)(7.56);

Eqs. (7.70)(7.72) must be inserted, respectively, at the right hand side of Eqs. (7.54)(7.56).
As in the previous sections the and denote the time derivatives of the polarizations with respect to the conformal time coordinate τ. This notation has been avoided in the previous equations of the present section since it could have been confused with the angular variables describing the polarizations of the outgoing photons. From now on this possible clash of notations does not arise.
Concerning the result obtained in Eqs. (7.77) and (7.78) few comments are in order:

in Eqs. (7.77)(7.78) denotes indifferently either or : the derivation reported in the case of h_{⊕} can be repeated in the case of h_{⊗} bearing in mind the differences in the angular structure (i.e. Eqs. (7.64)(7.66) should be used instead of Eqs. (7.61)(7.63));

by changing φ → φ and φ' → φ', Eqs. (7.61)(7.63) and Eqs. (7.64)(7.66) will be different because if the various sines appearing the various expressions; therefore the angular structure of the Stokes parameters will change but Eqs. (7.77) and (7.78) will keep their form.
In Eq. (7.96) h is the canonical amplitude of the graviton. Recalling the results of section 3 it has been defined that h⊕ = ℓ_{P} h and h_{⊗} = ℓ_{P} h. The factor simplifies once the brightness perturbations of Eqs. (7.94) are used. Using the common practice Eqs. (7.95)(7.97) have been written in units ℓ_{P} = 1.
In [98] the brightness perturbations for the poolarization have been assigned as in Eqs. (7.99)(7.100) and not as Eqs. (7.92) and (7.93). The ladder operators of Eqs. (2.52) and (2.53) are not invariant under the transformation φ → φ. More specifically the spin weight of a given function changes the sign whenever φ → φ. Equations (2.52) and (2.53) are consistent with the brightness perturbations written as in Eqs. (7.99) and (7.100). To be consistent with the customary notations we will therefore adhere to the conventions stipulated in [98].
7.4 Bmodes from relic gravitons
The analytic and numerical solutions of the evolution equations for the brightness perturbations allow for an explicit evaluation of the temperature and polarization observables. The results obtained so far imply that long wavelength gravitons will produce both temperature and polarization anisotropies. More specifically, following the terminology of section 2 (see, in particular, Eqs. (2.48)(2.50)) the relevant angular power spectra induced by the relic gravitons will be the TT, EE, TE and BB angular power spectra.
where the minus sign in the second term at the right hand side of Eq. (7.122) arises since the CondonShortley phase has been included in the definition of the associated Legendre functions (see Eqs. (2.54)(2.55) and discussion therein). The right hand side of Eq. (7.122) now be simplified by using the following steps:

Eq. (7.111) can be used to simplify the term ;

Eq. (7.112) can be used to simplify the term ;

the term (1  μ^{2})P_{j}(μ) can be simplified by using the equation of the Legendre polynomials and the recursion relations (7.111) and (7.112).
As in the case of the Bmode the presence of the factor (i)^{ j }in Eq. (7.110) leads to a sign difference in the terms Δ_{P ℓ ± 2} appearing in Eq. (7.125). The absolute values of the coefficients of the three terms appearing (inside the the square bracket) in Eq. (7.125) are independent on the conventions related to Eq. (7.110). The coefficient of the term Δ_{P ℓ} appears to be different from the homologue term reported in Eq. (4.41) of [100]. After many checks and crosschecks of the derivations leading to Eq. (7.125) it has been concluded that the difference (i.e. 6ℓ(ℓ + 1) instead of 6(ℓ  1)(ℓ + 2) in the numerator of the coefficient of Δ_{P ℓ}) is just the result of a typo and the correct expression is 6(ℓ  1)(ℓ + 2) as in Eq. (7.125).
7.5 Angular power spectra induced by long wavelength gravitons
where the second equality follows from the Fourier transform of Δ_{I}( ) and by using Eq. (7.98) in the obtained expression. The angular integrations can now be performed with the approach already exploited in previous subsection for the polarization observables. It is however useful to solve Eq. (7.131) with a slightly different method (see, e. g. [98]) where explicit use is made of the solutions (7.126) and (7.127).
where, as in Eq. (2.46), . Having discussed the temperature autocorrelation induced by the long wavelength gravitons it is now the moment of discussing the polarization observables. In Fig. 11 (plot at the left) the EE angular power spectra (full line) are compared with the Bmode autocorrelations (dashed line) induced by the tensor modes in the case r_{T} = 0.1 while the other cosmological parameters are fixed to the bestfit values of the WMAp 5yr data alone. Always in Fig. 11 (plot at the right) the full line illustrates the Bmode autocorrelation in the case r_{T} = 0.1 while the dashed line illustrated the EE angular power spectrum stemming from the standard adiabatic mode. In Fig. 12 (plot at the left) the BB angular power spectra are illustrated for different values of r_{T} compatible with current bounds on r_{T} (see Tabs. 1 and 2).
The Emode and the Bmode angular power spectra can be obtained from Eqs. (7.99) and (7.100) by carefully following all the steps leading to the temperature autocorrelation of Eq. (7.135). The crucial difference will be, of course, that and arise, respectively, in the expansion of Δ_{E}( ) and Δ_{B}( ) as reported in Eq. (2.46). To compute and in terms of S_{P}(k, τ) the steps are, in short, the following:
where, Δ_{Tℓ}(k, τ_{0}) and Δ_{Eℓ}(k, τ_{0}) are given, respectively, by Eqs. (7.135) and (7.147). It should finally be mentioned that the derivatives with respect to x appearing in Eqs. (7.147) and (7.148) can be transformed into derivatives with respect to τ by making use of integration by parts [98]. This step is carried on in full analogy with what happens with the scalar modes of the geometry [127].
Various experiments provided, so far, direct limits on the Bmode polarization.
The main polarization experiments, their typical observational frequencies, and the upper limits on the Bmode polarization are illustrated.
Experiments  Years  Pivot frequencies  Upper limits on BB 

20002003  2636 GHz  (2π) ≤ 50(μK)^{2}  
Boomerang [168]  2003  145 GHz  /(2π) ≤ 8.6(μK)^{2} 
2003  140 GHz  /(2π) ≤ 112.3(μK)^{2}  
2007  100/150 GHz  /(2π) ≤ 10(μK)^{2}  
20022005  2636 GHz  /(2π) ≤ 3.76(μK)^{2}  
20042005  3546/84100 GHz  /(2π) ≤ 4.8(μK)^{2}  
20012006  2394 GHz  /(2π) ≤ 0.5(μK)^{2} 
8 Highfrequency spikes in the relic graviton background
In the ΛCDM paradigm long wavelength gravitons can affect the CMB polarization. As the frequency increases towards the region accessible to wide band interferometers, the ΛCDM signal can only decrease for number of independent reasons:

in the ΛCDM paradigm the spectral energy density is only approximately scaleinvariant (see, e.g. Fig. 2) but the scaling violations always tend to make the spectral energy density smaller at high frequencies;

there are various secondary effects associated, for instance, with the variation of the effective number of relativistic degrees of freedom, with the neutrino anisotropic stress and with the transition to a phase dominated by the darkenergy contribution: all these features reduce the spectral energy density in different frequency regions;

the addition of supplementary scalar fields driving the inflationary phase does not change the two previous statements.
The previous conclusions are all based (either directly or indirectly) on the assumption that the thermal history of the Universe is minimal in the sense that, after inflation, the Universe soon becomes dominated by relativistic particles so that the sound speed of the plasma soon reaches values c_{st} = 1/ . The postinflationary thermal history might not be minimal. For instance it could happen that the transition to the radiationdominated regime is not instantaneous [82]. More specifically it can happen that, after inflation, the sound speed of the plasma is such that c_{st} > 1/ . When the sound p speed is larger than 1/ , the fluid is said to be stiff. In the system of units used in the present paper the speed of light is such that c = 1. A natural upper limit for the sound speed is exactly 1 which is the maximally stiff fluid compatible with causality [132] (see, however, also [133]). If the thermal history of the Universe contemplates a postinflationary phase stiffer than radiation, a spike in the relic graviton spectrum is expected at high frequencies [82] (see also [83, 84] as well as Eqs. (6.29)(6.30) and discussion therein). The possibility of having a postinflationary phase stiffer than radiation has been also investigated in different contexts such as in [151–153].
In the early Universe, the dominant energy condition might be violated and this observation will also produce scaling violations in the spectral energy density [154, 155]. If we assume the validity of the ΛCDM paradigm, a violation of the dominant energy condition implies that, during an early stage of the life of the Universe, the effective enthalpy density of the sources driving the geometry was negative and this may happen in the presence of bulk viscous stresses [154, 155] (see also [156, 157] for interesting reprises of this idea). In what follows the focus will be on the more mundane possibility that the thermal history of the plasma includes a phase where the speed of sound was close to the speed of light.
Absent any indirect tests on the thermal history of the Universe prior to the formation of light nuclear elements, it is legitimate to investigate situations where, before nucleosyntheis, the sound speed of the plasma was larger than 1/ , at most equalling the speed of light. In this plausible extension of the current cosmological paradigm, hereby dubbed TensorΛCDM (i.e. TΛCDM) scenario, highfrequency gravitons are copiously produced [46, 47]. Without conflicting with the bounds on the tensor to scalar ratio stemming from the combined analysis of the three standard cosmological data sets (i.e. cosmic microwave background anisotropies, largescale structure data and observations of type Ia supenovae), the spectral energy density of the relic gravitons in the TΛCDM scenario can be potentially observable by wideband interferometers (in their advanced version) operating in a frequency window which ranges between few Hz and few kHz.
The presence of a stiff phase increases the spectral energy density for frequencies larger than a pivotal frequency ν_{s} which is related to the total duration of the stiff phase. If the stiff phase takes place before BBN, then ν_{s} > 10^{2} nHz. If the stiff phase takes place for equivalent temperatures larger than 100 GeV, then ν_{s} ≥ μHz. If the stiff phase takes place for T ≥ 100 TeV, then ν_{s} > mHz. The frequency ν_{s} marks the beginning of a branch of the spectrum where the tilt of the spectral energy density is blue (i.e. increasing in slope) rather than nearly scale invariant or slightly red (as it is the case in the conventional scenario).
8.1 Scaling violations
On top of the standard parameters of the ΛCDM scenario (see, e. g., Eq. (2.7)) the minimal TΛCDM scenario demands two supplementary parameters

the frequency ν_{s} defining the region of the spectrum at which the scaling violations take place;

the slope of the spectrum arising during the stiff phase.
Using Eq. (8.11) and the definition of Σ (see Eq. (8.3)) it turns out that Σ = λ^{1/4}, which is always smaller than 1 and, at most, (1).
The latter requirement would imply that the stiff age is already finished by the time the Universe has a temperature of the order of 100 TeV when, presumably, the number of relativistic degrees of freedom was much larger than in the minimal standard model (in Eq. (8.14) the typical value of g_{ ρ }is the one arising in the minimal supersymmetric extension of the standard model).
8.2 Spectral energy density in the minimal TΛCDM scenario
In Fig. 14 (plot at the right) two examples of the scaling violations on the spectral energy density are illustrated. Both examples are compatible with the bounds illustrated in the plot at the left. Similar examples have been already illustrated in Fig. 4. These examples will now be discussed in greater detail. The spectral energy density has been computed by using the numerical approach presented in section 6.
In the two examples of Fig. 14 (plot at the right) the ΛCDM parameters are fixed to the values reported in Eq. (2.7) (see [3–7]). The spectral index has been allowed to run, i.e. α_{T} ≠ 0 (see Eqs. (1.6) and (6.53)). The two supplementary parameters should be identified with the sound speed during the stiff phase (i.e. c_{st}) and the threshold frequency (i.e. ν_{s}). Besides c_{st} and ν_{s}, there will also be r_{T} which controls, at once, the normalization and the slope of the lowfrequency branch of the spectral energy density. At the moment wide band interferometers have sensitivities which are insufficient for cutting through the phenomenologically interesting region [69–71]. In the near future, however, there is the hope of a dramatic improvement of the sensitivity: even 5 or 6 orders of magnitude at least heeding the original design (see e.g. [61]) together with the recent proposals for an advanced Ligo program.
As specifically discussed in Eqs. (8.5)(8.6) the frequency of the elbow, i.e. ν_{s}, is fully determined by Σ (see Eq. (8.3) and discussion therein). The two supplementary parameters ν_{s} and c_{st} can be traded for Σ and w_{t} as already done in Fig. 14 (plot at the left). In doing so there is also a potential advantage since, according to Eq. (8.4), Σ shifts the maximal frequency of the spectrum.
As soon as the frequency increases from the aHz up to the nHz (and even larger) the spectral energy density increases sharply in comparison with the nearly scaleinvariant case where the spectral energy density was, for ν > nHz, at most (10^{16}). In the case of Fig. 14 the spectral energy density is clearly much larger. The accuracy in the determination of the infrared branch of the spectrum is a condition for the correctness of the estimate of the spectral energy density of the highfrequency branch. The plots of Fig. 14 (see also Fig. 4) demonstrate that the lowfrequency bounds on r_{T} do not forbid a larger signal at higher frequencies.
A decrease of r_{T} implies a suppression of the nearly scaleinvariant plateau in the region ν_{eq} <ν <ν_{s}. At the same time the amplitude of the spectral energy density still increases for frequencies larger than the frequency of the elbow (i.e. ν_{s}). The latter trend can be simply understood since, at high frequency, the transfer function for the spectral energy density grows faster than the power spectrum of inflationary origin. For instance, in the case w_{t} = 1 and neglecting logarithmic corrections, for ν ≫ ν_{s}. Now, recall that n_{T} is given by Eq. (1.4). If r_{T} → 0, the combination (n_{T} + 1) will be much closer to 1 than in the case when, say, r_{T} ≃ 0.3. This aspect can be observed in Fig. 4 where different values of r_{T} have been reported. By decreasing the w_{t} from 1 to, say, 0.6 the extension of the nearly flat plateau gets narrower. This is also a general effect which is particularly evident by comparing the two curves of Fig. 14 (plot at the right).
up to logarithmic corrections. The result of Eq. (8.15) stems from the simultaneous integration of the background evolution equations and of the tensor mode functions according to the techniques described in section 3. The semianalytic estimate of the slope (see [82]) agrees with the results obtained by means of the transfer function of the spectral energy density.
To conclude the discussion it is appropriate to elaborate on the interplay between the stiff spectra and the phenomenological bounds mentioned in subsection 1.6. Let us start with millisecond pulsar bound of Eq. (1.14).
Since ν_{pulsar} ≃ 10^{3} ν_{bbn}, Eq. (8.16) implies that (ν_{pulsar}, τ_{0}) ≃ 10^{13}or even 10^{14} depending upon r_{T}. But this value is always much smaller than the constraint stemming from pulsar timing measurements (see Eq. (1.14)). If either ν_{s} ≫ ν_{bbn} or c_{st} < 1 the value of (ν_{pulsar}, τ_{0}) will be even smaller. This conclusion follows immediately from the hierarchy between ν_{pulsar} and ν_{bbn}. If either ν_{s} ≫ ν_{bbn}, can only grow very little and certainly much less than required to violate the bound of Eq. (1.14). Consequently, even in the extreme cases when the frequency of the elbow is close to ν_{bbn}, the spectral energy density is always much smaller than the requirement of Eq. (1.14).
As noticed in the past [82], the most significant constraint on the stiff spectra stems from BBN. The models illustrated in Fig. 14 (plot at the right) are on the verge of saturating the bounds of Eqs. (1.18)(1.19). This conclusion stems directly from the form of spectral energy density: the broad spike dominates the (total) energy density of relic gravitons which are inside the Hubble radius at the time of big bang nucleosynthesis. For consistency with the lowfrequency determinations of the tensor power spectrum r_{T} must be bounded from above according to the values reported, for instance, in Tabs. 1 and 2. Once the the value of r_{T} has been selected, the constraints of Eqs. (1.18)(1.19) can be imposed. From Fig. 14 a large signal is expected for ν_{LV} ≃ 100 Hz for 0.35 <w_{t} < 0.6. This range turns out to be compatible with the bounds of Eqs. (1.18)(1.19). In the opposite limit (e. g. w_{t} ≃ 1) the spike becomes narrower, the elbow frequency augments and the signal at the interferometer scale diminishes.
8.3 Detectability prospects
The compatibility with the phenomenological constraints demands that the parameters of the TΛCDM paradigm must lie above the full lines of Fig. 14 (plot at the left). The requirements of Fig. 14 suggest, therefore, that Σ should be raised a bit. In this case the frequency of the elbow gets shifted to the right but, at the same time, the overall amplitude of the spike diminishes. The putative amplitude remains still much larger than the conventional inflationary signal reported, for instance, in Fig. 3.