Searching for Primordial Magnetic Fields with CMB B-modes
Levon Pogosian Alex Zucca
astro-ph.CO astro-ph.GA












#1#2#3#4#1 #2, #3 (#4)

    Searching for Primordial Magnetic Fields with CMB B-modes]Searching for Primordial Magnetic Fields with CMB B-modesDepartment of Physics, Simon Fraser University, Burnaby, BC, V5A 1S6, Canada
        
        
        
        
        Was the primordial universe magnetized? The answer to this question would help explain the origin of micro-Gauss strength magnetic fields observed in galaxies. It is also of fundamental importance in developing a complete theory of the early universe.  While there can be other signatures of cosmological magnetic fields, a signature in the cosmic microwave background (CMB) would prove their primordial origin. The B-mode polarization of CMB is particularly promising in this regard because there are relatively few other sources of B-modes, and because the vortical modes sourced by the primordial magnetic field (PMF) survive diffusion damping up to a small fraction of the Silk length. At present, the Planck temperature and polarization spectra combined with the B-mode spectrum measured by the South Pole Telescope (SPT) constrain the PMF strength to be no more than ∼ 1 nano-Gauss (nG). Because of the quartic scaling of the CMB anisotropy spectra with the PMF strength, this bound will not change by much even with the significantly better measurements of the B-mode spectrum by the Stage III and Stage IV CMB experiments. On the other hand, being able to tighten the bound well below the 1 nG threshold is important for ruling out the purely primordial origin (requiring no dynamo action) of galactic fields. Considering Faraday rotation, which converts some of the E-modes into B-modes and scales linearly with the field strength, will help to achieve this goal. As we demonstrate, the upcoming experiments, such as SPT-3G and the Simons Observatory, will be sensitive to fields of ∼ 0.5 nG strength thanks to the mode-coupling signature induced by Faraday rotation. A future Stage IV ground based experiment or a Space Probe will be capable of probing fields below 0.1 nG, and would detect a scale-invariant PMF of 0.2 nG strength without de-lensing or subtracting the galactic rotation measure.
        [
            Levon Pogosian & Alex Zucca
            December 30, 2023
        ===============================
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        § INTRODUCTION
        
        
        The pre-recombination universe was fully ionized and could sustain a frozen-in magnetic field configuration over many epochs of its evolution <cit.>. There are observational and theoretical arguments favouring the existence of such a primordial magnetic field (PMF)<cit.>. Chief amongst them is the unexplained origin of the galactic magnetic fields, typically of micro-Gauss (μG) strength and coherent over the extent of the galaxy. Such magnetic fields could be generated through a dynamo mechanism, but it would still require a seed field of a certain minimum strength <cit.>. Moreover, μG strength fields are observed in high redshift galaxies that are too young to have gone through the number of revolutions necessary for the dynamo to work <cit.>. Alternatively, the galactic fields could trace their origin to turbulent magnetic field amplification during structure formation at high redshift <cit.>, with the small-scale dynamo converting turbulent energy into magnetic energy. Another possibility is that the supernova explosions in protogalaxies provided the magnetic seed fields that were later amplified by compression, shearing and stochastic motions <cit.>.
        A PMF of sufficient strength on appropriate scales could provide the seed, or eliminate the need for dynamical amplification altogether. Since a typical galactic halo is a few kilo-parsec (kpc) in size and forms in a collapse of a few Mega-parsec (Mpc) sized region, the latter possibility requires a PMF of nano-Gauss (nG) present day strength[It is conventional to use comoving units to describe the strength of cosmological magnetic fields,  i.e.𝐁=𝐁_ today = a^2𝐁(a), where a is the scale factor normalized to a=1 today.] coherent over a comoving Mpc scale in order to compress into a μG level galactic field.
        
        PMFs are expected to have been produced at some level in the electroweak (EW) and the QCD phase transitions <cit.>. Detailed studies of hydromagnetic and magnetohydrodynamic (MHD) turbulence <cit.> have demonstrated that magnetic power can transfer from smaller to larger scales through both inverse transfer and inverse cascade, but it is challenging to obtain a coherence scale and the amplitude of the field large enough to explain the galactic fields. For instance, simulations performed in <cit.> show that the PMF generated in the EW phase transition can at best reach a strength of 0.3 nG at a scale of 30 kpc if the initial field is maximally helical <cit.> and if its initial amplitude is the largest allowed by the big bang nucleosynthesis (BBN) bound <cit.>. Another possibility is for the PMF to have been generated during <cit.> and at the end of Inflation <cit.>, although there it is also challenging to generate a PMF of nG strength. As there are still significant gaps in our understanding of the early universe phase transitions and connecting inflation with particle physics, the possibility of an observationally relevant PMF cannot be ruled out <cit.>. In fact, constraints on the PMF can help us select among different theories of the early universe <cit.>. 
        
        The PMF contributes to CMB anisotropies through metric perturbations and through the Lorentz force felt by baryons in the pre-recombination plasma <cit.>. Of particular importance are the vortical perturbations, which lead to B-mode polarization patterns <cit.>. As they contribute to anisotropies on smaller angular scales  (l ∼ 1000), the magnetic vector modes are not as obscured by the galactic foregrounds as the inflationary B-mode. Another source of magnetic B-modes is Faraday rotation (FR) <cit.>, which scales linearly with B and, as we will show, offers a very promising probe of the PMF.
        
        Current CMB bounds on the PMF strength are at a few nG level, and have been at a few nG level for many years. The reason for this seemingly slow improvement, despite the significantly better CMB data, is the fact that the CMB anisotropy spectrum induced by the PMF scales as B^4. Thus, an order of magnitude stronger constraint on the PMF contribution to the CMB spectra translates into just a factor of ∼ 2 improvement in the bound on B.  Prior to the recent measurements of the B-mode polarization spectrum by the South Pole Telescope (SPT) <cit.>, the strongest CMB bounds on the PMF strength were derived by the Planck Collaboration <cit.> using the 2015 Planck data <cit.>. They placed a limit on the PMF strength smoothed over a 1 Mpc region of B_1 Mpc < 4.4 nG at the 95% confidence level (CL). Comparable bounds were obtained in <cit.> using the 7-year WMAP data <cit.> combined with the high-l temperature anisotropy spectrum from SPT <cit.>.  These constraints assumed a stochastic PMF described by a power-law spectrum with the spectral index treated as a free parameter. The bounds become tighter if one assumes a particular spectral index. For instance, the Planck 2015 bound on the scale-invariant PMF, favoured by the simplest models of inflationary magnetogenesis <cit.>, was B_ 1Mpc<2 nG <cit.>.
        
        Recent measurements of the B-mode polarization by POLARBEAR <cit.>, BICEP/Keck <cit.> and SPT <cit.> have opened a new opportunity for probing the PMF<cit.>. The POLARBEAR collaboration placed a bound of B_1 Mpc < 3.9 nG  <cit.> using only their B-mode spectrum with the cosmological parameters fixed to the best fit LCDM values. Later, a joint fit to the 2015 Planck data and the SPT B-mode spectrum <cit.>, performed in <cit.>, reduced the Planck-only bound by a factor of two, pushing it to B_ 1Mpc<2 nG at 95% CL after marginalizing over the spectral index and B_ 1Mpc<1.2 nG for the scale-invariant PMF spectrum. These early results have demonstrated the potential of B-mode measurements for constraining the PMF.  
        
        The bounds mentioned above are based on the analysis of the power spectra of CMB temperature and polarization. Even if the PMF strength B is Gaussian distributed, the PMF stress-energy sourcing the CMB anisotropies is not, as it is quadratic in B. Thus, there can be potentially interesting constraints derived from higher order CMB correlations <cit.>. Magnetically induced CMB bispectra and trispectra have been studied in <cit.>. It is generally found that they do not improve on the bounds derived from power spectra, with a possible exception of the trispectrum sourced by the  passive modes[Passive and compensated modes will be discussed in Sec. <ref>]. As reported in <cit.>, the bound on the trispectrum obtained from the 2013 Planck data <cit.> constrains the amplitude of a scale-invariant PMF to below 0.6 nG if the energy scale at the time of the PMF generation was T_B=10^14 GeV. While potentially promising, one should be aware of some factors that may weaken this bound.  The amplitude of the passive mode trispectrum scales as B^8 (ln T_B/T_ν)^4, where T_ν = 1 MeV is the neutrino decoupling scale, thus the bound on B will weaken if T_B, which is unknown, happens to be smaller. Also, most of the contribution comes from the largest angular scales and any increase in the magnetic spectral index will also weaken the bound. Finally, the result in <cit.> was based on the variance in the trispectrum calculated by Planck <cit.> without taking into account the PMF contribution, which must be included when constraining the PMF signal.
        
        In addition to sourcing anisotropies, the PMF dumps energy into the plasma, influencing the ionization history of the Universe <cit.>. The impact of this effect on the CMB anisotropy spectra offers a potentially more sensitive probe of the PMF than the magnetically induced anisotropies. However, there are remaining uncertainties in modelling the baryon heating due to ambipolar diffusion and decaying magnetic turbulence <cit.>. The PMF can also induce μ- and y-type distortions of the black-body spectrum of the CMB <cit.>, which are not well constrained by existing observations but can be promising future probes of the PMF <cit.>. Similarly, the PMF affects the thermal Sunyaev-Zeldovich contribution to the CMB through modification of the intergalactic medium evolution <cit.>. It has been suggested that PMF would induce small-scale baryonic density fluctuations at the epoch of cosmic recombination, leading to an inhomogeneous recombination process that alters the heights of the peaks in the CMB spectrum <cit.>. This process was recently simulated in <cit.> and very strong bounds (O(0.01-0.05) nG depending on the PMF spectrum) were derived. These results are extremely interesting and deserve follow up investigations.
        
        The search for the PMF is not limited to CMB observations. The magnetic pressure and the diffusion damping of the magnetic fields affect the growth of large scale structure on scales k > 1 h/Mpc <cit.>. Since non-linearities are important on these scales, accurate modelling of these effects requires numerical simulations of structure formation in the presence of the PMF that do not exist at present. There is also preliminary evidence of magnetic fields coherent over cosmological (O(Mpc)) distances coming from the non-observation of GeV gamma-ray halos around the TeV blazars <cit.>. Upcoming data, combined with improved modelling of plasma instabilities, should be able to conclusively confirm or deny this possibility <cit.>. However, even if cosmological scale magnetic fields were discovered, it would not necessarily prove their primordial origin, as they could, in principle, be produced through ejection of magnetic field lines in explosions of highly magnetized stars or via other astrophysical processes <cit.>. On the other hand, detecting a magnetic field signature in CMB would provide a conclusive proof of its primordial origin.
        
        The top scientific objective of the next generation CMB polarization experiments is the detection of inflationary gravitational waves <cit.>. They will also aim to provide an accurate measurement of the weak lensing contribution to the B-mode spectrum, which happens to be on angular scales where the PMF effects are also important. Their specifications will be well suited for reducing the bound on the PMF strength by an order of magnitude <cit.>. In particular, as we will show, measurements of the FR induced mode-coupling correlations by a Stage IV ground based experiment will be able to probe PMFs of 0.1 nG strength and finally push the sensitivity well below the critical 1 nG threshold.
        
        
        
        § PMF CONTRIBUTIONS TO CMB TEMPERATURE AND POLARIZATION ANISOTROPIES
        
        
        Under the ideal magnetohydrodynamic (MHD) approximation, which is known to hold well in the highly conducting primeval plasma of the early universe and on cosmological scales at smaller redshifts, the PMF configuration remains frozen in the plasma. Its strength evolves with time according to 𝐁(𝐱, τ) = 𝐁(𝐱, τ_0)/a^2(τ) (see <cit.> for a review), where τ_0 denotes the present conformal time and a is the scale factor normalized to a(τ_0) =1. Following the convention, we will quote all bounds on the PMF in terms of the “comoving” field strength 𝐁(𝐱, τ_0). 
        
        We assume that the PMF does not have a homogeneous component, as it would break the isotropy of the Universe and has already been strongly constrained by CMB <cit.>. Instead, we model the PMF as a statistically isotropic Gaussian distributed random field described by a power spectrum
        
            ⟨ B_i (𝐤) B_j^*(𝐤^') ⟩ = (2 π)^3 δ^(3)(𝐤 - 𝐤^') P_ij P_B(k),
        
        where P_ij = δ_ij - k̂_i k̂_j. In the equation above, we have neglected the parity-odd contribution from the helical component of the PMF. Including the helical component requires introducing additional free parameters and is largely degenerate with the overall normalization, weakening the bounds on the PMF strength by about 25%<cit.>. Given this, for the sake of simplicity, we omit the helical contribution in the remaining discussion. However, future data can, in principle, measure it through parity-odd CMB polarization spectra <cit.>. 
        
        The magnetic power spectrum is commonly taken to be a power law parametrized as
        
        P_B(k)  = 
        {[ S_0   k^n_B,        k<k_D;           0,      k ≥ k_D ].
        
        where n_B is the spectral index and 2π/k_D is the damping scale below which magnetic fields dissipate due to radiation viscosity <cit.>.  The spectral index depends on the generation mechanism of the PMF. The originally proposed simple inflationary models of magnetogenesis predicted a scale-invariant spectrum, corresponding to n_B=-3, although other values are possible in more complicated models <cit.>. The PMFs produced in phase transitions generically have n_B=2 on large scales probed by CMB <cit.>. The damping scale depends on the amplitude of the PMF spectrum and the spectral index as <cit.>
            k_DMpc^-1 = [ 5.5 × 10^4 h ( B_λ/nG)^-2( 2 π/λ / Mpc)^n_B+3Ω_b h^2/0.022]^1/n_B+5,
        
        where Ω_b is the baryon density fraction, h is the reduced Hubble constant, H_0 = 100   h  km/(s Mpc), and B_λ is obtained by smoothing the magnetic energy density over a comoving wavelength λ using a Gaussian filter:
        
        B^2_λ =  1/(2 π)^3∫_0^∞ d^3 k   e^-k^2 λ^2 P_B(k)
        =  2 S_0/(2 π)^21/λ^n_B+3Γ( n_B+3/2).
        
        It is common to take λ = 1  Mpc, as it corresponds to the size of a typical region at the time of last scattering that later collapses to form a galactic halo. Another way to quantify the strength of the PMF is in terms of its contribution to the radiation energy density, namely,
        
            ℰ_B = 1/(2 π)^3∫_0^k_D dk   k^2 P_B(k),
        
        and defining the effective PMF strength as B_eff = √(8 πℰ_B). Correspondingly, one can define the fractional PMF density as  Ω_Bγ = ℰ_B / ρ_γ where ρ_γ is the total radiation energy density.  The relation between B_eff, B_λ and Ω_Bγ is given by <cit.>
            B_eff = B_λ (k_D λ)^n_B+3/2/√(Γ ((n_B+5)/2)) = 3.3 × 10^3√(Ω_Bγ)  nG.
        
        The stress-energy tensor associated with the PMF can be written as <cit.>
        T^   0_B  0(𝐱, τ) = - 1  8 π a^4 B^2(𝐱) ≡ -ρ_γΔ_B;   
        T^  i_B  j(𝐱, τ) = 1/4 π a^4( 1/2 B^2(𝐱) δ^i_j - B^i(𝐱)B_j(𝐱) )
        ≡ p_γ(Δ_B δ^i_j + Π^   i_B   j),
        
        where B_i(𝐱) = B_i (𝐱, τ_0), ρ_γ and p_γ=ρ_γ/3 are the photon density and pressure, Δ_B is the magnetic contribution to the radiation density contrast and Π^   i_B   j is the dimensionless anisotropic stress. 
        
        CMB is studied through measurements of four Stokes parameters: the intensity I, functions Q and U describing linear polarization, and V which quantifies circular polarization. CMB is not expected to be circularly polarized, and none of the processes discussed in this paper can generate V at an appreciable level, hence we will ignore V from now on. A typical CMB experiment would measure I, Q and U at different points on the sky and several frequencies. The primary sources of anisotropies preserve the black body spectrum of the primordial CMB signal, while altering the temperature of the black body along different directions on the sky. For this reason, I, Q and U are quantified in terms of fluctuations of the black body temperature in a given line of sight . 
        
        Rather than interpreting the data in terms of Q and U, which depend on the choice of the coordinate basis, it is beneficial to project them onto a parity-even and a parity-odd basis polarization patterns, referred to as E (gradient) and B (curl) modes <cit.>. Namely, one can write 
        
        (Q ± iU)() = ∑_lm (E_lm± i B_lm) _± 2Y_lm()
        
        where _± 2Y_lm() are the spin-2 spherical harmonics, and E_lm and B_lm are the E and B mode multipoles. The E and B mode spectra are obtained by taking ensemble averages of products of the multipole coefficients. For instance, the B-mode spectrum, C_l^BB, is defined as ⟨ B^*_lm B_l'm'⟩ = δ_ll'δ_mm' C_l^BB. Separating the polarization signal into E and B modes helps to reveal the physical processes responsible for polarizing CMB. Intensity gradients present at last scattering generate parity-even E-mode polarization patterns, while B-modes require sources that have handedness, such as gravitational waves or vorticity in the metric. 
        
        The correspondence between the observed E and B modes and their primordial origin is partially obscured by the galactic foregrounds, primarily the polarized dust and the synchrotron emission, which generate both modes at roughly the same magnitude. They pose a particularly strong challenge for B-modes, since the target primordial signal is extremely faint in comparison. For example, the B-mode spectrum sourced by gravitational waves generated during inflation with tensor-to-scalar ratio r=0.001, which is the target of the next generation of CMB experiments, is at least two orders of magnitude weaker than the foreground contribution. Since foregrounds typically depend on frequency, they can be subtracted, in principle, by comparing CMB maps at different frequencies. However, doing it with the above-mentioned accuracy is an extremely challenging ongoing problem. We are specifically interested in B-modes sourced by the PMF.  As we will discuss below, most of the information containing the PMF signatures comes from multipoles in the 500<l<2000 range, where the EE spectrum peaks. This happens to correspond to the range of scales on which the foregrounds contribute the least. Also, it is sufficient to work with only 10-20% of the sky around the galactic poles, where the galactic contamination is smaller.
        
        The PMF stress-energy sources perturbations in the metric, which include vector (vorticity) and tensor (gravitational waves) modes. As we will discuss below, B-modes are primarily sourced by the “compensated” vector modes and the “passive” tensor modes. In addition, the polarization of CMB is rotated along the line of sight by FR, which is frequency-dependent and introduces a characteristic mode-coupling signature. We briefly review these effects below.
        
        
        
         §.§ Anisotropies from metric fluctuations sourced by a stochastic PMF
        
        
        Let us consider a stochastic PMF generated at some early epoch before last scattering. The PMF induced CMB effects are constrained to be small on cosmological scales, hence we can work at linear order in perturbation theory where the PMF is treated as a stiff source with no back-reaction from the metric or matter fluctuations. In Fourier space, the traceless symmetric tensor Π_B  j^ i in Eq. (<ref>) can be decomposed into scalar (S), vector (V) and tensor (T) components <cit.> as
        
            Π_Bj^ i (𝐤,τ) e^i 𝐤·𝐱= Π_Bj^  i(S) + Π_Bj^  i(V) + Π_Bj^  i(T),
        
        with the components given by
        Π_B  ij^(S) = Π_B^(0)Q_ij^(0);   Π_B  ij^(V) = Π_B^(+1)Q_ij^(+1) +  Π_B^(-1)Q_ij^(-1) ;  Π_B  ij^(T) = Π_B^(+2)Q_ij^(+2) +  Π_B^(-2)Q_ij^(-2), 
        
        where Q_ij^(0) = -(k̂_i k̂_j - 1/3 δ_ij) exp(i 𝐤·𝐱), Q_ij^(± 1) = i k̂_(   i e^±_j   )exp(i  𝐤·𝐱) and Q_ij^(± 2) = e_ij^(± 2)exp(i  𝐤·𝐱) are, respectively, the scalar, vector and tensor harmonic functions <cit.>. Here, 𝐞^(±) = -i/√(2) (𝐞^1 ± i 𝐞^2), 𝐞^1,2 are the unit vectors orthogonal to the wave vector 𝐤 and e_ij^(± 2) = √(3/2) e_i^(±) e_j^(±).
        
        At linear order in perturbation theory, the scalar, vector and tensor modes evolve independently. There is a subtlety in setting their initial conditions that stems from the presence of neutrinos. Namely, after neutrinos decouple from photons, which happens at energies below 1 MeV, they free stream and can develop a non-zero anisotropic stress that compensates the stiff anisotropic stress of the PMF. However, prior to their decoupling, neutrinos are bound in a tightly coupled fluid with photons and baryons and are unable to compensate for the magnetic anisotropic stress, which then acts as a source of adiabatic scalar and tensor mode perturbations <cit.>. On scales that remain outside the horizon, the amplitudes of these perturbations grow logarithmically with time until neutrinos decouple and remain constant afterwards. In practice, the initial conditions for perturbations are set on super-horizon scales at some time after the neutrino decoupling. Then, perturbations sourced prior to decoupling set the initial amplitude of the adiabatic mode, which subsequently evolves passively without further influence from the PMF stress-energy. Further, they are assumed to be uncorrelated with the adiabatic fluctuations generated by inflation and are treated as separate passive magnetic modes <cit.>. After the neutrino decoupling, the PMF stress-energy continues to actively source  compensated perturbations of isocurvature type  <cit.>, which evolve independently from the adiabatic mode. The adiabatic vector modes decay rapidly when not sourced and can be ignored, thus CMB spectra receive contributions from five magnetic modes in total: a passive and a compensated scalar mode, a passive and a compensated tensor mode, and a compensated vector mode. At linear order, contributions of the different modes to the CMB spectra are uncorrelated and can be computed separately. We will focus only on the two modes that can contribute to the B-mode spectrum at an observationally relevant level: the compensated vector mode and the passive tensor mode[A detailed discussion of the different modes and their CMB signatures can be found in <cit.>]. A publicly available code[The code can be downloaded at <https://alexzucca90.github.io/MagCAMB/>.] for computing PMF sourced CMB spectra was developed in <cit.>  based on CAMB <cit.> and an earlier code written by R. Shaw <cit.>.
        
        The PMF acts as an active source of vector perturbations with the vector component of the magnetic anisotropic stress, Π_B^(± 1), being compensated by the anisotropic stress of neutrinos, Π_ν^(± 1). The compensated tensor mode has a negligibly small amplitude and can be ignored. However, by the time of neutrino decoupling, the PMF sources an adiabatic tensor mode of amplitude <cit.>
            H^± 2≈ R_γΠ_B^(± 2)[ ln(τ_ν / τ_B) + (5/8 R_ν -1 ) ],
        
        where τ_ν and τ_B are the conformal times at the epochs of neutrino decoupling and the PMF generation, and R_γ and R_ν are the photon and neutrino fractions of the total radiation density. 
        
        
        
        In Fig. <ref> (left and centre panels), we show the compensated vector and the passive tensor  modes contributions to the B-mode spectrum (BB) for B_1  Mpc = 2 nG, β=log_10 (τ_ν/τ_B) = 17 and different values of the spectral index. Note the change in the dependence on n_B that occurs at n_B = -1.5. In the -3<n_B<-1.5 range, the shapes of the CMB spectra are determined by a Fourier spectrum that scales as k^2 n_B + 6<cit.>. Thus, in that regime, an increase in n_B causes a shift of power from lower to higher l. This reduces the CMB anisotropy power on scales within the observational window. In contrast, for n_B ≥ -1.5, the CMB spectra are independent on n_B because the Fourier spectrum is cutoff dominated. In that case, the amplitude of the spectrum increases with n_B, thus larger values of n_B lead to more CMB power for the same PMF strength on 1 Mpc scale. 
        
        The vector mode contribution to BB, being independent of τ_B, offers a less model-dependent constraint on the PMF strength. Furthermore, it contributes to scales on which foregrounds are less of a problem and where the lensing contribution to the B-mode spectrum has already been measured. In the future, when measurements at lower l become available, simultaneous constraints on the tensor and vector modes will allow to constrain τ_B, which is of fundamental interest. 
        
        The contribution of the magnetic vector mode to the B-mode polarization power spectrum can be well constrained by the current and future CMB experiments capable of detecting the B-modes from weak lensing. It is clear from the right panel of Fig. <ref> that the high-l measurements of the B-mode polarization performed by SPT can place competitive bounds on the amplitude of the PMF. As shown in <cit.>, they significantly reduce the upper bound obtained from Planck.
        
        
        
         §.§ Mode-coupling correlations induced by Faraday rotation
        
        
        The polarization vector of CMB photons passing through magnetized plasma is rotated by an angle <cit.>
            α(n̂) = 3c^2 ν^-2/16 π^2 e∫τ̇  B· d l 
            =c^2 ν^-2  RM(n̂),
        
        where τ̇ is the differential optical depth, 
        ν is the frequency at which CMB is observed, B is the comoving magnetic field, and d l is the comoving length element along the photon trajectory. The rotation measure (RM) is a commonly used frequency independent characteristic of the strength of FR. Eq. (<ref>) implies that a significant RM can be produced by a small PMF over a large distance, which is the case at recombination, or by a larger magnetic field over a smaller path, which is the case inside our galaxy. The latter acts as noise in the estimates of the RM due to a PMF and was studied in detail in <cit.>. The main contribution of the PMF to FR happened just after the polarization was generated at the time of last scattering, when the mean ionized fraction was still high. Subsequently, additional FR was produced in ionized regions along the line of sight that contain magnetic fields, such as clusters of galaxies and our own galaxy. We will not discuss FR from clusters, as they are expected to induce rotations of CMB polarization on relatively small scales. On the other hand, the Milky Way contribution to the RM near galactic poles can look very similar to that of a scale-invariant PMF <cit.>.
        
        When the polarization is rotated, the two Stokes parameters transform as
        
            Q(ν) + iU(ν) =(Q^(0) + iU^(0)) exp(2iα(ν))  ,
        
        where Q^(0) and U^(0) are the Stokes parameters at last scattering. Since the FR angle falls off as 1/ν^2, Q^(0) and U^(0) can also be taken to be the Stokes parameters measured at a very high frequency. In theory, it is possible to extract the FR angle by taking maps of Q and U at different frequencies and using Eq. (<ref>) 
        to solve for the rotation in each pixel. Each additional frequency channel provides a separate measurement of α(ν) thus reducing the error bar on the measurement of RM(). Such a direct measurement of FR would be interesting to attempt, but could be challenging given that Q and U maps require significant processing to subtract the foreground contributions.
        
        Another way to extract the rotation angle is based on quadratic estimators utilizing mode-coupling correlations between E and B modes induced by FR <cit.>. It is similar to the procedure used to isolate the weak lensing contribution to CMB anisotropies <cit.>. Unlike a direct extraction of FR from Eq. (<ref>), the quadratic estimator method does not utilize the frequency dependence. It is statistical in nature and formally involves summing over all pixels of Q and U in order to reconstruct α in a given direction on the sky. 
        
        In the small rotation angle approximation, which is well-justified since the PMF is constrained to be very small, the relation between the spherical expansion coefficients of the E, B and α fields can be written as
        
        B_lm=2∑_LM∑_l' m'α_LM E_l' m'ξ_lml'm'^LMH_ll'^L  ,
        
        where ξ_lml'm_2^LM and H_ll'^L are related to Wigner 3-j symbols via <cit.>ξ_lml'm'^LM   ≡    (-1)^m √( (2l+1)(2L+1)(2l'+1)  4π)(
        [  l  L l'; -m  M m' ]) 
        
        H_ll'^L    ≡   (
        [  l  L l';  2  0 -2 ])  ,
        
        with the summation restricted to even L+l'+l. Eq. (<ref>) assumes “instant last scattering”, namely, that all of the E-mode polarization was generated instantly at some early time, and subsequently rotated. In reality, during last scattering, the polarization is rotated while it is being produced. This approximation has been investigated in great detail in <cit.>, where the exact radiation transfer equations in the presence of magnetic fields were solved. It was shown that corrections introduced by taking into account the width of the last scattering are not negligible, but not significant enough to warrant complicating the discussion of the observability of the FR effects.
        
        Eq. (<ref>) implies correlations between multipoles of E and B modes. Since the primordial T and E are correlated, FR also correlates T and B. The quadratic estimator method assumes the primordial E mode to be a statistically isotropic Gaussian random field characterized by a known power spectrum, so that ⟨ E_l m^* E_l'm'⟩ =  δ_ll'δ_mm' C_l^EE. Then, given a CMB polarization map, one can calculate a quantity
        
        [α̂_BE,LM]_ll' = 2π∑_mm' B_lmE^*_l'm'ξ_lml'm'^LM (2l+1)(2l'+1)C_l^EEH_ll'^L
        which is an unbiased estimator of α_LM<cit.> from a given pair of multipoles ll'. 
        The optimal estimate of α_LM is obtained by summing over all available ll' pairs weighted appropriately to minimize the variance. An analogous estimator can be constructed from products of T and B. The variance is further reduced when combinations of E, B and T from different frequency channels are combined. In such a case, one works with the estimator of RM defined, in the case of EB, as <cit.>
        [r̂_B^iE^j,LM]_ll' = c^-2ν^2_i [α̂_B^iE^j,LM]_ll' ,
        
        where i and j label different frequency channels. 
        
        The variance in r̂_LM, for a statistically isotropic RM, is defined as
        ⟨r̂^*_LMr̂_L'M'⟩ =δ_LL'δ_MM' [C_L^ RM+σ^2_ RM,L],
        where C_L^ RM is the RM power spectrum that receives contributions from the PMF and the galaxy, while σ^2_ RM,L is the combined variance of individual estimators [r̂_B^iE^j,LM]_ll'. Using notation similar to that of <cit.>, the variance can be written as
        σ^-2_ RM,L = ∑_l, l' ≥ l G_l l'∑_A,A' [( C^l l')^-1]_AA' 
        Z_l l'^A Z_l l'^A' ,
        
        where the sum is restricted to even l+l'+L, G_l l'≡ (2l+1)(2l'+1) (H^L_l l')^2 / π, A and A' label all relevant quadratic combinations of E, B and T, 
        
        Z_l l'^X^iB^j = c^2ν_j^-2 W^ij_l l' C^XE_l  , 
        
        Z_l l'^B^iX^j = c^2ν_i^-2 W^ij_l l' C^EX_l' ,
        
        with X denoting either T or E, and W^ij_l l'≡exp[-(l^2+l'^2) θ^2_ij/16ln 2] accounts for the finite width of the beam. We take θ_ij = max[θ^i_ fwhm,θ^j_ fwhm], where θ^i_ fwhm is the full-width-at-half maximum (FWHM) of the Gaussian beam of the i-th channel. The covariance matrix elements, [ C^l l']_AA', are
        
        [ C^l l']_X^iB^j,Y^kB^n= C̃^X^iY^k_l C̃^B^jB^n_l' ;  
        [ C^l l']_B^iX^j,B^kY^n= C̃^B^iB^k_l C̃^X^jY^n_l'
        with X and Y standing for either E or T, and
        C̃^X^iY^j_l ≡ C^XY,  prim_l+f_ LC^XY,  WL_l+ δ_X^iY^jσ^2_P,i ,
        
        is the measured spectrum, that includes the primordial contribution C^XY,  prim_l, the WL contribution C^XY,  WL_l, and the detector noise σ^2_P,i which is taken to be uncorrelated between different maps. The de-lensing fraction f_ L is introduced to account for the partial subtraction of the WL contribution. According to <cit.>, the quadratic estimator method of <cit.> can reduce the WL contribution to C̃_l^BB by a factor of 7 (implying f_ L =0.14), with iterative methods promising a further reduction <cit.>.
        
        The signal-to-noise ratio of the detection of the primordial RM spectrum C_L^ RM,PMF is given by
        ( S  N )^2 = ∑_L=1^L_max(f_ sky/2) (2L+1) [C_L^ RM,PMF]^2  [C_L^ RM,PMF + f_ GC_L^ RM,G+σ^2_ RM,L]^2 ,
        
        where f_ G is the fraction of the galactic RM spectrum that contributes to the total measured signal. Estimates of C_L^ RM,G for different sky cuts were obtained in <cit.> based on the galactic RM map of Oppermann et al <cit.>. The RM map <cit.> was produced using a compilation of FR of extragalactic radio sources. They found that the RM is the strongest along the galactic plane, falling off as one moves towards the poles and becoming latitude independent at latitudes higher than 80^∘. Near the galactic poles, the RM is ∼ 20 rad/m^2 and has a nearly scale-invariant anisotropy spectrum for L>200. It was found that, because of the strong RM away from the poles, it is optimal to use f_ sky∼ 0.2 for constraining the PMF.
        
        The quadratic estimator works for the monopole (L=0) and the dipole (L=1) of the FR field which, in principle, should not be ignored. However, the FR angle is not expected to have a monopole, since it would imply a non-zero magnetic charge enclosed by the CMB surface, while the dipole, corresponding to a uniform magnetic field, is constrained to be negligible for the PMF and is relatively unimportant for the galactic RM.
        
        The variance of the rotation measure estimator σ_ RM,L is the smallest for large scale rotations and increases steadily with increasing L. This makes physical sense – a rotation that is uniform over a large part of the sky will couple many ll' pairs in exactly the same way, boosting the signal-to-noise of its detection. Thus, while one wants to measure E and B modes at l ∼ 1000, one is actually probing rotations at L < 100<cit.>. For this reason, the mode-coupling estimator is generally more sensitive to PMFs that have a large fraction of their power on large scales. 
        
        For a scale-invariant PMF, the dependence of the signal-to-noise on the survey parameters can be approximated by <cit.>S  N≈( B_λ 0.1   nG) ( f_ sky 0.2 )^1/2(28' - θ_ fwhm 10' ) (95  GHzν)^2 (1 μK-arcminσ_P)  ,
        
        which works reasonably well for experiments capable of detecting sub-nG level fields. One can deduce the expected bound on the PMF strength from (<ref>) by setting the S/N according to the desired significance level,  e.g.S/N=1 for the 1σ bound,  or S/N=2 for 2σ,  etc, and solving for B_λ. The reader should be cautioned that this expression does not account for the weak lensing contribution and the galactic RM and hence should be taken only as a very approximate “in principle” limit achievable with a given experiment. In the forecasts presented in the next Section, all bounds are calculated numerically using the exact expression (<ref>).
        
        
        
        In addition to the mode-coupling correlations of EB and TB type, FR contributes to the B-mode polarization spectrum, C_l^BB<cit.>. A detailed discussion of the FR induced B-mode spectrum can be found in <cit.>. As shown in <cit.>, for scale-invariant PMF, the FR signatures are always less visible in the B-mode spectrum compared to the mode-coupling quadratic estimator. This is due to two primary factors: a) the mode-coupling estimator is linear in the small rotation angle α, while BB scales as α^2, and b) in the mode-coupling estimator (<ref>), many ll' multipole pairs of E and B fields probe the same L mode of the rotation field. The mode-coupling estimator is primarily sensitive to large scale rotations and, for that reason, it is more sensitive to near-scale-invariant PMFs and less so to causally generated PMF, which have blue spectra with n_B=2. In the latter case, most of the FR rotation happens on very small scales and the S/N of detecting FR in BB is actually higher than that in the quadratic estimator <cit.>. 
        
        Fig. <ref> shows the FR induced BB spectra at 70 GHz for the scale-invariant field of 1 nG strength and a causal field with B_ eff=10 nG. Also shown is the BB spectrum from the galactic FR obtained in <cit.> using the rotation measure map of <cit.> after applying the Planck galactic cut. As one can see, the FR sourced BB is always well below the weak lensing signal. Note that the amplitude of the FR sourced B-mode spectrum scales as ν^-4, and the relatively foreground-free CMB frequencies are in the 90-150 GHz range. Thus, overall, the BB sourced by primordial or the galactic FR is not expected to be significant and we will not discuss it further.
        
        
        
        § CURRENT AND FUTURE B-MODE BOUNDS ON PRIMORDIAL MAGNETIC FIELDS
        
        
        CMB polarization was first detected by DASI in 2002 <cit.>, followed shortly by WMAP <cit.>. In the subsequent decade, the E-mode polarization measurements steadily improved, culminating in the high-resolution EE spectra obtained by Planck <cit.>. Until recently, only upper bounds on the B-mode spectrum were available. However, thanks to the pioneering work by the POLARBEAR <cit.>, BICEP/Keck <cit.> and SPT <cit.> collaborations, cosmology has entered the era of precision B-mode science. The ongoing and future experiments will improve the accuracy and resolution of B-mode measurements by orders of magnitude. In what follows, we will review the current PMF bounds derived from B-modes and present forecasts for the upcoming and future CMB experiments. We will consider both ways in which constraints on the PMF can be derived from CMB polarization, namely, the B-mode spectrum and the mode-coupling correlations induced by FR.
        
        
        
         §.§ Current bounds from the B-mode spectrum and Faraday rotation
        
        
        Prior to the B-mode data becoming available, CMB bounds on the PMF were derived from the available measurements of TT, TE and EE <cit.>. These spectra are dominated by the primary CMB anisotropies, with the PMF contribution being relatively small and partially degenerate with changes in cosmological parameters. The availability of the B-mode spectra significantly improved the bounds. This was first demonstrated by the POLARBEAR collaboration, who placed a bound of B_1 Mpc < 3.9 nG <cit.>, comparable to the one from Planck, using only their B-mode spectrum. One should note that Planck data was used indirectly to determine the cosmological parameters assumed in <cit.>, however the PMF contributions to TT, TE and EE were not used in deriving the constraint. 
        
        A comprehensive analysis was performed by Zucca et al <cit.>, who considered the PMF contributions to all CMB spectra and obtained constraints based on the Planck TT, TE and EE <cit.> in combination with the B-mode spectrum measured by SPT <cit.>. Fig. <ref>, based on the analysis in <cit.>, allows one to compare the PMF bounds obtained from the Planck data alone, and after adding the BB measurements. The PMF bounds obtained by the Planck collaboration in <cit.> are broadly consistent with <cit.>, with minor differences that have been addressed in <cit.>.
        
        To obtain the distributions in Fig. <ref>, in addition to the standard cosmological parameters, three PMF parameters were varied: the amplitude B_1 Mpc, the spectral index n_B, and the PMF generation time parameter β = log_10(τ_ν / τ_B). The left and the centre panels of Fig. <ref> show the joint probability obtained in <cit.> for the magnetic amplitude B_1  Mpc and the magnetic index n_B from Planck alone and after combining Planck with SPT. One can see that the bound on B_1  Mpc becomes weaker with increasing n_B in the -3<n_B<-1.5 range. However, the bound becomes stronger when n_B>-1.5. This is due the qualitative change in the dependence of the CMB spectra on n_B that occurs at n_B=-1.5.  Namely, as discussed in <cit.> and can be seen from Fig. <ref>, for -3<n_B<-1.5, an increase in n_B results in a shift of power from lower to higher l, reducing the signal inside the observational window, and thus allowing for larger PMF strengths. In contrast, when n_B>-1.5, the CMB spectra are cutoff dominated and no longer change in shape, with larger n_B giving larger amplitudes for the same PMF strength, leading to tighter constraints on B_1  Mpc.
        
        Since there was no detection of the PMF amplitude, the magnetic spectral index n_B and the PMF generation epoch parameter β are unconstrained. To see this, one must choose appropriate priors on B_1  Mpc. Fig. <ref> separately shows the cases with a uniform (left panel) and the logarithmic (center panel) prior on B_1  Mpc. One can see an apparent upper bound on n_B in the case of the uniform prior, which is clearly unphysical and is due to the fact that extremely small values of B_1  Mpc have effectively zero prior probability. The upper bound on n_B disappears if one uses a logarithmic prior on B_1  Mpc which treats values of different orders of magnitude with equal probability. Incidentally, the choice of the prior does not affect the upper bound on B_1  Mpc.
        
        There is a potential degeneracy between the tensor-to-scalar ratio r and the time of the generation of the PMF parameter β, since, for nearly scale-invariant PMFs, the passive tensor magnetic mode is similar in shape to the inflationary tensor mode. To address this, Zucca et al <cit.> separately considered the cases with a fixed r=0, and when the two parameters are co-varied. It was found that this degeneracy does not affect the current bounds since the B-mode data on the relevant scales is still very weak. This will change in the future, when data from experiments such as the Simons Observatory <cit.>, designed to detect r ∼ 0.001, becomes available. Having a reliable B-mode spectrum measured over the 50 < l < 3000 range would simultaneously constrain the amplitude and the generation epoch of the PMF, helping to narrow the search for the correct theory of magnetogenesis.
        
        As Fig. <ref> demonstrates, adding the B-mode data significantly tightens the constraints. In the case when the magnetic spectral index n_B is allowed to vary, the 95% CL bound on B_1  Mpc is reduced from 3.3 nG to 1.55 nG. The bounds are stronger when the spectral index is fixed to n_B ≈ -3, which is of particular interest because the simplest models of inflationary magnetogenesis <cit.> predict a nearly scale-invariant PMF. In the case of n_B=-2.9, the 95% CL bound is reduced from 2 nG for Planck to 1.2 nG for Planck+SPT. We note that the dependence on the smoothing scale disappears in the scale-invariant case, and B_1  Mpc = B_ eff.
        
        The other special cases is that of PMFs generated after inflation,  e.g. during preheating <cit.> or in phase transitions <cit.>. Such fields have very small coherence lengths and are uncorrelated on cosmological scales. This leads to a generic prediction, based on causality, for such fields to have n_B = 2 on scales relevant to CMB anisotropies <cit.>. Since these fields have most of their power concentrated on small scales, 1 Mpc is not a representative scale for assessing their strength. It is more appropriate to use Ω_Bγ or B_ eff, which quantify the total energy density in the PMF <cit.>. The Planck+SPT provides a constraint of Ω_Bγ < 10^-3 or B_ eff < 100 nG at 95% CL, which is two orders of magnitude stronger than the BBN bound of Ω_Bγ < 0.1<cit.>. The corresponding bound on B_1  Mpc happens to be B_1  Mpc< 0.002 nG at 95% CL and represents the fact that most of the PMF power is concentrated on much smaller scales.
        
        Future BB measurements will not seriously change the current  constraint on the PMF because of the quartic scaling of CMB spectra with the field strength. Instead, the mode-coupling induced by FR scales linearly and should offer a more sensitive probe. The method has already been applied to the POLARBEAR maps <cit.> to derive a 95% CL bound of 93 nG on the scale-invariant PMF[Note that the POLARBEAR forecast performed in <cit.> assumed a frequency of 90 GHz and f_ sky=0.024, while <cit.> used the actual 148 GHz channel data with f_ sky=0.0006. This accounts for the significant difference between the 2σ bound of 3 nG forecasted in <cit.> and the 93 nG bound found in <cit.>.]. As we will see in the next subsection, the upcoming experiments will improve this bound by up to three orders of magnitude. Such a dramatic improvement will be due to a combination of factors, such as utilizing 100 GHz channels in addition to just the 150 GHz map used in <cit.>, having a larger sky fraction, f_ sky∼ 0.2 vs 0.0006 used in <cit.>, and having lower noise levels, σ_P ∼ 2μK-arcmin vs the 6.9μK-arcmin in <cit.>. Thus, even though today the mode-coupling method is not competitive compared to using the B-mode spectrum, it holds a huge promise for the upcoming and future experiments.
        
        
        
         §.§ Forecasted bounds from SPT-3G, the Simons Observatory, CMB-S4 and Space Probe
        
        
        The SPT B-mode spectrum <cit.>, shown in Fig. <ref>, constrains the magnetic vector mode, which contributes roughly on the same angular scale as the weak lensing B-mode. Over the next few years, Stage III CMB experiments, such as SPT-3G <cit.> and the Simons Observatory (SO) <cit.>, will measure the B-mode spectrum with a significantly better accuracy. However, while they will reduce the uncertainty in the PMF contribution to BB by more than an order of magnitude <cit.>, it will translate into only a minor improvement in the bound on the PMF strength thanks to the quartic dependence of the spectrum on B_λ. It will remain at or just below 1 nG even after the data from a future Stage IV ground based experiment <cit.> or a future space mission become available. 
        
        The FR induced mode-coupling, on the other hand, scales linearly with B_λ, and the qualitative improvements in the quality of CMB polarization maps will translate into dramatically improved bounds on nearly scale-invariant PMFs. As we show below, the 1σ uncertainty in B_λ will decrease from about ∼ 50 nG today <cit.> to ∼ 0.5 nG for Stage III and to ∼ 0.1 nG for Stage IV experiments. Crossing the 1 nG threshold is significant for ruling out a purely primordial origin (requiring no dynamo action) of the galactic magnetic field, and the FR induced mode-coupling will achieve this milestone.
        
        In what follows, we will present forecasts for three representative classes of experiments:
        
          *  Stage III ground based experiments, such as the planned SPT-3G <cit.> and the Simons Observatory (SO) <cit.>. Both are expected to measure the Stokes parameters at three frequencies: 95, 150 and 220 GHz. SPT-3G will offer ∼ 1 arcmin resolution covering 6% of the sky and will achieve noise levels below the weak lensing B-mode signal (< 5  μK-arcmin). The SO will cover a larger fraction of the sky with a resolution of a few arcmin and noise levels comparable to the weak lensing B-mode.
        
          *  A proposed Stage IV ground based experiment (CMB-S4) <cit.>, that will cover 40% of the sky at 100 and 150 GHz with ∼ 1 arcmin resolution and noise levels of ∼ 1μK-arcmin.
        
          *  A space mission, similar to the proposed CMBPol EPIC-IM <cit.> and CORE <cit.>, mapping the full sky using multiple channels covering a wide range of frequencies at a resolution of a few arcmin and noise levels comparable to the weak lensing B-mode. We will consider the channels in the 45-150 GHz range, where the mode-coupling induced by FR can be significant.
        
        The parameters assumed for each of the experiments are summarized in Table <ref>. Since the mode-coupling estimator is generally more sensitive to PMFs that have a large fraction of their power on large scales, we only consider the case of a scale-invariant PMF, with n_B =-2.9, with our conclusions being valid for a range of n_B close to that value.
        
        The forecasted bounds on the scale-invariant PMF amplitude are shown in Table <ref> for different choices of parameters f_L and f_G quantifying the fraction of the weak lensing contribution to BB and the fraction of the galactic rotation measure, respectively. The case of no de-lensing corresponds to f_L=1, while a perfect de-lensing corresponds to f_L=0. Similarly, f_G=1 and 0 correspond to no galaxy subtraction and a full subtraction, respectively. While a full subtraction of either weak lensing or the galactic rotation measure is impossible in practice, we opt to include these possibilities to demonstrate the ultimate sensitivity of the quadratic estimator to the PMF.
        
        We find that the two Stage III experiments will have comparable sensitivities to the PMF, promising to achieve ∼ 0.6-0.7 nG for SPT-3G and ∼ 0.4 nG for the SO at 68% CL. While SPT-3G will offer a higher resolution and lower noise levels that the SO, it will cover a smaller fraction of the sky. Because the SPT-3G detector noise is expected to be securely below the weak lensing BB, the SPT-3G bounds are more sensitive to the degree of de-lensing. On the other hand, the expected noise levels of the SO are comparable to the weak lensing contribution and subtracting it does not make as much of a difference. Instead, the stronger bounds on the PMF from the SO are due to a larger sample of ll' pairs available with the larger sky coverage. As expected, the superior accuracy with which Stage III experiments will measure the B-mode spectrum compared to the existing data from SPT will not lead to a notable improvement in the BB-based bound on the PMF strength. Instead, the mode-coupling based bound will be reduced by 2 orders of magnitude, from ∼ 50 to ∼ 0.5 nG.
        
        The SO is considered to be a path-finder experiment towards a future Stage IV experiment. Both will cover about half the sky, but CMB-S4 will do it at a much higher resolution and lower noise levels. It can constrain the PMF to below 0.2 nG without any de-lensing and can, in principle, achieve bounds below 0.1 nG. The bounds would get even tighter if channels at frequencies lower than 95 GHz were added. However, at some point, the galactic RM, which is roughly equivalent to the RM produced by a ∼ 0.1 nG scale-invariant PMF <cit.>, will pose an obstacle to pushing the bound lower.
        
        
        
        As one can see from Table <ref>, a realistic space mission would not necessarily offer stronger bounds than CMB-S4. A space based experiment would generally be better for studying large scale CMB anisotropies and for measurements at higher frequencies needed to better understand the galactic foregrounds. Since the mode-coupling estimator relies on measurements at l∼ 1000, it does not benefit significantly from a larger sky fraction. Also, the galactic rotation measure is very strong near the galactic plane, restricting the useful fraction of the sky to f_ sky∼ 0.2 around the galactic poles. Nevertheless, we find that a space probe such as EPIC-IM or CORE would offer bounds of ∼ 0.1 nG, comparable to CMB-S4.
        
        Overall, our forecasts show that one can achieve bounds close to 0.1 nG without worrying about FR due to our galaxy. Going below 0.1 nG would only be possible <cit.> if there was an independent measurement of the galactic RM. A full sky RM map was constructed in Oppermann et al <cit.> based on a compilation of RMs of known radio sources. Their map was used in <cit.>, as well as in the present analysis, to estimate the obstruction to the PMF posed by the galaxy. While the map of <cit.> has large uncertainties, it is plausible to expect better maps, based on many more radio sources, to become available in the future. Thus, in principle, it should be possible to use the FR mode-coupling estimator to push the bound on the PMF below 0.1 nG.
        
        
        
        § OUTLOOK
        
        
        Cosmology has entered the era of precision B-mode science thanks to the breakthrough results achieved by POLARBEAR <cit.>, BICEP/Keck <cit.> and SPT <cit.>. The next generation of CMB polarization experiments will significantly improve the accuracy of B-mode measurements. In addition to searching for signatures of inflationary gravitational waves <cit.> and PMFs, they will probe neutrino masses <cit.>, modifications of gravity <cit.>, cosmic (super)strings <cit.> and other fundamental physics <cit.>.
         
        In this article, we have reviewed the most relevant B-mode signatures of a PMF. At present, and in the near term, the tightest bounds on the PMF strength will come from the vector and tensor mode magnetic contributions to the BB spectrum. When accurate measurements of BB at lower l become available, it will be possible to simultaneously constrain the amplitude and the generation epoch of the PMF, which is of great relevance to theories of the early universe. Also, having accurate estimates of parity-odd TB and EB spectra would allow to constrain the helical component of the PMF<cit.>, which has been neglected in the present work. 
        
        In the future, looking for mode-coupling correlations induced by FR will produce tighter bounds. For the mode-coupling estimator approach to become competitive, one needs an experiment with a sufficiently high resolution to collect as many modes as possible near the peak of the EE spectrum, or in the 500<l<3000 range. One also needs noise levels below the weak lensing B-mode, which is around 5μK-arcmin. Achieving a higher resolution and lower noise is, in principle, a matter of cost – it can be accomplished by building larger telescopes with more detectors. While still challenging, it is arguably a simpler task than subtracting the large galactic foregrounds that obscure the primordial B-modes at lower l. FR estimators would also be aided by having more channels at frequencies below 100 GHz, but they do not benefit from increasing the sky coverage beyond 20%, as the galactic FR becomes too strong away from the galactic poles. 
        
        B-modes offer a relatively clean probe of the PMF, as there are relatively few other sources contributing on the relevant scales. The most relevant competing contribution to the BB spectrum comes from weak lensing, which can be partially subtracted using de-lensing techniques <cit.>. Fortunately, the mode-coupling correlations induced by lensing are of opposite parity to those induced by FR, and do not interfere with reconstruction of the FR angle. The lensing does, however, contribute to the variance of the quadratic estimator of FR, which will be important for Stage III and Stage IV experiments with noise levels below 5μK-arcmin. The frequency dependence of the FR signal also makes it distinct from other potential sources of birefringence,  e.g. caused by pseudo-scalar fields <cit.>.
        
        Our analysis did not account for instrumental systematic effects, such as the beam asymmetry and imperfect scanning of the sky, which can play an important role. Such systematics generate non-zero parity-odd spectra of EB and TB type, in addition to contributing to all parity-even spectra, including BB <cit.>. The impact of instrumental systematics on the ability to constrain the PMF is currently under investigation <cit.>.
        
        The other important obstruction comes from the galactic RM, which can hinder the detection of FR from the PMF. To detect PMF strengths below 0.1 nG one would need a more accurate galactic RM map <cit.>. However, even with no de-lensing and no subtraction of the galactic RM, Stage III experiments will be able to decrease the upper bound on the PMF strength to ∼ 0.5 nG, while a Stage IV experiment or an EPIC-IM class space mission would reduce it to below 0.2 nG. This will conclusively rule out or confirm the existence of ∼ nG strength PMFs required to generate μG galactic fields with no help from a dynamo mechanism.
        
        
        We benefited from earlier collaborations with Soma De, Yun Li, Bess Ng, Tanmay Vachaspati and Amit Yadav, as well as communications and discussion with Camille Bonvin, Fabio Finelli, Tina Kahniashvili, Brian Keating, Antony Lewis, Daniela Paoletti, Richard Shaw and Meir Shimon. We acknowledge support by the Natural Sciences and Engineering Research Council of Canada (NSERC). This research was enabled in part by support provided by WestGrid <cit.> and Compute Canada <cit.>.
        
        99
        Widrow:2002ud 
          L. M. Widrow,
          
          Rev. Mod. Phys.  74, 775 (2002)
          
          [astro-ph/0207240].
          
          
        Grasso:2000wj 
          D. Grasso and H. R. Rubinstein,
          
          Phys. Rept.  348, 163 (2001)
          [astro-ph/0009061].
          
          
        
        Widrow:2011hs 
          L. M. Widrow, D. Ryu, D. R. G. Schleicher, K. Subramanian, C. G. Tsagas and R. A. Treumann,
          
          Space Sci. Rev.  166, 37 (2012)
          
          [arXiv:1109.4052 [astro-ph.CO]].
          
          Athreya:1998
        R. M. Athreya, V. K. Kapahi, P. J. McCarthy, & W. van Breugel, Astronomy & Astrophysics  329, 809 (1998)
        
        
        Ryu:2008hi 
          D. Ryu, H. Kang, J. Cho and S. Das,
          
        Science  320, 909 (2008)
        doi:10.1126/science.1154923
        [arXiv:0805.2466 [astro-ph]].
        
          
        Ryu:2011hu 
          D. Ryu, D. R. G. Schleicher, R. A. Treumann, C. G. Tsagas and L. M. Widrow,
          
        Space Sci. Rev.  166, 1 (2012)
        doi:10.1007/s11214-011-9839-z
        [arXiv:1109.4055 [astro-ph.CO]].
        
          
        Schleicher:2012fx 
          D. R. G. Schleicher  et al.,
          
        Astron. Nachr.  334, 531 (2013)
        doi:10.1002/asna.201211898
        [arXiv:1211.4356 [astro-ph.CO]].
        
          
        Latif:2012aq 
          M. A. Latif, D. R. G. Schleicher, W. Schmidt and J. Niemeyer,
          
        Mon. Not. Roy. Astron. Soc.  432, 668 (2013)
        doi:10.1093/mnras/stt503
        [arXiv:1212.1619 [astro-ph.CO]].
        
          
        Wagstaff:2013yna 
          J. M. Wagstaff, R. Banerjee, D. Schleicher and G. Sigl,
          
        Phys. Rev. D  89, no. 10, 103001 (2014)
        doi:10.1103/PhysRevD.89.103001
        [arXiv:1304.4723 [astro-ph.CO]].
        
          
          
        
        Beck:2013gca 
          A. M. Beck, K. Dolag, H. Lesch and P. P. Kronberg,
          
        Mon. Not. Roy. Astron. Soc.  435, 3575 (2013)
        doi:10.1093/mnras/stt1549
        [arXiv:1308.3440 [astro-ph.GA]].
        
          
        Seifried:2013uta 
          D. Seifried, R. Banerjee and D. Schleicher,
          
        Mon. Not. Roy. Astron. Soc.  440, no. 1, 24 (2014)
        doi:10.1093/mnras/stu294
        [arXiv:1311.4991 [astro-ph.GA]].
        
          
            
        
        Vachaspati:1991nm 
          T. Vachaspati,
          
          Phys. Lett. B  265, 258 (1991).
          
          
          
        Brandenburg:2014mwa 
          A. Brandenburg, T. Kahniashvili and A. G. Tevzadze,
          
        Phys. Rev. Lett.  114, no. 7, 075001 (2015)
        doi:10.1103/PhysRevLett.114.075001
        [arXiv:1404.2238 [astro-ph.CO]].
        
          
          
        
        Brandenburg:2016odr 
          A. Brandenburg and T. Kahniashvili,
          
        Phys. Rev. Lett.  118, no. 5, 055102 (2017)
        doi:10.1103/PhysRevLett.118.055102
        [arXiv:1607.01360 [physics.flu-dyn]].
        
          
        Reppin:2017uud 
          J. Reppin and R. Banerjee,
          
        Phys. Rev. E  96, no. 5, 053105 (2017)
        doi:10.1103/PhysRevE.96.053105
        [arXiv:1708.07717 [astro-ph.CO]].
        
          
          
        
        Brandenburg:2017neh 
          A. Brandenburg, T. Kahniashvili, S. Mandal, A. R. Pol, A. G. Tevzadze and T. Vachaspati,
          
        Phys. Rev. D  96, no. 12, 123528 (2017)
        doi:10.1103/PhysRevD.96.123528
        [arXiv:1711.03804 [astro-ph.CO]].
        
          
        Brandenburg:2018ptt 
          A. Brandenburg, R. Durrer, T. Kahniashvili, S. Mandal and W. W. Yin,
          
        arXiv:1804.01177 [astro-ph.CO].
        
        Kernan:1995bz 
          P. J. Kernan, G. D. Starkman and T. Vachaspati,
          
          Phys. Rev. D  54, 7207 (1996)
          
          [astro-ph/9509126].
          
          
        Grasso:1996kk 
          D. Grasso and H. R. Rubinstein,
          
          Phys. Lett. B  379, 73 (1996)
          
          [astro-ph/9602055].
          
          
        Cheng:1996yi 
          B. l. Cheng, A. V. Olinto, D. N. Schramm and J. W. Truran,
          
          Phys. Rev. D  54, 4714 (1996)
          
          [astro-ph/9606163].
          
          
        Yamazaki:2012jd 
          D. G. Yamazaki and M. Kusakabe,
          
          Phys. Rev. D  86, 123006 (2012)
          doi:10.1103/PhysRevD.86.123006
          [arXiv:1212.2968 [astro-ph.CO]].
          
          
          
          
        Turner:1987bw 
          M. S. Turner and L. M. Widrow,
          
          Phys. Rev. D  37, 2743 (1988).
          
          
          
        Ratra:1991bn 
          B. Ratra,
          
          Astrophys. J.  391, L1 (1992).
          
          
          
        DiazGil:2007dy 
          A. Diaz-Gil, J. Garcia-Bellido, M. Garcia Perez and A. Gonzalez-Arroyo,
          
          Phys. Rev. Lett.  100, 241301 (2008)
          doi:10.1103/PhysRevLett.100.241301
          [arXiv:0712.4263 [hep-ph]].
          
          
        DiazGil:2008tf 
          A. Diaz-Gil, J. Garcia-Bellido, M. Garcia Perez and A. Gonzalez-Arroyo,
          
          JHEP  0807, 043 (2008)
          doi:10.1088/1126-6708/2008/07/043
          [arXiv:0805.4159 [hep-ph]].
          
          
        Adshead:2016iae 
          P. Adshead, J. T. Giblin, T. R. Scully and E. I. Sfakianakis,
          
          JCAP  1610, 039 (2016)
          doi:10.1088/1475-7516/2016/10/039
          [arXiv:1606.08474 [astro-ph.CO]].
          
          
         
        
        Durrer:2013pga 
          R. Durrer and A. Neronov,
          
          Astron. Astrophys. Rev.  21, 62 (2013)
          
          [arXiv:1303.7121 [astro-ph.CO]].
          
          
        Barnaby:2012tk 
          N. Barnaby, R. Namba and M. Peloso,
          
          Phys. Rev. D  85, 123523 (2012)
          doi:10.1103/PhysRevD.85.123523
          [arXiv:1202.1469 [astro-ph.CO]].
          
          
        
        Long:2013tha 
          A. J. Long, E. Sabancilar and T. Vachaspati,
          
          JCAP  1402, 036 (2014)
          doi:10.1088/1475-7516/2014/02/036
          [arXiv:1309.2315 [astro-ph.CO]].
          
          
        Mack:2001gc 
          A. Mack, T. Kahniashvili and A. Kosowsky,
          
          Phys. Rev. D  65, 123004 (2002)
          
          [astro-ph/0105504].
          
          
        
        Lewis:2004ef 
          A. Lewis,
          
          Phys. Rev. D  70, 043011 (2004)
          
          [astro-ph/0406096].
          
          
          
          
        Finelli:2008xh 
          F. Finelli, F. Paci and D. Paoletti,
          
          Phys. Rev. D  78, 023510 (2008)
          
          [arXiv:0803.1246 [astro-ph]].
          
          
          
          
        Paoletti:2008ck 
          D. Paoletti, F. Finelli and F. Paci,
          
          Mon. Not. Roy. Astron. Soc.  396, 523 (2009)
          
          [arXiv:0811.0230 [astro-ph]].
          
          
          
          
          Shaw:2010
          	R. Shaw and A. Lewis,
          	Phys. Rev. D  81, 043517 (2010)
          	 [arXiv:0911.2714 [astro-ph]].
          	 
         
        Subramanian:1997gi 
          K. Subramanian and J. D. Barrow,
          
          Phys. Rev. D  58, 083502 (1998)
          
          [astro-ph/9712083].
          
          
        Subramanian:1998fn 
          K. Subramanian and J. D. Barrow,
          
          Phys. Rev. Lett.  81, 3575 (1998)
          
          [astro-ph/9803261].
          
          
          
          
        Seshadri:2000ky 
          T. R. Seshadri and K. Subramanian,
          
          Phys. Rev. Lett.  87, 101301 (2001)
          
          [astro-ph/0012056].
          
          
        Subramanian:2002nh 
          K. Subramanian and J. D. Barrow,
          
          Mon. Not. Roy. Astron. Soc.  335, L57 (2002)
          doi:10.1046/j.1365-8711.2002.05854.x
          [astro-ph/0205312].
          
          
        Subramanian:2003sh 
          K. Subramanian, T. R. Seshadri and J. D. Barrow,
          
          Mon. Not. Roy. Astron. Soc.  344, L31 (2003)
          doi:10.1046/j.1365-8711.2003.07003.x
          [astro-ph/0303014].
          
          
          	 
        
        Kosowsky:1996yc 
          A. Kosowsky and A. Loeb,
          
          Astrophys. J.  469, 1 (1996)
          
          [astro-ph/9601055].
          
          
        Kosowsky:2004zh 
          A. Kosowsky, T. Kahniashvili, G. Lavrelashvili and B. Ratra,
          
          Phys. Rev. D  71, 043006 (2005)
          
          [astro-ph/0409767].
          
          
          
        
        Pogosian:2011qv 
          L. Pogosian, A. P. S. Yadav, Y. F. Ng and T. Vachaspati,
          
          Phys. Rev. D  84, 043530 (2011)
          Erratum: [Phys. Rev. D  84, 089903 (2011)]
        
          [arXiv:1106.1438 [astro-ph.CO]].
          
          
          
        
        
        Keisler:2015hfa 
          R. Keisler  et al. [SPT Collaboration],
          
          Astrophys. J.  807, no. 2, 151 (2015)
          
          [arXiv:1503.02315 [astro-ph.CO]].
          
          
        Ade:2015cva 
          P. A. R. Ade  et al. [Planck Collaboration],
          
          arXiv:1502.01594 [astro-ph.CO].
          
          
        Adam:2015rua 
          R. Adam  et al. [Planck Collaboration],
          
          arXiv:1502.01582 [astro-ph.CO].
          
          
        Paoletti:2012bb 
          D. Paoletti and F. Finelli,
          
          Phys. Lett. B  726, 45 (2013)
          doi:10.1016/j.physletb.2013.08.065
          [arXiv:1208.2625 [astro-ph.CO]].
          
          
        Larson:2010gs 
          D. Larson  et al.,
          
          Astrophys. J. Suppl.  192, 16 (2011)
          doi:10.1088/0067-0049/192/2/16
          [arXiv:1001.4635 [astro-ph.CO]].
          
          
        Keisler:2011aw 
          R. Keisler  et al.,
          
          Astrophys. J.  743, 28 (2011)
          doi:10.1088/0004-637X/743/1/28
          [arXiv:1105.3182 [astro-ph.CO]].
          
          
        Ade:2014afa 
          P. A. R. Ade  et al. [POLARBEAR Collaboration],
          
          Astrophys. J.  794, no. 2, 171 (2014)
        
          [arXiv:1403.2369 [astro-ph.CO]].
          
          
        Ade:2014xna 
          P. A. R. Ade  et al. [BICEP2 Collaboration],
          
          Phys. Rev. Lett.  112, no. 24, 241101 (2014)
        
          [arXiv:1403.3985 [astro-ph.CO]].
          
          
        Array:2015xqh 
          P. A. R. Ade  et al. [BICEP2 and Keck Array Collaborations],
          
          Phys. Rev. Lett.  116, 031302 (2016)
        
          [arXiv:1510.09217 [astro-ph.CO]].
          
          
        Ade:2015cao 
          P. A. R. Ade  et al. [POLARBEAR Collaboration],
          
          Phys. Rev. D  92, 123509 (2015)
        
          [arXiv:1509.02461 [astro-ph.CO]].
          
          
        Bonvin:2014xia 
          C. Bonvin, R. Durrer and R. Maartens,
          
          Phys. Rev. Lett.  112, no. 19, 191303 (2014)
          doi:10.1103/PhysRevLett.112.191303
          [arXiv:1403.6768 [astro-ph.CO]].
          
          
          
        
        Zucca:2016iur 
          A. Zucca, Y. Li and L. Pogosian,
          
          Phys. Rev. D  95, no. 6, 063506 (2017)
          doi:10.1103/PhysRevD.95.063506
          [arXiv:1611.00757 [astro-ph.CO]].
          
          
        Brown:2005kr 
          I. Brown and R. Crittenden,
          
          Phys. Rev. D  72, 063002 (2005)
          doi:10.1103/PhysRevD.72.063002
          [astro-ph/0506570].
          
          
          
        
        Seshadri:2009sy 
          T. R. Seshadri and K. Subramanian,
          
          Phys. Rev. Lett.  103, 081303 (2009)
          doi:10.1103/PhysRevLett.103.081303
          [arXiv:0902.4066 [astro-ph.CO]].
          
          
         
          
        Caprini:2009vk 
          C. Caprini, F. Finelli, D. Paoletti and A. Riotto,
          
          JCAP  0906, 021 (2009)
          doi:10.1088/1475-7516/2009/06/021
          [arXiv:0903.1420 [astro-ph.CO]].
          
          
          
          
        Cai:2010uw 
          R. G. Cai, B. Hu and H. B. Zhang,
          
          JCAP  1008, 025 (2010)
          doi:10.1088/1475-7516/2010/08/025
          [arXiv:1006.2985 [astro-ph.CO]].
          
          
        Trivedi:2010gi 
          P. Trivedi, K. Subramanian and T. R. Seshadri,
          
          Phys. Rev. D  82, 123006 (2010)
          doi:10.1103/PhysRevD.82.123006
          [arXiv:1009.2724 [astro-ph.CO]].
          
          
        Shiraishi:2010yk 
          M. Shiraishi, D. Nitta, S. Yokoyama, K. Ichiki and K. Takahashi,
          
          Phys. Rev. D  82, 121302 (2010)
          Erratum: [Phys. Rev. D  83, 029901 (2011)]
          doi:10.1103/PhysRevD.83.029901, 10.1103/PhysRevD.82.121302
          [arXiv:1009.3632 [astro-ph.CO]].
          
          
         
         
        Brown:2010jd 
          I. A. Brown,
          
          Astrophys. J.  733, 83 (2011)
          doi:10.1088/0004-637X/733/2/83
          [arXiv:1012.2892 [astro-ph.CO]].
          
          
         
         
        Shiraishi:2011dh 
          M. Shiraishi, D. Nitta, S. Yokoyama, K. Ichiki and K. Takahashi,
          
          Phys. Rev. D  83, 123003 (2011)
          doi:10.1103/PhysRevD.83.123003
          [arXiv:1103.4103 [astro-ph.CO]].
          
          
          
        
        Shiraishi:2012rm 
          M. Shiraishi, D. Nitta, S. Yokoyama and K. Ichiki,
          
          JCAP  1203, 041 (2012)
          doi:10.1088/1475-7516/2012/03/041
          [arXiv:1201.0376 [astro-ph.CO]].
          
          
        Shiraishi:2013wua 
          M. Shiraishi and T. Sekiguchi,
          
          Phys. Rev. D  90, no. 10, 103002 (2014)
          doi:10.1103/PhysRevD.90.103002
          [arXiv:1304.7277 [astro-ph.CO]].
          
          
         
         
        Shiraishi:2013vha 
          M. Shiraishi,
          
          JCAP  1311, 006 (2013)
          doi:10.1088/1475-7516/2013/11/006
          [arXiv:1308.2531 [astro-ph.CO]].
          
          
                   
          
        Trivedi:2011vt 
          P. Trivedi, T. R. Seshadri and K. Subramanian,
          
          Phys. Rev. Lett.  108, 231301 (2012)
          doi:10.1103/PhysRevLett.108.231301
          [arXiv:1111.0744 [astro-ph.CO]].
          
          
          
          
        Trivedi:2013wqa 
          P. Trivedi, K. Subramanian and T. R. Seshadri,
          
          Phys. Rev. D  89, no. 4, 043523 (2014)
          doi:10.1103/PhysRevD.89.043523
          [arXiv:1312.5308 [astro-ph.CO]].
          
          
        Ade:2013ydc 
          P. A. R. Ade  et al. [Planck Collaboration],
          
          Astron. Astrophys.  571, A24 (2014)
          doi:10.1051/0004-6361/201321554
          [arXiv:1303.5084 [astro-ph.CO]].
          
          
        Kunze:2014eka 
          K. E. Kunze and E. Komatsu,
          
          JCAP  1506, no. 06, 027 (2015)
        
          [arXiv:1501.00142 [astro-ph.CO]].
          
          
        Kunze:2017rsb 
          K. E. Kunze,
          
          Phys. Rev. D  96, no. 6, 063526 (2017)
          doi:10.1103/PhysRevD.96.063526
          [arXiv:1703.03650 [astro-ph.CO]].
          
        Chluba:2015lpa 
          J. Chluba, D. Paoletti, F. Finelli and J. A. Rubio-Martn,
          
          Mon. Not. Roy. Astron. Soc.  451, no. 2, 2244 (2015)
        
          [arXiv:1503.04827 [astro-ph.CO]].
          
          
              
        
          
        
        Jedamzik:1996wp 
          K. Jedamzik, V. Katalinic and A. V. Olinto,
          
          Phys. Rev. D  57, 3264 (1998)
          doi:10.1103/PhysRevD.57.3264
          [astro-ph/9606080].
          
          
        Jedamzik:1999bm 
          K. Jedamzik, V. Katalinic and A. V. Olinto,
          
          Phys. Rev. Lett.  85, 700 (2000)
          [astro-ph/9911100].
          
          
          
        
        Kunze:2013uja 
          K. E. Kunze and E. Komatsu,
          
          JCAP  1401, 009 (2014)
        
          [arXiv:1309.7994 [astro-ph.CO]].
          
          
        Miyamoto:2013oua 
          K. Miyamoto, T. Sekiguchi, H. Tashiro and S. Yokoyama,
          
          Phys. Rev. D  89, no. 6, 063508 (2014)
          doi:10.1103/PhysRevD.89.063508
          [arXiv:1310.3886 [astro-ph.CO]].
          
          
         
         
        Ganc:2014wia 
          J. Ganc and M. S. Sloth,
          
          JCAP  1408, 018 (2014)
          doi:10.1088/1475-7516/2014/08/018
          [arXiv:1404.5957 [astro-ph.CO]].
          
          
           
        
        Wagstaff:2015jaa 
          J. M. Wagstaff and R. Banerjee,
          
          Phys. Rev. D  92, no. 12, 123004 (2015)
          doi:10.1103/PhysRevD.92.123004
          [arXiv:1508.01683 [astro-ph.CO]].
          
          
        
        Minoda:2017iob 
          T. Minoda, K. Hasegawa, H. Tashiro, K. Ichiki and N. Sugiyama,
          
          arXiv:1705.10054 [astro-ph.CO].
          
            
        
        Jedamzik:2013gua 
          K. Jedamzik and T. Abel,
          
          JCAP  1310, 050 (2013).
          doi:10.1088/1475-7516/2013/10/050
          
          
        Jedamzik:2018itu 
          K. Jedamzik and A. Saveliev,
          
        arXiv:1804.06115 [astro-ph.CO].
        
            
        
        Shaw:2010ea 
          J. R. Shaw and A. Lewis,
          
          Phys. Rev. D  86, 043510 (2012)
          doi:10.1103/PhysRevD.86.043510
          [arXiv:1006.4242 [astro-ph.CO]].
          
          
          
          
        Kahniashvili:2010wm 
          T. Kahniashvili, A. G. Tevzadze, S. K. Sethi, K. Pandey and B. Ratra,
          
          Phys. Rev. D  82, 083005 (2010)
          doi:10.1103/PhysRevD.82.083005
          [arXiv:1009.2094 [astro-ph.CO]].
          
          
          
        
        Kahniashvili:2012dy 
          T. Kahniashvili, Y. Maravin, A. Natarajan, N. Battaglia and A. G. Tevzadze,
          
          Astrophys. J.  770, 47 (2013)
          doi:10.1088/0004-637X/770/1/47
          [arXiv:1211.2769 [astro-ph.CO]].
          
          
        Yamazaki:2016sov 
          D. G. Yamazaki,
          
          Phys. Rev. D  93, no. 4, 043004 (2016)
          doi:10.1103/PhysRevD.93.043004
          [arXiv:1602.02750 [astro-ph.CO]].
          
          
          
        
        Neronov:1900zz 
          A. Neronov and I. Vovk,
          
          Science  328, 73 (2010)
          [arXiv:1006.3504 [astro-ph.HE]].
          
          
        	
        Tashiro:2013ita 
          H. Tashiro, W. Chen, F. Ferrer and T. Vachaspati,
          
          Mon. Not. Roy. Astron. Soc.  445, no. 1, L41 (2014)
          
          [arXiv:1310.4826 [astro-ph.CO]].
          
          
        Berezhiani:2003ik 
          Z. Berezhiani and A. D. Dolgov,
          
          Astropart. Phys.  21, 59 (2004)
          
          [astro-ph/0305595].
          
          
        Gnedin:2000ax 
          N. Y. Gnedin, A. Ferrara and E. G. Zweibel,
          
          Astrophys. J.  539, 505 (2000)
          
          [astro-ph/0001066].
          
          
        Langer:2005gw 
          M. Langer, N. Aghanim and J. L. Puget,
          
          Astron. Astrophys.  443, 367 (2005)
          
          [astro-ph/0508173].
          
          
        Ando:2010ry 
          M. Ando, K. Doi and H. Susa,
          
          Astrophys. J.  716, 1566 (2010)
          
          [arXiv:1005.1123 [astro-ph.CO]].
          
          
        Durrive:2015cja 
          J. B. Durrive and M. Langer,
          
          Mon. Not. Roy. Astron. Soc.  453, no. 1, 345 (2015)
          
          [arXiv:1506.08177 [astro-ph.CO]].
          
          
        Abazajian:2016yjj 
          K. N. Abazajian  et al. [CMB-S4 Collaboration],
          
          arXiv:1610.02743 [astro-ph.CO].
          
          
          
          
        Yadav:2012uz 
          A. Yadav, L. Pogosian and T. Vachaspati,
          
          Phys. Rev. D  86, 123009 (2012)
        
          [arXiv:1207.3356 [astro-ph.CO]].
          
          
          
        
        De:2013dra 
          S. De, L. Pogosian and T. Vachaspati,
          
          Phys. Rev. D  88, no. 6, 063527 (2013)
        
          [arXiv:1305.7225 [astro-ph.CO]].
          
          
          
          
        Pogosian:2013dya 
          L. Pogosian,
          
          Mon. Not. Roy. Astron. Soc.  438, no. 3, 2508 (2014)
         
          [arXiv:1311.2926 [astro-ph.CO]].
          
          
          
          
        Subramanian:2015lua 
          K. Subramanian,
          
          Rept. Prog. Phys.  79, no. 7, 076901 (2016)
          
          [arXiv:1504.02311 [astro-ph.CO]].
          
          
        Barrow:1997mj 
          J. D. Barrow, P. G. Ferreira and J. Silk,
          
          Phys. Rev. Lett.  78, 3610 (1997)
          
          [astro-ph/9701063].
          
          
        Barrow:1997sy 
          J. D. Barrow,
          
          Phys. Rev. D  55, 7451 (1997)
          doi:10.1103/PhysRevD.55.7451
          [gr-qc/9701038].
          
          
          
          
          
        Kahniashvili:2005xe 
          T. Kahniashvili and B. Ratra,
          
          Phys. Rev. D  71, 103006 (2005)
        
          [astro-ph/0503709].
          
          
        Kunze:2011bp 
          K. E. Kunze,
          
          Phys. Rev. D  85, 083004 (2012)
          
          [arXiv:1112.4797 [astro-ph.CO]].
          
          
        Ballardini:2014jta 
          M. Ballardini, F. Finelli and D. Paoletti,
          
          JCAP  1510, no. 10, 031 (2015)
          
          [arXiv:1412.1836 [astro-ph.CO]].
          
          
        
        Bonvin:2013tba 
          C. Bonvin, C. Caprini and R. Durrer,
          
          Phys. Rev. D  88, 083515 (2013)
          
          [arXiv:1308.3348 [astro-ph.CO]].
          
          
        Durrer:2003ja 
          R. Durrer and C. Caprini,
          
          JCAP  0311, 010 (2003)
        
          [astro-ph/0305059].
          
          
        Jedamzik:2010cy 
          K. Jedamzik and G. Sigl,
          
          Phys. Rev. D  83, 103005 (2011)
         
          [arXiv:1012.4794 [astro-ph.CO]].
          
          
        Kahniashvili:2009qi 
          T. Kahniashvili, A. G. Tevzadze and B. Ratra,
          
          Astrophys. J.  726, 78 (2011)
          
          [arXiv:0907.0197 [astro-ph.CO]].
          
          
        Kamionkowski:1996zd 
          M. Kamionkowski, A. Kosowsky and A. Stebbins,
          
          Phys. Rev. Lett.  78, 2058 (1997)
          [astro-ph/9609132].
          
        Seljak:1996gy 
          U. Seljak and M. Zaldarriaga,
          
          Phys. Rev. Lett.  78, 2054 (1997)
          [astro-ph/9609169].
          
          
        
        Kodama:1985bj 
          H. Kodama and M. Sasaki,
          
          Prog. Theor. Phys. Suppl.  78, 1 (1984).
        
          
          
        Hu:1997hp 
          W. Hu and M. J. White,
          
          Phys. Rev. D  56, 596 (1997)
        
          [astro-ph/9702170].
          
          
          
          
        
        Giovannini:2004aw 
          M. Giovannini,
          
          Phys. Rev. D  70, 123507 (2004)
          
          [astro-ph/0409594].
          
          
        Giovannini:2006kc 
          M. Giovannini,
          
          Phys. Rev. D  74, 063002 (2006)
          
          [astro-ph/0606759].
          
          
          
          
        
        Lewis:1999bs 
          A. Lewis, A. Challinor and A. Lasenby,
          
          Astrophys. J.  538, 473 (2000)
          
          [astro-ph/9911177].
          
          
          
          
          
        Harari:1996ac 
          D. D. Harari, J. D. Hayward and M. Zaldarriaga,
          
          Phys. Rev. D  55, 1841 (1997)
          [astro-ph/9608098].
          
          
          
          
          
        Kamionkowski:2008fp
          M. Kamionkowski,
          
          Phys. Rev. Lett.  102, 111302 (2009)
          [arXiv:0810.1286 [astro-ph]].
          
          
        
        Yadav:2009eb 
          A. P. S. Yadav, R. Biswas, M. Su and M. Zaldarriaga,
          
          Phys. Rev. D  79, 123009 (2009)
          [arXiv:0902.4466 [astro-ph.CO]].
          
          
        Gluscevic:2009mm
          V. Gluscevic, M. Kamionkowski and A. Cooray,
          
          
          Phys. Rev.  D  80, 023510 (2009)
          [arXiv:0905.1687 [astro-ph.CO]].
          
        Gluscevic:2012me 
          V. Gluscevic, D. Hanson, M. Kamionkowski and C. M. Hirata,
          
          arXiv:1206.5546 [astro-ph.CO].
          
          
         
        Hu:2001kj 
          W. Hu and T. Okamoto,
          
          Astrophys. J.  574, 566 (2002)
          [astro-ph/0111606].
          
          
          
        
        Seljak:2003pn 
          U. Seljak and C. M. Hirata,
          
          Phys. Rev. D  69, 043005 (2004)
          [astro-ph/0310163].
          
          
         
         
        Oppermann:2011td 
          N. Oppermann, H. Junklewitz, G. Robbers, M. R. Bell, T. A. Ensslin, A. Bonafede, R. Braun and J. C. Brown  et al.,
          
          arXiv:1111.6186 [astro-ph.GA].
          
          
         
        Kahniashvili:2008hx 
          T. Kahniashvili, Y. Maravin and A. Kosowsky,
          
          Phys. Rev. D  80, 023009 (2009)
          [arXiv:0806.1876 [astro-ph]].
          
          
           
        
        Kovac:2002fg 
          J. Kovac, E. M. Leitch, C. Pryke, J. E. Carlstrom, N. W. Halverson and W. L. Holzapfel,
          
          Nature  420, 772 (2002)
          [astro-ph/0209478].
          
          
        Kogut:2003et 
          A. Kogut  et al. [WMAP Collaboration],
          
          Astrophys. J. Suppl.  148, 161 (2003)
          doi:10.1086/377219
          [astro-ph/0302213].
          
          
        Page:2006hz 
          L. Page  et al. [WMAP Collaboration],
          
          Astrophys. J. Suppl.  170, 335 (2007)
          doi:10.1086/513699
          [astro-ph/0603450].
          
          
            
          
        
        Aghanim:2015xee 
          N. Aghanim  et al. [Planck Collaboration],
          
          
          [arXiv:1507.02704 [astro-ph.CO]].
          
          
          
          
        Suzuki:2015zzg 
          A. Suzuki  et al. [POLARBEAR Collaboration],
          
          J. Low. Temp. Phys.  184, no. 3-4, 805 (2016)
          doi:10.1007/s10909-015-1425-4
          [arXiv:1512.07299 [astro-ph.IM]].
          
          
        Benson:2014qhw 
          B. A. Benson  et al. [SPT-3G Collaboration],
          
          Proc. SPIE Int. Soc. Opt. Eng.  9153, 91531P (2014)
          doi:10.1117/12.2057305
          [arXiv:1407.2973 [astro-ph.IM]].
          
          
        Sutton:2017jgr 
          D. R. Sutton, C. Feng and C. L. Reichardt,
          
          Astrophys. J.  846, no. 2, 164 (2017)
          doi:10.3847/1538-4357/aa85e2
          [arXiv:1702.01871 [astro-ph.CO]].
          
          
          
        
        Baumann:2008aq 
          D. Baumann  et al. [CMBPol Study Team],
          
          AIP Conf. Proc.  1141, 10 (2009)
          doi:10.1063/1.3160885
          [arXiv:0811.3919 [astro-ph]].
          
          
          
        
        Delabrouille:2017rct 
          J. Delabrouille  et al. [CORE Collaboration],
          
          arXiv:1706.04516 [astro-ph.IM].
          
          
        Crittenden:1993ni 
          R. Crittenden, J. R. Bond, R. L. Davis, G. Efstathiou and P. J. Steinhardt,
          
          Phys. Rev. Lett.  71, 324 (1993)
          [astro-ph/9303014].
          
          
            
        
        Abazajian:2013oma 
          K. N. Abazajian  et al. [Topical Conveners: K.N. Abazajian, J.E. Carlstrom, A.T. Lee Collaboration],
          
          Astropart. Phys.  63, 66 (2015)
        
          [arXiv:1309.5383 [astro-ph.CO]].
          
          
        Amendola:2014wma 
          L. Amendola, G. Ballesteros and V. Pettorino,
          
          Phys. Rev. D  90, 043009 (2014)
        
          [arXiv:1405.7004 [astro-ph.CO]].
          
          
        Raveri:2014eea 
          M. Raveri, C. Baccigalupi, A. Silvestri and S. Y. Zhou,
          
          Phys. Rev. D  91, no. 6, 061501 (2015)
        
          [arXiv:1405.7974 [astro-ph.CO]].
          
          
        Seljak:1997ii 
          U. Seljak, U. -L. Pen and N. Turok,
          
          Phys. Rev. Lett.  79, 1615 (1997)
          [astro-ph/9704231].
          
          
        
        Avgoustidis:2011ax 
          A. Avgoustidis, E. J. Copeland, A. Moss, L. Pogosian, A. Pourtsidou and D. A. Steer,
          
          Phys. Rev. Lett.  107, 121301 (2011)
        
          [arXiv:1105.6198 [astro-ph.CO]].
          
          
        Moss:2014cra 
          A. Moss and L. Pogosian,
          
          Phys. Rev. Lett.  112, 171302 (2014)
         
          [arXiv:1403.6105 [astro-ph.CO]].
          
          
         
         
        Lizarraga:2016onn 
          J. Lizarraga, J. Urrestilla, D. Daverio, M. Hindmarsh and M. Kunz,
          
          JCAP  1610, no. 10, 042 (2016)
          doi:10.1088/1475-7516/2016/10/042
          [arXiv:1609.03386 [astro-ph.CO]].
          
          
        
        Pospelov:2008gg 
          M. Pospelov, A. Ritz, C. Skordis, A. Ritz and C. Skordis,
          
          Phys. Rev. Lett.  103, 051302 (2009)
          doi:10.1103/PhysRevLett.103.051302
          [arXiv:0808.0673 [astro-ph]].
          
          
        Shimon:2007au 
          M. Shimon, B. Keating, N. Ponthieu and E. Hivon,
          
          Phys. Rev. D  77, 083003 (2008)
          doi:10.1103/PhysRevD.77.083003
          [arXiv:0709.1513 [astro-ph]].
          
          CPR_paper
        L. Pogosian, M. Shimon, M. Mewes and B. Keating, in preparation.
            
        westgridcompute
