Radiative transport of relativistic species in cosmology

TL;DR

This paper reviews the construction of distribution functions for fermions and photons in cosmology, emphasizing their similarities and differences, and introduces methods to handle polarization, anisotropy, and spectral distortions in the Boltzmann equations.

Contribution

It provides a unified framework for fermion and photon distribution functions, including polarization and anisotropy, and extends the Kompaneets equation to more general cases.

Findings

Unified treatment of fermion and photon distribution functions.

Extension of the Kompaneets equation to anisotropic and polarized photons.

Parameterization of photon spectra using logarithmic moments.

Abstract

We review the general construction of distribution functions for gases of fermions and bosons (photons), emphasizing the similarities and differences between both cases. The central object which describes polarization for photons is a tensor-valued distribution function, whereas for fermions it is a vector-valued one. The collision terms of Boltzmann equations for fermions and bosons also possess the same general structure and differ only in the quantum effects associated with the final state of the reactions described. In particular, neutron-proton conversions in the early universe, which set the primordial Helium abundance, enjoy many similarities with Compton scattering which shapes the cosmic microwave background and we show that both can be handled with a Fokker-Planck type expansion. For neutron-proton conversions, this allows to obtain the finite nucleon mass corrections, required for precise theoretical predictions, whereas for Compton scattering it leads to the thermal and recoil effects which enter the Kompaneets equation. We generalize the latter to the general case of anisotropic and polarized photon distribution functions. Finally we discuss a parameterization of the photon spectrum based on logarithmic moments which allows for a neat separation between temperature shifts and spectral distortions.

Figures

Tables

Table 1. Table 1: Main two-body reactions for the collision term of neutrinos. The electronic neutrino is noted ν 𝜈 \nu and the muonic neutrino is noted ν μ subscript 𝜈 𝜇 \nu_{\mu} . Similar reactions for antineutrinos can be deduced with a global charge conjugation on all these reactions, and thus on all the collision Kernels. The squared amplitudes | M | 2 superscript 𝑀 2 |M|^{2} are expressed with Eq. ( III.4 ). The integration on momenta is defined in Eq. ( 144 ).

Reaction	Particles names	Chiral couplings	$2^{-9}{G}_F^{-2}E{I}_{C\left[A\right]}=𝒦|M{|}^2\times$
$\nu +{\nu }_{\mu }\leftrightarrow \nu +{\nu }_{\mu }$	$a+b\leftrightarrow c+d$	${ϵ}_{-}^{ac}=\frac{1}{2}{ϵ}_{+}^{ac}=0{ϵ}_{-}^{bd}=\frac{1}{2}{ϵ}_{+}^{bd}=0$	$\left({\hat{f}}_a{\hat{f}}_b{f}_c{f}_d-f_a{f}_b{\hat{f}}_c{\hat{f}}_d\right)$
$\nu +{\overline{\nu }}_{\mu }\leftrightarrow \nu +{\overline{\nu }}_{\mu }$	$a+\overline{d}\leftrightarrow c+\overline{b}$	${ϵ}_{-}^{ac}=\frac{1}{2}{ϵ}_{+}^{ac}=0{ϵ}_{-}^{bd}=\frac{1}{2}{ϵ}_{+}^{bd}=0$	$\left({\hat{f}}_a{\overline{f}}_b{f}_c{\hat{\overline{f}}}_d-f_a{\hat{\overline{f}}}_b{\hat{f}}_c{\overline{f}}_d\right)$
$\nu +\overline{\nu }\leftrightarrow {\overline{\nu }}_{\mu }+{\nu }_{\mu }$	$a+\overline{c}\leftrightarrow \overline{b}+d$	${ϵ}_{-}^{ac}=\frac{1}{2}{ϵ}_{+}^{ac}=0{ϵ}_{-}^{bd}=\frac{1}{2}{ϵ}_{+}^{bd}=0$	$\left({\hat{f}}_a{\overline{f}}_b{\hat{\overline{f}}}_c{f}_d-f_a{\hat{\overline{f}}}_b{\overline{f}}_c{\hat{f}}_d\right)$
$\nu +\nu \leftrightarrow \nu +\nu$	$a+b\leftrightarrow c+d$	${ϵ}_{-}^{ac}=\frac{1}{2}{ϵ}_{+}^{ac}=0{ϵ}_{-}^{bd}=\frac{1}{2}{ϵ}_{+}^{bd}=0$	$2\left({\hat{f}}_a{\hat{f}}_b{f}_c{f}_d-f_a{f}_b{\hat{f}}_c{\hat{f}}_d\right)$
$\nu +\overline{\nu }\leftrightarrow \nu +\overline{\nu }$	$a+\overline{d}\leftrightarrow c+\overline{b}$	${ϵ}_{-}^{ac}=\frac{1}{2}{ϵ}_{+}^{ac}=0{ϵ}_{-}^{bd}=\frac{1}{2}{ϵ}_{+}^{bd}=0$	$4\left({\hat{f}}_a{\overline{f}}_b{f}_c{\hat{\overline{f}}}_d-f_a{\hat{\overline{f}}}_b{\hat{f}}_c{\overline{f}}_d\right)$
$\nu +e^{-}\leftrightarrow \nu +e^{-}$	$a+b\leftrightarrow c+d$	${ϵ}_{-}^{ac}=\frac{1}{2}{ϵ}_{+}^{ac}=0{ϵ}_{-}^{bd}={ϵ}_{-}^e+1{ϵ}_{+}^{bd}={ϵ}_{+}^e$	$\left({\hat{f}}_a{\hat{f}}_b{f}_c{f}_d-f_a{f}_b{\hat{f}}_c{\hat{f}}_d\right)$
$\nu +e^{+}\leftrightarrow \nu +e^{+}$	$a+\overline{d}\leftrightarrow c+\overline{b}$	${ϵ}_{-}^{ac}=\frac{1}{2}{ϵ}_{+}^{ac}=0{ϵ}_{-}^{bd}={ϵ}_{-}^e+1{ϵ}_{+}^{bd}={ϵ}_{+}^e$	$\left({\hat{f}}_a{\overline{f}}_b{f}_c{\hat{\overline{f}}}_d-f_a{\hat{\overline{f}}}_b{\hat{f}}_c{\overline{f}}_d\right)$
$\nu +\overline{\nu }\leftrightarrow e^{+}+e^{-}$	$a+\overline{c}\leftrightarrow \overline{b}+d$	${ϵ}_{-}^{ac}=\frac{1}{2}{ϵ}_{+}^{ac}=0{ϵ}_{-}^{bd}={ϵ}_{-}^e+1{ϵ}_{+}^{bd}={ϵ}_{+}^e$	$\left({\hat{f}}_a{\overline{f}}_b{\hat{\overline{f}}}_c{f}_d-f_a{\hat{\overline{f}}}_b{\overline{f}}_c{\hat{f}}_d\right)$

Equations

$$E=p^0,\Rightarrow E^2=m^2+|\mathbf{p}{|}^2.$$

$$[\mathrm{d}p]\equiv \frac{{\mathrm{d}}^3\mathbf{p}}{(2\pi {)}^{3}2p^0},$$

$$\underline{\delta }(p-p^{\prime })=(2\pi {)}^{3}2p^0{\delta }^3(\mathbf{p}-{\mathbf{p}}^{\prime }),$$

$$\eta =\mathrm{diag}(-,+,+,+).$$

$${ϵ}_{0123}=-{ϵ}^{0123}=1.$$

$$N_{rs}(\mathbf{p},{\mathbf{p}}^{\prime })\equiv a_r^{†}(\mathbf{p})a_s({\mathbf{p}}^{\prime }).$$

$$N_{rs}\equiv \int [\mathrm{d}p]N_{rs}(p,p).$$

$$⟨\Psi |N_{rs}(\mathbf{p},{\mathbf{p}}^{\prime })|\Psi ⟩=\underline{\delta }(p-p^{\prime })f_{rs}(\mathbf{p}).$$

$$⟨\Psi |N_{rs}|\Psi ⟩=\mathcal{V}\int \frac{{\mathrm{d}}^3\mathbf{p}}{(2\pi {)}^3}{f}_{rs}(p),$$

$$N_{rs}^{\star }(p,p^{\prime })=N_{sr}(p^{\prime },p)\Rightarrow f_{rs}^{\star }(p)=f_{sr}(p).$$

$$\mathbf{n}\equiv \frac{\mathbf{p}}{|\mathbf{p}|}.$$

$$n^i$$

$$e_{\theta }^i$$

$$S^{\mu }(u^{\nu },p^{\nu })=-\frac{m}{|\mathbf{p}|}{u}^{\mu }+\frac{E}{m|\mathbf{p}|}{p}^{\mu },$$

$$S^0=m^{-1}|\mathbf{p}|,S^i=m^{-1}E{n}^i.$$

$${ϵ}_{\pm }^{\mu }\equiv \frac{1}{\sqrt{2}}(e_{\theta }^{\mu }\mp \mathrm{i}{e}_{\varphi }^{\mu }),$$

$$\displaystyle{\epsilon}^{0}_{\pm}=0\,,\quad{\epsilon}^{i}_{\pm}=\frac{1}{\sqrt{2}}\left(\begin{array}[]{c}\text{cos}(\theta)\text{cos}(\phi)\pm i\text{sin}(\phi)\\ \text{cos}(\theta)\text{sin}(\phi)\mp i\text{cos}(\phi)\\ -\text{sin}(\theta)\end{array}\right)\,,$$

$$\psi =\sum _{s=\pm \frac{1}{2}}\int [\mathrm{d}p]\left[{\mathrm{e}}^{-\mathrm{i}p\cdot x}{\bar{a}}_s^{†}(p)v_s(p)+{\mathrm{e}}^{\mathrm{i}p\cdot x}{a}_s(p)u_s(p)\right],$$

$$\{a_{r}(p),a^{\dagger}_{s}(p^{\prime})\}&=\delta^{\rm K}_{rs}\underline{\delta}(p-p^{\prime})\,,\\ \{{\bar{a}}_{r}(p),{\bar{a}}^{\dagger}_{s}(p^{\prime})\}=\delta^{\rm K}_{rs}\underline{\delta}(p-p^{\prime})\,,$$

$$F_{\mathfrak{a}}^{\mathfrak{b}}(p)\equiv \{\begin{array}{l}\sum _{rs}{f}_{rs}(p)u_{s,\mathfrak{a}}(p){\bar{u}}_r^{\mathfrak{b}}(p)part.,\\ \sum _{rs}{f}_{rs}(p)v_{r,\mathfrak{a}}(p){\bar{v}}_s^{\mathfrak{b}}(p)antipart.\end{array}$$

$$(/p+m)u_s(p)=0,(/p-m)v_s(p)=0,$$

$$\mathcal{O}\equiv \{\text{mathds}1,{\gamma }^{\mu },{\Sigma }^{\mu \nu },{\gamma }^{\mu }{\gamma }^5,{\gamma }^5\},$$

$$F_{\mathfrak{a}}^{\mathfrak{b}}=\{\begin{array}{l}\sum _{X\in \mathcal{O}}{c}_X\sum _{rs}{f}_{rs}({\bar{u}}_{r}X{u}_s)X_{\mathfrak{a}}^{\mathfrak{b}}part.\\ \sum _{X\in \mathcal{O}}{c}_X\sum _{rs}{f}_{sr}({\bar{v}}_{r}X{v}_s)X_{\mathfrak{a}}^{\mathfrak{b}}antipart.\end{array}$$

$$\displaystyle u_{s}=\left(\begin{array}[]{c}\sqrt{E+2s|\vec{p}|}\chi_{s}\\ \sqrt{E-2s|\vec{p}|}\chi_{s}\end{array}\right)\quad v_{s}=\left(\begin{array}[]{c}2s\sqrt{E-2s|\vec{p}|}\chi_{-s}\\ -2s\sqrt{E+2s|\vec{p}|}\chi_{-s}\end{array}\right),$$

$$\chi_{\tfrac{1}{2}}=\left(\begin{array}[]{c}{\rm e}^{-{\rm i}\phi/2}{\rm cos}(\theta/2)\\ {\rm e}^{{\rm i}\phi/2}{\rm sin}(\theta/2)\end{array}\right)\quad\chi_{-\tfrac{1}{2}}=\left(\begin{array}[]{c}-{\rm e}^{-{\rm i}\phi/2}{\rm sin}(\theta/2)\\ {\rm e}^{{\rm i}\phi/2}{\rm cos}(\theta/2)\end{array}\right),$$

$$u_{s}\bar{u}_{s}&=\frac{1}{2}\left(\mathds{1}-2s\gamma_{5}\not{S}\right)\left(-\not{p}+m\right)\,,\\ v_{s}\bar{v}_{s}=\frac{1}{2}\left(\mathds{1}-2s\gamma_{5}\not{S}\right)\left(-\not{p}-m\right)\,.$$

$$\sum _s{u}_s{\bar{u}}_s=-/p+m,\sum _s{v}_s{\bar{v}}_s=-/p-m.$$

$$u_{s}\bar{u}_{r}&=\frac{1}{\sqrt{2}}\gamma^{5}\not{\epsilon}_{r-s}\left(\not{p}-m\right)\quad{\rm if}\quad r\neq s\,,\\ v_{s}\bar{v}_{r}=\frac{1}{\sqrt{2}}\gamma^{5}\not{\epsilon}_{s-r}\left(\not{p}+m\right)\quad{\rm if}\quad r\neq s\,.$$

$$I$$

$$Q_{\pm }$$

$$Q^{\mu }\equiv Q_{+}{ϵ}_{+}^{\mu }+Q_{-}{ϵ}_{-}^{\mu },{\mathcal{Q}}^{\mu }\equiv Q^{\mu }+V{S}^{\mu }.$$

$$\mathbf{F}=\frac{1}{2}\left(I-{\gamma }^5{\gamma }^{\mu }{\mathcal{Q}}_{\mu }\right)\left(M-/p\right),$$

Full text

1. Introduction

2. Theoretical framework

3. I Distribution functions

4. I.1 General construction

5. I.1.1 Notation

6. I.1.2 Number operator

7. I.1.3 Adapted orthonormal basis

8. I.1.4 Fermions

9. I.1.5 Massless fermions

10. I.1.6 Massless bosons

11. I.2 Multipolar decomposition

12. I.2.1 Fermions

13. I.2.2 Massless bosons

14. I.3 Observer independence

15. I.3.1 Massive fermions

16. I.3.2 Massless fermions

17. I.3.3 Massless bosons

18. I.3.4 Relation to abstract tensor indices

19. II Boltzmann equation

20. II.1 Liouville equation in curved space-time

21. II.1.1 Massive fermions

22. II.1.2 Massless fermions

23. II.1.3 Massless bosons

24. II.2 Quantum evolution in the interaction picture

25. II.3 Molecular chaos

26. Weak interactions

27. III General collision term

28. III.1 Fermi theory of weak interactions

29. III.1.1 Neutral currents

30. III.1.2 Charged currents

31. III.2 Two-body processes

32. III.2.1 Notation

33. III.2.2 Collision term structure

34. III.3 General collision term

35. III.4 Standard reactions with neutrinos

36. III.4.1 Muon decay

37. III.4.2 neutrino/muon neutrino scattering

38. III.4.3 neutrino/neutrino scattering

39. III.4.4 neutrino/electron scattering

40. IV Neutrons-protons conversions

41. IV.1 General expression of the rates

42. IV.2 Isotropy of distributions

43. IV.3 Expansion in the energy transfer

44. IV.4 Lowest order $n\leftrightarrow p$ reaction rate

45. Compton scattering

46. V Compton collision term

47. V.1 Extended Klein-Nishina formula

48. V.2 Kinematics

49. V.3 Electron velocity distribution

50. V.4 Fokker-Planck expansion

51. V.4.1 Energy shift

52. V.4.2 Electron distribution

53. V.4.3 Squared amplitude

54. V.4.4 Integration measure

55. V.5 Structure of the expansion

56. VI Evolution of isotropic distributions

57. VI.1 Kompaneets equations

58. VI.2 Thermalization on electrons

59. VI.2.1 Test distribution

60. VI.2.2 Dominant distribution

61. VI.2.3 Compton cooling

62. VII Anisotropic distribution functions

63. VII.1 Choice of frame

64. VII.2 Thomson term

65. VII.3 Thermal effects

66. VII.4 Recoil effects

67. Spectral distortions

68. VIII Spectrum parameterization

69. VIII.1 Distribution of Planck spectra

70. VIII.1.1 Temperature transform

71. VIII.1.2 Spectral moments

72. VIII.1.3 Spectral moment of a polarized spectrum

73. VIII.1.4 Discussion on the choice of a set of variables

74. VIII.2 Spectral moments evolution

75. VIII.2.1 General form of the Boltzmann equation

76. VIII.2.2 Doppler, SZ effect and $y$-type distortion

77. VIII.2.3 Structure of the numerics

78. VIII.3 Summary and notation

79. IX Angular correlations of spectral distortions

80. IX.1 Collision term of spectral moments

81. IX.2 Non-linear kSZ effect during reionization

82. Appendices

83. A Spinor valued operators

84. B Construction of the classical collision term

85. C Finite nucleon mass corrections

86. D Symmetric trace-free (STF) tensors

87. D.1 Notation

88. D.2 Relation to spherical harmonics

89. D.3 Relation to spin-weighted spherical harmonics

90. D.4 Products and contractions of STF tensor

91. D.5 The ${\cal Y}_{\ell m}^{L}$ in the literature

92. E Collision term in a different frame

93. E.1 Intensity

94. E.2 Linear polarization

95. E.3 Circular polarization

Introduction

We review theoretical and practical aspects of the radiative transport of relativistic species. Our emphasis is on cosmological applications, hence we focus on fermions (neutrons, protons, electrons, positrons and neutrinos) during big-bang nucleosynthesis (BBN) and photons of the cosmic microwave background (CMB). Relativistic species cannot be described by perfect fluids and one must account for the distribution of particles using a distribution function $f({\bm{x}},{\bm{p}},t)$, whose evolution is dictated by a Boltzmann equation $L[f]=C[f]$. The left hand side is the Liouville term and describes the free streaming of particles. In a curved space-time, this requires the use of cosmological perturbation theory, that is general relativity. The emphasis of this article is on the right hand side which is the collision term, and describes the evolution of the distributions under the influence of collisions, that is because of the micro-physics. Hence all the results presented here are independent of any perturbation theory, as they are derived from the basic principles of particle physics. From the equivalence principle, they are formulated in a local orthonormal frame, that is in the context of special relativity.

It is instructive to consider the cases of fermions and bosons side by side as their description by distribution functions have numerous similarities. In fact, the case of massless fermions is simpler than the case of photons in many respects, essentially because the spin of fermions ($1/2$) is smaller than the spin of photons ($1$). Hence the paper is organized to allow for a detailed comparison of these two cases. We show that the collision term for weak interactions during BBN has a structure which is extremely similar to the structure of the collision term for photons due to Compton scattering. Furthermore, we can use common techniques to express in practice these collision terms in functions of the distribution functions moments. Even though one could follow the analogy for anisotropic distribution functions, it is not useful for the case of fermions in the context of BBN. Hence we study in details the angular structure only for Compton scattering and derive the extended Kompaneets equation, valid for anisotropic and polarized photon distributions. Since the emphasis is on the derivation and the structure of the equations, we only summarize how the equations must be applied in the cosmological context, overviewing briefly the main physical effects.

In § I a general procedure to build a classical distribution function out of the quantum number operator is summarized. We then detail in § II how a classical Boltzmann equation can be derived, given a set of suitable approximation and assumptions, from the quantum evolution of the number operator. Sections III is then dedicated to the collision terms of weak interactions processes for fermions in the early universe. In order to compute in practice these collision terms, we review the Fokker-Planck expansion in § IV applied to BBN weak interactions. The similar treatment of Compton interactions between electrons and photons and its Fokker-Planck expansion are reviewed in § V. Applications for the evolution of isotropic photon distributions under Compton interactions are presented in § VI, with a brief discussion on its implications for cosmology. The detailed form of the collision term for anisotropic distributions, including polarization is subsequently exposed in § VII, using symmetric trace-free tensors to decompose the angular dependence. It is the first general derivation of the Compton collision term in the literature which includes thermal and recoil effects while describing consistently polarization. Finally we present in § VIII a parameterization for spectral distortions and we collect in § IX the equations governing the generation of distortions from the Thomson part of the collision term, with plots of the associated angular power spectrum generated during the reionization era.

Theoretical framework

I Distribution functions

In this section, based on Fidler and Pitrou (2017), we review how the distribution function is built from the quantum expectations of the number operator, and how its covariant components can be extracted. We also show that for each spin there is an adapted expansion in spin-weighted spherical harmonics for the dependence on the spatial momentum direction. The case of fermions is presented first, even though it is less known, as it allows to understand better the photon case.

I.1 General construction

I.1.1 Notation

Before considering the kinetic theory in curved spacetime, we build the formalism in a flat space-time (that is the Minkowski space-time of special relativity) in which the quantum theory of particles is very well established. An inertial frame is defined by a tetrad field, that is by a timelike vector field $e_{0}$ and three spacelike vector fields $e_{i}$, together with the associated co-tetrad $e^{0},e^{i}$. Latin indices such as $i,j,\dots$ indicate spatial components in the tetrad basis. A four-vector is written as $V^{\mu}$ where Greek indices $\mu,\nu,\dots$ denote components in the tetrad basis. In particular, the components of the tetrad vectors and co-vectors in the tetrad basis are by definition $[e_{\mu}]^{\nu}=\delta_{\mu}^{\nu}$ and $[e^{\mu}]_{\nu}=\delta^{\mu}_{\nu}$. If gravity can be ignored, that is in the context of special relativity, the inertial frame is global. Later, when including the effect of gravity in the context of general relativity, the inertial frame is local and one must employ general coordinates whose indices are labelled by $\alpha,\beta,\dots$. For a given vector this implies $V^{\alpha}=V^{\mu}[e_{\mu}]^{\alpha}$.

The momentum vector $p^{\mu}$ will often simply be denoted as $p$ and its spatial components $p^{i}$ allow to build the spatial momentum ${\bm{p}}=p^{i}{\bf e}_{i}$. More generally, we reserve boldface notation to spatial vectors. The energy associated with the momentum is given by the time component

[TABLE]

When a quantity depends on the spatial momentum, we use indifferently ${\bm{p}}$ or $p$ when no ambiguity can arise. The (special) relativistic (and Lorentz covariant) integration measure is defined as

[TABLE]

and its associated (special) relativistic Dirac function is defined accordingly as

[TABLE]

such that $\int[{\rm d}p]\underline{\delta}(p-p^{\prime})=1$. Our metric convention follows the standard notation employed in cosmology, which is the opposite of the metric commonly used in particle physics. In the tetrad basis, the metric $g_{\mu\nu}$ reduces to the Minkowski metric

[TABLE]

The Levi-Civita tensor is fully antisymmetric and in the tetrad basis all its components are deduced from the choice

[TABLE]

We identify the time-like vector of a tetrad $e_{0}$ with the velocity $u$ of an observer and its spatial Levi-Civita tensor is obtained from $\epsilon_{ijk}\equiv u^{\mu}\epsilon_{\mu ijk}$, such that $\epsilon_{123}=1$.

I.1.2 Number operator

Creation and annihilation operators, $a^{\dagger}_{r}({\bm{p}})$ and $a_{r}({\bm{p}})$ respectively, where the index $r$ refers to a helicity basis and ${\bm{p}}$ to the particle momentum, are defined for each particle type from its corresponding quantum field. It allows to define a quantum number operator as

[TABLE]

The total occupation operator is then obtained from a sum over all possible momenta of the diagonal part as

[TABLE]

When considering a given quantum state $|\Psi\rangle$, the average of the number operators allows to define a distribution function with helicity indices $f_{rs}$ as

[TABLE]

Hence the total number of particles is given by

[TABLE]

where we introduced the total volume $(2\pi)^{3}\delta^{3}(0)={\cal V}$. In this expression, $f_{rs}(p)$ corresponds exactly to the definition of a classical one-particle distribution function. By construction $N_{rs}$ and $f_{rs}$ are Hermitian, that is

[TABLE]

So far we have not specialized to particles nor antiparticles, not even to a special spin type (fermions or bosons), and this construction is very general. In the next two sections we study separately fermions and bosons, and we show how the distribution function with helicity indices ($f_{rs}({\bm{p}})$) can be decomposed into covariant components.

I.1.3 Adapted orthonormal basis

For a given observer with four-velocity $u^{\alpha}$ which is chosen to be aligned with the time-like tetrad vector $e_{0}$, we define the unit spatial vector of momentum direction by

[TABLE]

In spherical coordinates the momentum direction is given by $\theta,\phi$ and defines a radial unit vector. We then also consider the usual basis in spherical coordinates ${\bm{e}}_{\theta}$ and ${\bm{e}}_{\phi}$, which are purely spatial unit vectors. In tetrad components these are given by

[TABLE]

Let us introduce the helicity vector

[TABLE]

which is a unit vector in the direction of the spatial momentum that is transverse to $p^{\mu}$ in the sense $S^{\mu}p_{\mu}=0$, and is thus spacelike. Since the space of vectors orthogonal to $p^{\mu}$ is three-dimensional, the transverse property is not enough to specify the helicity vector and the definition (23) depends explicitly on the observer which is used to define the spatial part of the momentum. When no ambiguity can arise we write simply $S^{\mu}$. In components the helicity vector is given by

[TABLE]

Geometrically (see Fig. 1), the helicity vector corresponds to the spatial direction unit vector ${{\bm{n}}}$ boosted in its direction by the same boost needed to obtain $p^{\mu}/m$ from $u^{\mu}$.

Finally, we define the polarization basis ${\bm{\epsilon}}_{\pm}(u,p)$

[TABLE]

where again the dependence on $(u,p)$ can be omitted whenever no ambiguity can arise.

The set of vectors $S^{\mu},\epsilon_{\pm}^{\mu}$, and $p^{\mu}/m$ constitute an adapted orthonormal basis to a given observer and a given momentum.

I.1.4 Fermions

The quantum fermion field is

[TABLE]

and satisfies the Dirac equation $({\rm i}\gamma^{\mu}\partial_{\mu}+m)\psi=0$. In this expression, the creation and annihilation operators of the particles ($a_{s},a^{\dagger}_{s}$) and antiparticles (${\bar{a}}_{s},{\bar{a}}^{\dagger}_{s}$) satisfy the anti-commutation rules

[TABLE]

with all other anti-commutators vanishing and where $\delta^{\rm K}$ is the Kronecker function.

We then define the operator in spinor space (beware of position of helicity indices for antiparticles)

[TABLE]

For the sake of clarity we use a notation where components of operators in spinor space and Dirac spinors are explicit and are denoted by indices of the type $\mathfrak{a},\mathfrak{b},\dots$. The plane waves solutions $u_{s,\mathfrak{a}}$ and $v_{s,\mathfrak{a}}$ are the positive and negative frequency solutions and satisfy

[TABLE]

with the standard Dirac slashed notation $\not{p}\equiv p_{\mu}\gamma^{\mu}$ and the Dirac matrices satisfying the algebra $\{\gamma^{\mu},\gamma^{\nu}\}=-2g^{\mu\nu}$. As detailed in appendix A, all spinor space operators can be decomposed on the complete set

[TABLE]

and in particular the operators (31) are decomposed as

[TABLE]

Note also how the indices are in reverse order for antiparticles ($f_{sr}$ instead of $f_{rs}$) echoing a similar placement of indices in Eqs. (31).

The operators of the type (309) take a simple form in the adapted orthonormal basis defined in § I.1.3. Using the helicity basis

[TABLE]

where the right and left chiral parts are defined as

[TABLE]

we obtain

[TABLE]

In particular, we recover when summing on helicities the standard result

[TABLE]

Furthermore, when helicities are different, we obtain the so-called Bouchiat-Michel formulae (Bouchiat and Michel, 1958) [see also Dreiner et al. (2010, App. H.4) or Langenfeld (2007, App. E.3)]. Using the polarization basis (25) we show (Fidler and Pitrou, 2017, App. D7) that it is cast in the compact form

[TABLE]

Let us define111We use the obvious abuse of notation $f_{++}$ for e.g. $f_{+\tfrac{1}{2}\,+\tfrac{1}{2}}$. from the distribution function with helicity indices

[TABLE]

together with

[TABLE]

The functions $(I,Q,U,V)$ are the Stokes parameters. In detail, $I$ is the total intensity, $V$ the circular polarization, and $Q^{\mu}$ is the purely linear polarization vector. ${\cal Q}^{\mu}$ is the total polarization vector, taking into account both circular and linear polarization. By construction the total polarization ${\cal Q}^{\mu}$ is transverse to the momentum (${\cal Q}^{\mu}p_{\mu}=0$). The linear polarization $Q^{\mu}$ is transverse both to the momentum and to the observer velocity $u^{\mu}$, that is it is a purely spatial vector.

The covariant parts are defined from the decomposition

[TABLE]

[TABLE]

One degree of freedom corresponds to the total intensity $I({\bm{p}})$ while the three remaining degrees correspond to the state of polarization and are covariantly contained in a vector ${\cal Q}^{\mu}({\bm{p}})$ because of its transverse property.

The decomposition (45) can be understood from group representations. Indeed, the total polarization vector is a spin-$1$ representation of ${\rm SO}(3)\simeq{\rm SU}(2)$ and the intensity is a spin-[math] representation. When forming the number operator (6), and thus $f_{rs}$, we are building the tensor product of spin-$1/2$ representations and what we have achieved is a decomposition of the reducible representation ${\bm{2}}\otimes{\bm{2}}$ in irreducible components ${\bm{3}}\oplus{\bm{1}}$, where we have denoted ${\bm{1}},{\bm{2}},{\bm{3}}$ the spin-$0,1/2,1$ representations of ${\rm SU}(2)$.

I.1.5 Massless fermions

In the massless limit, the previous decomposition is slightly modified to

[TABLE]

with $\widetilde{\Sigma}^{\mu\nu}$ defined in Eq. (305). Note that the linear polarization $Q^{\mu}$ and the circular polarization $V$ enter separately, and not as a total polarization vector ${\cal Q}^{\mu}$ as is the case for massive fermions. Using $\gamma^{5}\not{Q}\not{p}=\not{p}\gamma^{5}\not{Q}=2p^{\mu}\widetilde{\Sigma}_{\mu\nu}{Q}^{\nu}$ it can also be rewritten as

[TABLE]

In the massless case, the little group of the Lorentz group (Weinberg, 1995) is not ${\rm SO}(3)$ but ${\rm SO}(2)\simeq{\rm U}(1)$. Hence the decomposition in irreducible representations is of the form ${\bm{2}}_{1}\oplus{\bm{1}}\oplus{\bm{1}}$ where the purely linear polarization is in the spin-$1$ representation of ${\rm SO}(2)$ (noted ${\bm{2}}_{1}$) and circular polarization is in the representation ${\bm{1}}$.

Also, it is no longer possible to overtake the particles as they move at the speed of light in any coordinate system. This leads to both, the circular and linear polarizations $V$ and $Q^{\mu}$ to be individually observer independent. More rigorously, linear polarization is described by the coset of

[TABLE]

Indeed, since the polarization basis satisfies $\epsilon_{\pm}^{\mu}p_{\mu}=0$, but we also have $p_{\mu}p^{\mu}=0$, there is a gauge freedom in the definition of the polarization basis. The choice (25) corresponds to the particular choice which is also transverse to the observer velocity ($\epsilon_{\pm}^{\mu}u_{\mu}=0$), which selects unique representatives of polarization vectors. Therefore the polarization vector representative $\epsilon_{\pm}^{\mu}$ are observer dependent, but not the associated cosets $[\epsilon_{\pm}^{\mu}]$. As a consequence $Q^{\mu}$ is observer dependent but not its coset $[Q^{\mu}]$.

Given a representative of the coset, the one associated with a given observer (that is such that it is transverse to that observer velocity) is obtained by projection with a screen projector ${\cal H}_{\mu\nu}(u,p)$, which is abbreviated as ${\cal H}_{\mu\nu}$ when no ambiguity can arise. Using the decomposition of the null momentum into energy and unit direction

[TABLE]

where $E\equiv-u_{\mu}p^{\mu}=p^{0}$, the screen projector is built from the equivalent definitions

[TABLE]

with $k^{\mu}$ is a future directed null vector in the plane spanned by $(u^{\mu},p^{\mu})$ such that $k^{\mu}p_{\mu}=-1$, $k_{\mu}\epsilon_{\pm}^{\mu}=0$. It can be checked that the screen projector satisfies the expected properties ${{\cal H}_{\mu}}^{\sigma}{\cal H}_{\sigma\nu}={\cal H}_{\mu\nu}$ and ${\cal H}_{\mu\nu}p^{\nu}={\cal H}_{\mu\nu}u^{\nu}=0$. If the observer used in the definition is the natural observer associated with the tetrad with which components are taken (that is if $u^{\mu}=e_{o}^{\,\mu}$), the non-vanishing components of the screen projector are only ${\cal H}_{ij}=\delta_{ij}-n_{i}n_{j}$.

For two screen projectors associated with two observers $u$ and $\tilde{u}$ related by a boost

[TABLE]

but for the same momentum $p$, we find that they are related by

[TABLE]

where $\tilde{E}\equiv-\tilde{u}_{\mu}p^{\mu}$. In particular this implies

[TABLE]

Using the screen projector, another definition of the linear polarization coset is that two polarization vectors $Q_{1}^{\mu}$ and $Q_{2}^{\mu}$ describe the same state ($[Q_{1}^{\mu}]=[Q_{2}^{\mu}]$) if

[TABLE]

Note that for a transverse vector ($X_{\mu}p^{\mu}=0$) it is obvious from the decomposition (51) or the transformation rule (53) that

[TABLE]

implying that the definition (55) is unambiguous.

For photons, that is massless bosons, on which we focus in the next section, the structure is exactly similar and arises from the electromagnetic gauge freedom.

I.1.6 Massless bosons

The null mass bosonic vector field of quantum electrodynamics is

[TABLE]

where the creation and annihilation operators satisfy the commutation rule

[TABLE]

If vectors are massive, then the null helicity ($s=0$) must also be considered, see Fidler and Pitrou (2017, App. A).

A covariant distribution tensor is obtained by considering

[TABLE]

and by construction it is transverse to the momentum and the observer’s velocity ($f_{\mu\nu}p^{\nu}=f_{\mu\nu}u^{\mu}=0$). When no ambiguity can arise, we omit the dependence on the observer’s velocity $u^{\mu}$ used in its definition.

We define

[TABLE]

as the usual Stokes parameters222It is sometimes customary in the cosmic microwave background context to define the distribution function as Durrer (2008) $f^{\rm CMB}_{r\,s}\equiv f_{-r\,s}$. Accordingly, the tensor valued function (58) is defined as $f^{\mu\nu}=\sum_{rs}f^{\rm CMB}_{rs}\epsilon_{r}^{\mu}\epsilon_{s}^{\nu}$. With this definition the Stokes parameters are $I\equiv f^{\rm CMB}_{-+}+f^{\rm CMB}_{+-}$, $V\equiv f^{\rm CMB}_{-+}-f^{\rm CMB}_{+-}$ and $Q\pm{\rm i}U\equiv 2f^{\rm CMB}_{\pm\pm}$. corresponding to intensity, circular polarization and linear polarization. For a given observer with four-velocity $u^{\mu}$, we use as in the massless fermion case the spatial momentum direction unit vector ${\bm{n}}$ defined in the decomposition (50). Let us also define the two-dimensional Levi-Civita tensor

[TABLE]

We usually omit the dependence on $(u,p)$ and write simply $\epsilon_{\mu\nu}$. The tensor-valued distribution function is decomposed as

[TABLE]

where the screen projector is defined exactly as for massless fermions in Eqs. (51). The distribution tensor is doubly transverse, that is transverse to the momentum $p^{\mu}$ and also to the observer velocity $u^{\mu}$. ${\cal P}_{\mu\nu}$ is the linear polarization tensor and it is doubly transverse and traceless (it satisfies ${\cal P}_{\mu\nu}u^{\mu}={\cal P}_{\mu\nu}p^{\mu}={{\cal P}_{\mu}}^{\mu}=0$). It is defined as

[TABLE]

and its dependence on $u$ is often omitted. It can be extracted thanks to the transverse traceless projector

[TABLE]

In the ${\bm{e}}_{\theta},{\bm{e}}_{\phi}$ basis the components of the distribution tensor (61) form a $2\times 2$ Hermitian matrix Hu and White (1997); Tsagas et al. (2008); Durrer (2008)

[TABLE]

whereas in the ${\bm{\epsilon}}_{-},{\bm{\epsilon}}_{+}$ basis we obtain the Hermitian matrix

[TABLE]

For massless bosons, the structure of the decomposition can also be understood exactly like in the discussion following Eq. (48) for massless fermions. The difference is that for massless bosons, we decompose ${\bm{2}}_{1}\otimes{\bm{2}}_{1}$ into ${\bm{2}}_{2}\oplus{\bm{1}}\oplus{\bm{1}}$, where ${\bm{2}}_{1}$ (resp. ${\bm{2}}_{2}$) is the spin-$1$ (resp. spin-$2$) representation of ${\rm SO}(2)$.

As in the case of massless fermions, the definition (58) and the decomposition (61) of the distribution tensor is observer dependent, for exactly the same reasons that the polarization vectors are defined up to factors of $p^{\mu}$. Hence, we should rather consider the coset $[f_{\mu\nu}]$. Two polarization states $f^{1}_{\mu\nu}$ and $f^{2}_{\mu\nu}$ are in the same coset if

[TABLE]

In particular, the linear polarization parts ${\cal P}^{1}_{\mu\nu}$ and ${\cal P}^{2}_{\mu\nu}$ are equivalent if

[TABLE]

and one should rather consider the coset $[{\cal P}_{\mu\nu}]$ of linear polarization. With arguments similar to Eq. (56), this definition of equivalence (and its associated cosets) is observer independent.

I.2 Multipolar decomposition

I.2.1 Fermions

The intensity $I({\bm{p}})$ is easily decomposed into spherical harmonics. Indeed, once an observer choice is made, that is its four-velocity $u^{\mu}$ is identified with the time-like vector of the tetrad $[e_{0}]^{\mu}$, we can define the spatial momentum ${\bm{p}}$ and its direction unit vector ${{\bm{n}}}$ (see § I.1.3). We then perform the usual spherical harmonics decomposition

[TABLE]

Using Eq. (325) the multipoles are extracted as $I_{\ell m}=\{Y_{\ell m}|I\}$.

Alternatively one could use a decomposition based on symmetric trace-free (STF) tensors $I_{i_{1}\dots i_{\ell}}$ which is equivalent  Thorne (1980); Blanchet and Damour (1986); Pitrou (2009a, b)

[TABLE]

where we use the tools and notation summarized in appendix D. From Eq. (325) the STF tensors are extracted as

[TABLE]

The relation between both expansions is obtained from Eqs. (D.2) as

[TABLE]

For the polarization vector ${\cal Q}^{\mu}$ of fermions defined in Eq. (44), we have to pay attention to the transformation properties when performing a spatial rotation of the coordinate system around the direction of ${{\bm{n}}}$. The ordinary spherical harmonics, when evaluated at ${{\bm{n}}}$ do not transform under this rotation and are thus not suitable to decompose objects which have a non-trivial transformation under this rotation. The polarization vector $\cal Q$ transforms as an ordinary 4-vector (we have shown that it is observer-independent). However, this is not the case for the observer-dependent vectors and distribution functions used to build $\cal Q$. The vector in direction of the spatial momentum $S^{\mu}$ is invariant under this particular rotation as it points in the direction ${{\bm{n}}}$. Employing the observer-independence of $\cal Q$ which is discussed in the next section, we therefore conclude that $V$ must be invariant under this rotation and may be decomposed into ordinary spherical harmonics.

[TABLE]

Again an expansion in STF tensors of the type (69) is possible and is obtained by relations exactly similar to Eqs. (71).

The polarization vectors ${\epsilon}_{\pm}({{\bm{n}}})$ however transform with an additional spin $\mp 1$ complex rotation. To generate an observer-independent $\cal Q$ the corresponding $Q_{\pm}$ must transform with the opposite spin and they are decomposed into spin-weighted spherical harmonics $Y^{s}_{lm}$ Goldberg et al. (1967) as

[TABLE]

Note that this discussion only concerns the observer dependence under a specific spatial rotation and that due to the definition of helicity an additional dependence mixing $V$ and $Q_{\pm}$ exists for more general rotations and boosts. ${\cal E}$ and ${\cal B}$ modes multipoles can be defined from

[TABLE]

The ${\cal E}_{\ell m}$ have even parity (they get a factor $(-1)^{\ell}$ under parity transformation) whereas the ${\cal B}_{\ell m}$ have odd parity (they get a factor $(-1)^{\ell+1}$ under parity transformation) since spin-weighted spherical harmonics transform as $Y^{s}_{\ell m}\to(-1)^{\ell}Y^{-s}_{\ell m}$ and the polarization basis transforms as ${\bm{\epsilon}}_{\pm}\to-{\bm{\epsilon}}_{\mp}$. Equivalently since ${\bm{Q}}({\bm{p}})$ is a vector field on the unit sphere in momentum space, it can be decomposed as the gradient and the curl of two scalar functions as

[TABLE]

where $D_{i}$ is the covariant derivative on the unit sphere and $\epsilon_{ij}\equiv{\epsilon_{ijk}}n^{k}$ is the Levi-Civita tensor on the unit sphere already defined in Eq. (60). Decomposing the scalar functions ${E}$ and ${B}$ in multipoles ${E}_{\ell m}$ and ${B}_{\ell m}$ as in the expansion (68) and using Durrer (2008)

[TABLE]

the two possible definitions for the $E$ and $B$ modes multipoles are related by ${\cal E}_{\ell m}=\sqrt{\ell(\ell+1)/2}{E}_{\ell m}$ and ${\cal B}_{\ell m}=\sqrt{\ell(\ell+1)/2}{B}_{\ell m}$. Again a similar expansion can be obtained by using symmetric trace-free tensors to expand the scalar functions ${E}$ and ${B}$ directly in Eq. (75).

I.2.2 Massless bosons

The decomposition of intensity and circular polarization is performed with spherical harmonics as in Eqs. (68) and (72) or with STF tensors as detailed in § I.2.1 for fermions. However, the linear polarization part ${\cal P}_{\mu\nu}$ must be decomposed in spin-$2$ spherical harmonics. We decompose polarization as

[TABLE]

and the angular decomposition is

[TABLE]

Note that the factor $2$ in Eq. (77) is purely conventional. ${\cal E}$ and ${\cal B}$ modes are defined by

[TABLE]

Equivalently linear polarization can be decomposed with two potentials on the unit sphere in momentum space as (Tsagas et al., 2008, Eq. 4.3.8)

[TABLE]

and the associated multipoles ${E}_{\ell m}$ and ${B}_{\ell m}$ can be related to the ${\cal E}_{\ell m}$ and ${\cal B}_{\ell m}$ by some factors. Instead, if we use an expansion of ${E}$ and ${B}$ in STF tensors of the type (69), we can decompose ${\cal P}_{ij}$ with them. However, it is customary to remove the $\ell(\ell-1)$ factors brought by the covariant derivatives $D_{i}$ and use the expansion (Dautcourt and Rose, 1978) [see also Tsagas et al. (2008, Eq. 4.3.9) or Pitrou (2009a, Eq. 1.33)]

[TABLE]

The exponent ${\cal T}$ indicates that free indices are to be projected on the transverse traceless part with the operator (63). From the definition (77) of ${\cal P}_{\pm\pm}$, and using the notation (337), this expansion is equivalent to

[TABLE]

The STF tensors of the decomposition (81) are extracted thanks to (Tsagas et al., 2008; Pitrou, 2009a)

[TABLE]

where

[TABLE]

If we now associate to these STF tensors $E_{K_{\ell}}$ and $B_{K_{\ell}}$ the ${E}_{\ell m}$ and ${B}_{\ell m}$, using a relation of the type (71), these are related to the ${\cal E}_{\ell m}$ and ${\cal B}_{\ell m}$ defined in (79) by

[TABLE]

This is obtained using Eqs. (340) in Eq. (82) and comparing with Eq. (78).

I.3 Observer independence

In this section, we detail the transformation property of the distribution function under a general Lorentz transformation $\Lambda\in{\rm SO}^{+}(1,3)$. It is more appropriate to take the passive point of view and consider a transformed tetrad basis related to the initial one by

[TABLE]

The new observer’s velocity $\tilde{{\bm{u}}}$ is identified with the time-like vector of the new tetrad $\tilde{{\bm{e}}}_{\tilde{0}}$. That is we take the point of view that when considering a change of frame we also consider the associated change of observer, such a that any observer is not moving in its own frame. In that sense, the observer’s velocity is not observer independent.

The new components of the momentum $p^{\tilde{\mu}}\equiv p\cdot\tilde{{\bm{e}}}^{\tilde{\mu}}$ are related to the previous ones $p^{\mu}\equiv p\cdot{\bm{e}}^{\mu}$ by

[TABLE]

and we abbreviate $p^{\tilde{\mu}}$ as $\tilde{p}$.

I.3.1 Massive fermions

In Fidler and Pitrou (2017), we showed that the spinor valued operator transforms under the Lorentz transformation $\Lambda\in{\rm SO}^{+}(1,3)$ defined by Eqs. (86) as

[TABLE]

where $D(\Lambda)$ is the spinor-space representation of $\Lambda$. Using the property for Dirac matrices

[TABLE]

it implies that the covariant components for massive fermions transform as

[TABLE]

This means that they transform exactly as a scalar and vector field, and they are therefore observer independent.

The observer independence is important as it allows to build a statistical description of the fluid without the need to specify an observer first. This is particularly useful for deriving simple transport equations in general relativity.

The scalar $I$ describes the total intensity of the field and is observer independent since the local number of particles is identical for each observer. The information of the polarization of the fluid is contained in the observer independent vector ${\cal Q}^{\mu}$.

On the other hand the parameters $V$ and $Q_{\pm}$, describing individually the circular and linear polarizations are not observer independent. The circular polarization $V$, for example, changes if the observer is boosted and overtakes the momentum considered. We have defined

[TABLE]

where ${\cal Q}^{\mu}$ combines multiple observer dependent quantities into one observer independent vector. In the example of the observer overtaking a particle momentum, we change all left-helical $-\tfrac{1}{2}$ states into right-helical $+\tfrac{1}{2}$ states. This means that the boosted observer will find $\widetilde{V}=-V$. At the same time the vector $S^{\mu}$ is also observer dependent and the new observer will define the spatial momentum of the particles with the opposite sign. Therefore the combination $VS^{\mu}$ is invariant under this boost. At the same time the off-diagonal distributions are swapped: $\widetilde{f}_{+-}=f_{-+}$. However these are combined with the polarization vectors $\epsilon_{\pm}$ to form $Q^{\mu}$, which are also interchanged for the new observer, leading to $Q^{\mu}$ being invariant.

In a more general case $Q^{\mu}$ and $V$ cannot be disentangled in an observer independent manner and there always exists a subset of observers, all related by boosts along the momentum direction and rotations around the momentum direction, that will perceive the field to be entirely circularly polarised without any linear polarization. For this reason we will work with the observer independent polarization vector ${\cal Q}^{\mu}$ and only refer to the circular and linear polarizations when we have specified an observer. Only in the case of massless fermions, considered in § (I.1.5), the linear and circular polarization can be disentangled, and are observer independent, the latter in the sense of the polarization coset (49) as detailed in the next section.

I.3.2 Massless fermions

Using the decomposition (47) for massless fermions, we deduce that the covariant components transform as

[TABLE]

The screen projector [see Def. (51)] associated with the new observer and the new momentum components,

[TABLE]

(with $\tilde{E}\equiv-\tilde{u}_{\tilde{\mu}}p^{\tilde{\mu}}=p^{\tilde{0}}$) ensures that the linear polarization remains spatial for the new observer. Hence in the massless case, the linear polarization part is not strictly observer independent, but since this dependence introduced by the screen projector is there only as the result of a choice to remove a non physical degree of freedom, we can still conclude that in that sense the covariant components are observer independent. More rigorously, it is the coset of linear polarization [see definition (49)] which is observer independent and only the special choice of its representative element is observer dependent. Hence we should rather write the transformation rule of linear polarization cosets which is

[TABLE]

for which the observer independence is manifest.

I.3.3 Massless bosons

For massless bosons, the tensor-valued distribution function transforms as

[TABLE]

with the definition $(\widetilde{\cal H}\Lambda)^{\tilde{\mu}}_{\,\,\sigma}\equiv{{\widetilde{\cal H}}^{\tilde{\mu}}}_{\,\,\tilde{\nu}}\Lambda^{\tilde{\nu}}_{\,\,\sigma}$. Since the screen projector satisfies333This is exactly Eq. (53) but expressed with components associated to different tetrads.

[TABLE]

and the two-dimensional Levi-Civita tensor (60) satisfies a similar property, then we deduce that the covariant components transform as

[TABLE]

where $\widetilde{\cal T}$ is the transverse-traceless projector associated with $\widetilde{\cal H}_{\tilde{\mu}\tilde{\nu}}$, using the definitions (63) and (93). As in the case of massless fermions, it is the coset of linear polarization [see Eq. (66)] which is observer independent, and only the special choice of its representative element is observer dependent. Hence we should rather write the transformation rule as $\left[\widetilde{f}^{\tilde{\mu}\tilde{\nu}}(\tilde{p})\right]=\left[{\Lambda^{\tilde{\mu}}}_{\sigma}{\Lambda^{\tilde{\nu}}}_{\tau}{f}^{\sigma\tau}(p)\right]$, with the cosets defined by the equivalence relation (66), and for which the observer independence is manifest.

I.3.4 Relation to abstract tensor indices

Since we have shown that all components of the vectors or tensors associated to fermions and photons have the expected transformation properties, we could decide to work with abstract indices as in Challinor et al. (2000); Challinor (2000a); Tsagas et al. (2008) instead of working with indices referring to a particular tetrad. In most cases this reinterpretation is straightforward. However both approaches differ when it comes to expressing in practice the transformation of the STF tensors presented in § I.2 for the angular decomposition of the distribution functions. With abstract indices, projectors still appear in the transformation rules, as e.g. in Eqs. (4.3.31-4.3.33) of Tsagas et al. (2008), whereas with indices referring to components in tetrads the transformations relate STF tensors which all are purely spatial in their associated tetrad, that is we relate only spatial indices, as in Eqs. (1.56-1.58) of Pitrou (2009a). However, this subtlety only shows when the transformation of the multipoles is performed at least at second order in the boost velocity. In the remainder of this article, we use a method where no change of frame is needed, hence we do not detail any further the procedure to obtain the multipoles transformation rules. More details can be found in Pitrou (2009a).

II Boltzmann equation

II.1 Liouville equation in curved space-time

The previous construction was restricted to a homogeneous system, hence the functions appearing ($I$ and $Q^{\mu}$) depended only on $(t,p)$. In order to describe a gas of particles classically, one must assume that this construction is in fact valid only locally. That is we assume that there is a mesoscopic scale and that our previous construction was restricted to scales much smaller. The functions, which were dependent on $(t,p)$ must depend now on $(t,{\bm{x}},p)$. In order to derive a Liouville equation in curved space-time which describes the evolution of the covariant components, we must also distinguish between the massive and the massless cases.

II.1.1 Massive fermions

In the previous sections we have shown that $I$ and ${\cal Q}$ are observer independent. In addition, in the local Minkowski frame, they are also parallel transported in the absence of collisions. The helicity of particles does not change in free propagation and, considering that the momentum $p^{\mu}$ is conserved, the vectors $\epsilon_{\pm}^{\mu}$ and $S^{\mu}$ used to build the quantities $I$ and ${\cal Q}$ remain unchanged. Hence, in the local Minkowski space we obtain the equations of motion

[TABLE]

From the point of view of general relativity, these equations are only valid locally and neglect entirely the impact of the relativistic space-time. The intensity $I$ describes the total number of particles. The conservation of $I$ in the absence of collisions in Eq. (98) is equivalent to mass or particle number conservation. The geometrical impact of general relativity does not change the number of particles and we may generalise the equation of motion by requiring the conservation of $I$ along a full geodesic

[TABLE]

where $\frac{D}{D\lambda}$ is the derivative along the particle trajectory parameterized by $\lambda$.

The vector ${\cal Q}$ is parallel transported in the local space-time and describes the polarization of particles in an observer-independent way. Again, the geometrical nature of general relativity does not change the polarization of particles and we require that ${\cal Q}$ is parallel transported along the non-trivial trajectory of the particles. Note that the observer dependent linear and circular polarization may change non-trivially during the transport and require a specification of the dynamics of the observer.

Using the observer-independence, we are able to uniquely define the vector ${\cal Q}$ on our full space-time by employing the tetrads

[TABLE]

where we remind that the index $\mu$ is a tetrad component index, but the index $\alpha$ is a general coordinate index. Assuming parallel transport, we obtain the equation of motion

[TABLE]

Note that ${\cal Q}$ is by definition orthogonal to the momentum. This property is automatically conserved in the relativistic evolution as both the momentum and $\cal{Q}$ are parallel-transported along the geodesic of a free particle.

II.1.2 Massless fermions

In the massless case, linear polarization and circular polarization must be considered separately. Circular polarization $V$ is transported exactly like the intensity $I$ in Eqs. (99) because the direction of the helicity vector is identical to the momentum and therefore parallel-transported. However the linear polarization vector (considered in general coordinates with $Q^{\alpha}=Q^{\mu}[e_{\mu}]^{\alpha}$) cannot be parallel transported because it is transverse to both the momentum and the observer velocity $u^{\alpha}$, and the latter is not (necessarily) parallel transported. However, in the process of free streaming, any variation of $Q^{\alpha}$ in the direction of the momentum is not physical. Hence this unphysical degree of freedom must be eliminated by an appropriate projection so as to obtain an unambiguous equation for parallel transport. To that purpose, we use the screen projector (51) in general coordinates and write

[TABLE]

The transport of linear polarization in the massless case is the same as the transport of the full polarization vector in the massive case [Eq. (101)], up to an additional screen projection which ensures that the double transverse property holds. It can be equivalently formulated by saying that the coset $[Q^{\alpha}]$ is parallel transported, that is

[TABLE]

II.1.3 Massless bosons

The parallel transport of linear polarization for massless bosons, that is photons is very similar to massless fermions, except that instead of projecting a vector we must project a tensor (Challinor, 2000b, a; Tsagas et al., 2008; Pitrou, 2009a, b). We define polarization on the full spacetime as in Eq. (100), that is

[TABLE]

Similarly, the non-physical degree of freedom must be projected and the evolution of linear polarization is dictated by

[TABLE]

Note that for massless bosons, we need not postulate this equation as it is obtained from the eikonal approximation of electromagnetism, see e.g. Fleury (2015) for a detailed account on the procedure.

II.2 Quantum evolution in the interaction picture

So far we have discussed the free propagation of fermions or bosons. When in addition considering collisions, we will employ a separation of scales. We assume that the relativistic evolution is dominant on macroscopic scales, while individual collisions act on microscopic scales. We therefore may compute the collision term in the local tangent space corresponding to special relativity. Then averaging over the local Minkowski space-time of the observer we will provide an effective collision term for the relativistic evolution of the distribution functions.

We therefore introduce three separate scales, the microscopic scale of individual interactions, typically the Compton timescale of interacting particles. Then a mesoscopic scale over which we average the individual collisions, define our local distribution functions and describe the impact of the collisions on the averaged fluid. Finally, the macroscopic scale on which particles free stream on general relativistic geodesics. This separation of scales is illustrated in Fig. 2.

We begin with the description of collisions in the local frame of our observer. The full Hamiltonian $H$ can be separated into a free part $H_{0}$ and an interaction part $H_{\rm I}$. We employ the Heisenberg picture in which the states are time-independent. The time evolution of our distribution function is given by (omitting to specify the momentum dependence of $f_{rs}$ and $N_{rs}$ for simplicity)

[TABLE]

We find a differential equation for the operator $N_{rs}$ and are able to write an approximate solution as closed integration if we restrict ourselves to a given order in the interaction Hamiltonian. The details can be found in Fidler and Pitrou (2017) and are also summarized in appendix B. They require to separate between the microscopic scales of the quantum collisions and the macroscopic scales of the classical Boltzmann transport description.

Eventually, defining the collision term as

[TABLE]

the evolution of $f_{rs}$ is then ruled by the Boltzmann equation

[TABLE]

In the case of fermions, a spinor space operator associated with this collision term is obtained by contraction with $u_{s^{\prime}}(p)\bar{u}_{s}(p)$ (or $v_{s}(p)\bar{v}_{s^{\prime}}(p)$ for antiparticles) as in Eq. (34), and we define

[TABLE]

The covariant parts of this spinor space collision operator, $I_{C}(p)$ and ${\cal Q}_{C}^{\mu}(p)$ are obtained exactly like in Eq. (45). In the massless case the covariant parts are $I_{C}(p)$, $V_{C}(p)$ and ${Q}_{C}^{\mu}(p)$ and are obtained as in Eq. (47). For massless bosons (photons), a tensor valued collision function is built as in Eq. (58).

The classical Boltzmann equation is obtained when considering that this derivation, which has been made for a homogeneous system is in fact valid locally. That is in the derivation we assumed that the distribution function depends on time and momentum only $f_{rs}(t,p)$, but we now assume that it also depends on the position and employ $f_{rs}(t,{\bm{x}},p)$. This amounts to considering that under the mesoscopic length scale the system can be considered as homogeneous (see also Fig. 2), such that the volume integral in the Hamiltonian $H_{I}=\int{\rm d}^{3}{\bm{x}}{\cal H}_{I}$ can be extended to infinity in the computation of the local collision term $C[f_{ss^{\prime}}(t,{\bm{x}},p)]$. Expressed in terms of spinor valued or tensor valued operators the classical Boltzmann equation reads

[TABLE]

Finally, in order to include this collision term in the right hand side of the Liouville equation in curved space-time discussed in § II.1, one must multiply it by ${\rm d}t/{\rm d}\lambda=p^{0}=E$. This converts the collision term, seen as a rate of change of the distribution function per unit of proper time of the observer in the tetrad frame, to a collision term which is a change of the distribution function per unit of the affine parameter $\lambda$. In practice the Boltzmann equation obtained needs to be converted again to an equation giving the change of the distribution function per unit of a generalized time coordinate, and this final step requires a specific form of the metric.

II.3 Molecular chaos

In principle, when considering an interacting system, the one-particle distribution function is not enough to describe it statistically, because $n$-particle correlation functions are generated by collisions. In order to obtain a description only in terms of a one-particle distribution function, we must assume that the connected part of the $n$-particle functions vanishes and thus that $n$-particle functions are expressed only in terms of one-particle functions, corresponding to the assumption of molecular chaos. We review how this assumption is implemented in this section.

Let us introduce a multi-index notation which encodes both the helicity index and the momentum, and which consists in using $s$ for $(s,p)$ or $s^{\prime}$ for $(s^{\prime},p^{\prime})$. With this notation we write for instance $a_{s^{\prime}}$ instead of $a_{s^{\prime}}(p^{\prime})$. We also introduce a generalized delta function on both helicities and momenta which is

[TABLE]

In particular the number operator (6) is noted

[TABLE]

For fermions, we get from anticommuting rules

[TABLE]

which defines the Pauli blocking operator $\widehat{N}_{ss^{\prime}}$. Similarly for bosons, we get from commutation rules

[TABLE]

which defines the stimulated emission operator $\widehat{N}_{ss^{\prime}}$. The $n$-particle number operators for species $a$ are defined as

[TABLE]

Under the molecular chaos assumption, their expectation value for fermions in a quantum state is related to the expectation value of the number operator as444For bosons we remove all minus signs, that is the factor $(-1)^{n+1}\epsilon(\sigma)$ so we would for instance get $\langle N^{(2)}_{r_{1}r_{2}s_{1}s_{2}}\rangle=\langle N^{(1)}_{r_{1}s_{1}}\rangle\langle N^{(1)}_{r_{2}s_{2}}\rangle+\langle N^{(1)}_{r_{1}s_{2}}\rangle\langle N^{(1)}_{r_{2}s_{1}}\rangle$.

[TABLE]

where the sum is on the group of permutation $S_{n}$ and $\epsilon(\sigma)$ is the signature of the permutation. This approximation is exactly similar to the Boltzmann approximation of the BBGKY hierarchy (Volpe, 2015). In practice this assumption of molecular chaos is used to obtain the following property for the expectation in a quantum state of a product of one-particle number operators

[TABLE]

[TABLE]

where $\widehat{N}_{rs}$ is either the Pauli blocking operator (for fermions) or the stimulated emission operator (for bosons) defined by Eqs. (112) and (113). The expectation value of a product of one-particle number operators is the sum of products of expectation values of all possible pairings between creation and annihilation operators. For each pair, if the indices $(r_{a},s_{b})$ correspond to operators which were initially in the creation-annihilation order, that is with $a\leq b$ (resp. annihilation-creation order, that is $a>b$) we use $N_{r_{a}s_{b}}$ (resp. $\widehat{N}_{r_{a}s_{b}}$). For instance the expectation value for a product of two one-particle number operators is simply

[TABLE]

Finally, we also assume that species are uncorrelated such that the expectation value for operators of various species is the product of expectation values of the operators of each species. For instance for two species $a$ and $b$ we assume $\langle a_{r}^{\dagger}a_{s}b_{p}^{\dagger}b_{q}\rangle=\langle a_{r}^{\dagger}a_{s}\rangle\langle b_{p}^{\dagger}b_{q}\rangle$.

Weak interactions

The Fermi theory of weak interactions is a contact interaction between four fermions. All reactions in that approximation are of the type

[TABLE]

with other reactions involving antiparticles ($\bar{a},\bar{b},\bar{c},\bar{d}$) deduced from charge conjugation or crossing symmetry. In the next section we derive the general collision term for general weak currents and apply it to the case of neutrino interactions which is relevant for the early universe. In § IV we apply it to the case of neutron-proton conversions by weak interactions which controls the primordial Helium abundance.

III General collision term

III.1 Fermi theory of weak interactions

All weak interaction take the form of current-current interactions Nachtmann and Halzen (1991) at low energy (low compared to the $W^{\pm}$ and $Z$ masses), that is they are given by

[TABLE]

where $G_{F}\simeq 1.1663787\,\times\,10^{-5}\,{\rm GeV}^{-2}$ is the Fermi constant of weak interactions.

III.1.1 Neutral currents

Neutral currents describe the exchange of $Z$ bosons and as these are not charged they mediate elastic scatterings that do not alter the involved types of particles and only transfer momentum, spin and energy.

The neutral current is simply the sum of the neutral currents of all particles undergoing weak interactions

[TABLE]

For neutrinos, the neutral current couples only the left chiralities and, noting ${\bm{\nu}}$ the neutrino quantum field, it is simply

[TABLE]

with similar expressions for other flavors. However for electrons and (similarly pions and taus) the neutral currents must be further decomposed into left and right chiral interactions as

[TABLE]

where we noted ${\bm{e}}$ the electronic quantum field. The chiral coupling constants are for electrons

[TABLE]

with $\theta_{W}$ the Weinberg angle ($\sin^{2}\theta_{W}\simeq 0.23$).

III.1.2 Charged currents

Opposed to the neutral currents, the charged currents describe the exchange of charged $W$-bosons and therefore are inelastic. The structure of the charged current is more complex since it couples eigenmass states of different flavors, thanks to the Cabbibo-Kobayashi-Maskawa (CKM) matrix for quarks or the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) matrix for massive neutrinos. We ignore these complications for the examples that we shall consider and employ effective charged currents for the neutron/proton pair which is involved in beta decays and related processes, and the charged currents of the first two lepton flavors, that is of the electron/neutrino and muon/muon neutrino pairs. We use

[TABLE]

where $V_{ud}=0.97420(20)$ is a Cabbibo-Kobayashi-Maskawa (CKM) angle Patrignani and Particle Data Group (2016 and 2017 update).

The charged currents for electron/neutrino and muon/muon neutrino pairs are coupling only the left chiralities

[TABLE]

However, due to internal QCD effects, the coupling in the neutron/proton pair is not purely left chiral. The deviation from left chirality of the coupling is parameterized by the $g_{A}$ parameter whose measured value is approximately $1.2723(23)$ (Patrignani and Particle Data Group, 2016 and 2017 update) and the corresponding charged current reads

[TABLE]

When considering the cumulative effect of neutral currents and charged currents, we can use the Fierz identities which for anticommuting fields give Sarantakos et al. (1983); Sigl and Raffelt (1993)

[TABLE]

This means that the effect of multiple charged currents can be replaced by equivalent neutral currents. In the collision term we may therefore replace the charged currents by modifying the neutral chiral coupling factors (121), yielding

[TABLE]

III.2 Two-body processes

III.2.1 Notation

We use the compact notation introduced in § II.3 that we adapt to also account for the fact that we have several different species in the reaction (116). We introduce

[TABLE]

These multi-indices contain all information characterising one single particle (its momentum and helicity). We will typically label ingoing states as unprimed and outgoing states with primed indices. For species $a$ we employ the multi-index $\alpha$ and similarly for species $b$ (resp. $c$ and $d$) we use the multi-indices $\beta$ (resp. $\gamma$ and $\delta$). The plane wave solutions are written in a compact form in this notation. For instance for the species $a$ we write $u_{\alpha}\equiv u_{s_{\alpha}}(p_{\alpha})$ and $v_{\alpha}\equiv v_{s_{\alpha}}(p_{\alpha})$. Furthermore this allows to write a compact relativistic Dirac delta function which acts both on helicities and momenta as

[TABLE]

We denote the number operator associated with species $a$ as

[TABLE]

We also define the Pauli blocking operator

[TABLE]

The expectation value of these operators is denoted as

[TABLE]

where we introduce the short-hand notation ${\cal A}_{\alpha\alpha^{\prime}}(p_{\alpha})={\cal A}_{s_{\alpha}s^{\prime}_{\alpha}}(p_{\alpha})$. We recall that this quantity is exactly the one-particle distribution function associated with species $a$ [see Def. 8]. Note that for the Pauli blocking factor, $\widehat{\cal A}_{\alpha\alpha^{\prime}}(p_{\alpha})$ is a shorthand notation for $\delta^{\rm K}_{s_{\alpha}s^{\prime}_{\alpha}}-{\cal A}_{s_{\alpha}s^{\prime}_{\alpha}}(p_{\alpha})$. We associate to the one-particle distribution function (resp. the Pauli blocking function) a spinor valued operator following the procedure (34) that we note ${A_{\mathfrak{a}}}^{\mathfrak{b}}$ (resp. ${\widehat{A}_{\mathfrak{a}}}^{\,\,\mathfrak{b}}$) in component notation or simply ${\bm{A}}$ (resp. $\widehat{{\bm{A}}}$) in operator notation. Having defined for species $a$ the number operator $A_{\alpha\alpha^{\prime}}$, the distribution function ${\cal A}_{\alpha\alpha^{\prime}}$ and the spinor-valued (observer-independent) operator ${\bm{A}}$, we proceed identically for species $b$ (resp. $c$ , $d$) and we use $B_{\beta\beta^{\prime}}$, ${\cal B}_{\beta\beta^{\prime}}$ and ${\bm{B}}$ (resp. $C_{\gamma\gamma^{\prime}}$, ${\cal C}_{\gamma\gamma}^{\prime}$ and ${\bm{C}}$, $D_{\delta\delta^{\prime}}$, ${\cal D}_{\delta\delta^{\prime}}$ and ${\bm{D}}$), and associated hatted notations for Pauli blocking factors. Furthermore, for the antiparticles species $\bar{a},\bar{b},\bar{c},\bar{d}$ related to the species $a,b,c,d$, we use barred notation for number operators (e.g. $\overline{A}_{\alpha\alpha^{\prime}}$), distribution function (e.g. $\overline{\cal A}_{\alpha\alpha^{\prime}}$) and spinor valued operators (e.g. $\overline{{\bm{A}}}$), along with their hatted versions for Pauli blocking terms. Finally we define the collision term as in Eq. (107), that is

[TABLE]

such that the quantum Boltzmann equation (315) for species $a$ is written as (when neglecting forward scattering)

[TABLE]

III.2.2 Collision term structure

Following the previous discussion, our goal is to compute the collision term $C[{\cal A}_{ss^{\prime}}(p)]$ corresponding to the reaction (116) due to weak interactions. It is mediated by an Hamiltonian density of the form

[TABLE]

where, depending on the interaction, the same species may be represented by multiple indices. The chiral contributions of these currents are parameterized by $\epsilon^{ac}_{\pm}$ and $\epsilon^{bd}_{\pm}$ as

[TABLE]

with the notation

[TABLE]

The interaction Hamiltonian associated to the Hamiltonian density (134) is explicitly given by

[TABLE]

where we used the scattering operator for this reaction

[TABLE]

The $M$ matrices are defined with the multi-index notation (127)

[TABLE]

and for weak interactions they are of the general form

[TABLE]

To compute the collision term we first need to compute the operator $[{\cal M},[{\cal M},A_{ss^{\prime}}]]$. Using the commutation rules of Appendix B in Fidler and Pitrou (2017), and using the molecular chaos assumption described in § II.3, we get

[TABLE]

We now employ this result in Eqs. (III.2.2) and (132). We integrate a total of five momentum integrals (each one being itself three-dimensional in momentum space) using the Dirac distributions. Of these, four Dirac functions are contained in the expectation values of the number operators associated to the four species, and there is an extra Dirac function ($\underline{\delta}_{s\alpha^{\prime}}$ or $\underline{\delta}_{\alpha s^{\prime}}$ in Eq. (142)) from the collision term ensuring local energy and momentum conservation. Eventually, taking the expectation in the quantum state, we get

[TABLE]

with the integration on momenta

[TABLE]

We note that:

• •

The collision term is made of two types of terms. The first terms on the second and the third line of Eq. (143) correspond to scattering out processes, that is collisions which due to the minus sign deplete the distribution function associated with species $a$ and they correspond to $a+b\to c+d$. The second term on the second and third line correspond conversely to scattering in processes, which increase the distribution function of species $a$, and they are due to the reaction $c+d\to a+b$.

• •

For scattering out processes, the collision term is proportional to the distribution function of the initial states (species $a$ and $b$), but also to the Pauli blocking function of the final states (species $c$ and $d$), and the reverse is true for the scattering in processes.

• •

The distribution functions are Hermitian, that is ${\cal A}^{\star}_{ss^{\prime}}(p)={\cal A}_{s^{\prime}s}(p)$ as in Eq. (10). Let us now consider $C[{\cal A}_{ss^{\prime}}(p)]^{\star}$. Given the Hermiticity of the distribution functions and thus of the Pauli blocking functions, with a simple renaming of all primed indices as unprimed indices (and also of unprimed indices as primed indices), it is straightforward to show that this is equal to $C[{\cal A}_{s^{\prime}s}(p)]$, hence the collision term is also Hermitian as expected.

• •

In the previous computation when checking the Hermiticity, the second and third line of Eq. (143) are interchanged. Terms of the second line are proportional to $\delta^{\rm K}_{s\alpha^{\prime}}$ and correspond physically to the scattering of the helicity index $s^{\prime}$, and conversely in the third line the terms are proportional to $\delta^{\rm K}_{\alpha s^{\prime}}$ and it corresponds to the scattering of the helicity index $s$. Hence we see that the collision term possesses four terms corresponding to the in/out contributions and the $s/s^{\prime}$ contributions.

• •

Finally even though we computed the collision term for a homogenous system in a Minkowski space-time, the total volume, which appears as $\delta^{(3)}(0)$, drops out from both the left and the right hand side of Eq. (133). Hence, as argued before Eq. (110), we can consider that this collision term is valid locally, allowing us to consider in a classical macroscopic description that all distribution functions should be considered with a dependence on the point of space-time. We started a computation with total number of particles in a quantum system, but we end up using it with number densities of particles, considering that the collisions are point-like.

The procedure to follow is now transparent. The helicity indices of the distribution functions (or the related Pauli blocking functions) are contracted with the plane waves solutions contained in the $M$ matrices. From Eqs. (31) this is exactly what is needed to build the spinor space operators related to each species. Since only the indices $s$ and $s^{\prime}$ remain uncontracted in Eq. (143), we contract them with $u_{s^{\prime}}(p)\bar{u}_{s}(p)$ (or $v_{s}(p)\bar{v}_{s^{\prime}}(p)$ for antiparticles) so as to form a spinor space collision operator $C[{\bm{A}}(p)]$ as specified in the definition (109).

Note that the contraction of $\delta^{\rm K}_{ss^{\prime}}$ with $u_{s^{\prime}}(p)\bar{u}_{s}(p)$ or $v_{s}(p)\bar{v}_{s^{\prime}}(p)$ gives simply

[TABLE]

with the notation (46), as can be seen from Eqs. (41). We finally obtain the structure of the collision term

[TABLE]

where ${\bm{K}}={\bm{K}}[{\bm{B}},\widehat{{\bm{C}}},\widehat{{\bm{D}}}]$ is an operator depending on other species distribution functions integrated over momenta, and $\widehat{\bm{K}}={\bm{K}}[\widehat{{\bm{B}}},{{\bm{C}}},{{\bm{D}}}]$ is its hatted version. Its expression is

[TABLE]

where the momentum dependence ${\bm{A}}(p)$ and ${\bm{B}}(p_{\beta})$, ${\bm{C}}(p_{\gamma})$, ${\bm{D}}(p_{\delta})$ (and similarly for Pauli blocking operators) are omitted for a more compact notation. Since${\bm{K}}$, ${\bm{A}}$ and $\mathds{1}$ are all Hermitian, it is obvious from Eq. (III.2.2) that so is the collision term. Furthermore, its structure is again manifest. The first line corresponds to scattering out processes. As for the second line, it corresponds to the scattering in processes, and differs only by an overall sign and the exchange of the distribution and Pauli blocking functions.

This collision term $C[{\bm{A}}]$, being itself an operator in spinor space, can be decomposed into its covariant parts $I_{C[A]}$ and ${\cal Q}_{C[A]}^{\mu}$ as in the decomposition (45). These components can be found by multiplying by the appropriate $X\in{\cal O}$ and taking the trace, that is using the extraction (308). Since all operators involved in the collision term are made of $\gamma^{\mu}$ or $\gamma^{5}$ matrices, the problem is reduced to taking traces of products of these operators (Fidler and Pitrou, 2017, App. C). This systematic computation can be handled by a computer algebra package such as xAct Martín-García (2004) and this is particularly powerful since it also takes care of all simplifications involving space-time indices.

In particular, when using Eqs. (308) to extract the intensity part of the collision term (III.2.2), we find

[TABLE]

which is compactly written as

[TABLE]

Reactions related to the reaction (116) by crossing symmetry are deduced by replacing the operators describing the distributions by those of the antiparticle, and changing distribution operators for Pauli-blocking operator. For instance the collision term for $a+\bar{c}\leftrightarrow\bar{b}+d$ is deduced by ${{\bm{B}}}\to\widehat{\overline{{\bm{B}}}}$ and $\widehat{{\bm{C}}}\to{\overline{{\bm{C}}}}$, where the bar indicates that we consider the operator associated to the antiparticles [see Eq. (45)]. Similarly the reaction $\bar{a}+\bar{b}\leftrightarrow\bar{c}+\bar{d}$ is obtained by a global charge conjugation, where all operators are replaced by the one associated to the antiparticle. From the decomposition (45) it is obviously equivalent to $m\to-m$ for all masses. Finally the intensity part of the collision term for the species $a$ in the reaction (116) is the same as the intensity part fo the species $b$ 555When focusing on the polarization part of the collision term, this is no longer the case (Fidler and Pitrou, 2017)., and if we are to compute the collision term for $c$ or $d$ we need only to change the global sign.

III.3 General collision term

Let us now restrict to the case where all particles are unpolarized, the general case being detailed in Fidler and Pitrou (2017). More specifically, we assume that massive particles (such as electrons, positrons neutrons or protons) are unpolarized, that is for these species ${\cal Q}^{\mu}=0$. For these particles we define666The notation $f_{\mu}^{\pm}$ is obviously useless but we keep it as it is a particular case of the general case when species are polarized, which is considered in detail in Fidler and Pitrou (2017). Furthermore it allows to write the general collision term (156).

[TABLE]

However, for neutrinos, when considered as strictly massless, circular and linear polarization are separate concepts. We still assume that they do not have linear polarization. However, neutrinos have circular polarization since there are only left-helical neutrino and right-helical antineutrino states. We define for neutrinos

[TABLE]

where $\lambda=1$ for particles (neutrinos) and $\lambda=-1$ for antiparticles (antineutrinos). In fact given the left-chirality of weak interactions for neutrinos, we have $V=-\lambda I$ such that the previous definition reduces to

[TABLE]

Hence for neutrinos we also define

[TABLE]

$f$ is the distribution function per helicity state, which has a clear meaning if the distribution is unpolarized, and in thermal equilibrium it reduces to a Fermi-Dirac distribution. Since massless neutrinos exist only with left chiralities, that is left helicities $f=I$, whereas for other fermionic massive species, $f=I/2$ since they exist in two different helicities.

Under all these restrictions and with these definitions, the intensity part of the collision term is reduced to

[TABLE]

where the Kernel takes the general form (using the generic notation (46) for masses)

[TABLE]

The Kernel can be separated into a squared amplitude and a phase space in the form

[TABLE]

such that the collision term (155) is reduced to

[TABLE]

III.4 Standard reactions with neutrinos

Let us review the standard two-body reactions for neutrinos. These are required to describe the decoupling of neutrinos in the early universe (Dolgov et al., 1997; Mangano et al., 2005; Grohs et al., 2016; Froustey and Pitrou, 2020). We consider the various type of reactions one by one, and we summarize the results in table 1.

In the particular case that the species $a$ and $c$ are neutrinos or antineutrinos, that is can be considered as massless, and their coupling is only left-chiral ($\epsilon_{+}^{ac}=0$), we find

[TABLE]

III.4.1 Muon decay

The muon decay is due to the interaction between the muon ($\mu^{-}$)/muon neutrino ($\nu_{\mu}$) charged current and the electron ($e^{-}$)/neutrino ($\nu$) charged current. Furthermore it involves only left-chiral couplings. It thus corresponds to the case

[TABLE]

We remind that for the decay reaction $a\leftrightarrow\bar{b}+c+d$, the collision term is deduced from the reaction $a+b\leftrightarrow c+d$ by crossing symmetry. The collision term deduced from the general form (158) is therefore

[TABLE]

and where it is stressed by a barred notation that the covariant quantities related to the species $a$, $c$ and $d$ are those of particles, and those for the species $\bar{b}$ are those of antiparticles. The muon lifetime is recovered from this collision term evaluated at null spatial momentum of (${\bm{p}}=0$), and ignoring Pauli blocking effects, thanks to the definition $I_{C[A]}({\bm{p}}=0)={\rm d}I_{a}/{\rm d}t({\bm{p}}=0)\equiv-\Gamma_{a}I_{a}({\bm{p}}=0)$. We get

[TABLE]

and this is exactly the expression that would be obtained from the Fermi golden rule.

III.4.2 neutrino/muon neutrino scattering

The interactions between neutrinos of different types (e.g. electronic neutrinos and muonic neutrinos) are only due to neutral currents with a pure left chiral coupling. The effect of the reaction $\nu+\nu_{\mu}\leftrightarrow\nu+\nu_{\mu}$ thus corresponds to the case

[TABLE]

Using Eq. (III.4), the covariant parts of the collision term take the form

[TABLE]

The effect of the reaction $\nu+\bar{\nu}_{\mu}\leftrightarrow\nu+\bar{\nu}_{\mu}$, which in our general notation is $a+\bar{d}\leftrightarrow c+\bar{b}$, is obtained by a simple crossing symmetry. For instance the intensity part of the Kernel would be for that process

[TABLE]

and $|M|^{2}$ given by (III.4.2).

For completeness, we must stress again that the effect of antineutrino-muonic antineutrino reactions ($\bar{\nu}+\bar{\nu}_{\mu}\leftrightarrow\bar{\nu}+\bar{\nu}_{\mu}$) on antineutrinos is obtained by charge conjugation, that is by considering the case

[TABLE]

This means that the collision term takes the same form as Eqs. (III.4.2) but where all covariant components should now refer to antiparticle species. For instance the intensity part takes the form

[TABLE]

III.4.3 neutrino/neutrino scattering

Neutrino-neutrino scattering ($\nu+\nu\leftrightarrow\nu+\nu$) and neutrino-antineutrino scattering ($\nu+\bar{\nu}\leftrightarrow\nu+\bar{\nu}$) are special cases of the previous electronic neutrino-muonic neutrino scattering but there are a few crucial differences in the derivation of the collision term which are detailed in Fidler and Pitrou (2017).

To summarize, when considering interactions between neutrinos ($\nu+\nu\leftrightarrow\nu+\nu$) one must consider the two-body case (III.4.2) in the particular case $a=b=c=d=\nu$ and multiply the result by a factor $2$ [this point was omitted in Hannestad and Madsen (1995)]. And when considering interactions between neutrinos and antineutrinos of the same flavor ($\nu+\bar{\nu}\leftrightarrow\nu+\bar{\nu}$) one must consider the two-body interaction in the particular case $a=c=\nu$, $b=d=\bar{\nu}$ and multiply the result by a factor $4$ in agreement with Dolgov et al. (1997). In particular, a simple crossing symmetry allows to get the former reactions $\nu+\nu\leftrightarrow\nu+\nu$ from the $\nu+\bar{\nu}\leftrightarrow\nu+\bar{\nu}$ only up to a factor $1/2$. We can interpret this reduction by a factor two using the fact that outgoing particles are identical and one must not double count the outgoing states.

III.4.4 neutrino/electron scattering

Contrary to neutrino-neutrino scattering, electron-neutrino scattering is due to both charged and neutral currents. However the Fierz reordering reduces the problem to an interaction of neutral currents with modified chiral couplings. Using Eqs. (121) and (126), the effect of $\nu+e^{-}\leftrightarrow\nu+e^{-}$ on neutrinos corresponds to the case

[TABLE]

which must be used in Eq. (III.4).

The effect of $\nu+e^{+}\leftrightarrow\nu+e^{+}$ is obtained by a crossing symmetry. The effect of $\bar{\nu}+e^{+}\leftrightarrow\bar{\nu}+e^{+}$ on antineutrinos is obtained from charge conjugation of (166), that is it corresponds to the case

[TABLE]

and the effect of $\bar{\nu}+e^{-}\leftrightarrow\bar{\nu}+e^{-}$ is obtained from crossing symmetry.

Finally, we can check that in the unpolarized case, these results for neutrino/electrons interaction and those for neutrino/neutrinos interactions obtained in § III.4.2 and III.4.3 are exactly the results of Grohs et al. (2016). However note that as mentionned in this reference, there is a typo in the annihilation of neutrino and antineutrinos into electrons and positron in tables $1$ and $2$ of Dolgov et al. (1997), and thus Tables 1.5 and 1.6 of Lesgourgues et al. (2013). The process described in these tables should be of the form $\nu+\bar{\nu}\leftrightarrow e^{-}+e^{+}$ and not $\nu+\bar{\nu}\leftrightarrow e^{+}+e^{-}$. Up to this typographical correction our results agree also with Dolgov et al. (1997); Lesgourgues et al. (2013) and we gather all reactions in Table 1.

IV Neutrons-protons conversions

Neutron-proton conversions are controlled by weak interactions in the early universe. As they enforce statistical equilibrium, and since the neutron is more massive and thus less likely statistically, the frozen neutron abundance depends directly on the reaction rates. For larger reaction rates, the frozen abundance is smaller and thus it leads to less primordial Helium production (Pitrou et al., 2018). We now review the general form of these rates and we detail how a Fokker-Planck expansion can be used to compute them in practice.

IV.1 General expression of the rates

Let us first consider the reactions

[TABLE]

They are mediated by the coupling of the neutron ($n$)/proton ($p$) charged current and the neutrino/electron charged current. While the latter is purely left chiral, the former has both chiral couplings due to the effective constant $g_{A}$ defined in Eq. (124). Hence we must consider the case

[TABLE]

With unpolarized species, the collision term for the forward reaction (168) takes the simpler form

[TABLE]

where the coupling constants are

[TABLE]

and the left-left right-right and left-right chiral couplings are

[TABLE]

All other reactions are related by crossing symmetry or time reversal, which affect only the phase space, but not $|M|^{2}$, that is we only need to make sure to put the distribution function $f$ for initial particles and the Pauli-blocking factor $1-f$ for final particles.

The number density of nucleons $N=n,p$ is related to the distribution function $f_{N}$ by

[TABLE]

Hence from Eq. (169) we can define reaction rates for the densities of neutrons and protons. The forward rates $\Gamma_{n\to p}$ are of the form

[TABLE]

where $\alpha_{e}=1$ (resp. $\alpha_{e}=-1$) if the electron/positron is in the initial (resp. final) state, and with a similar definition for the neutrino/antineutrino coefficient $\alpha_{\nu}$. Hence, Eq. (IV.1) describes all reactions (168). Note, that we have neglected Pauli-blocking effects of the final proton, since the baryon-to-photon ratio is very low. However we have correctly included Pauli-blocking effects of electrons/positrons and neutrinos/antineutrinos since for a Fermi-Dirac (FD) distribution without chemical potential

[TABLE]

The vanishing of the electron/positron chemical potential is enforced by the very low baryon-to-photon number ratio (Pitrou et al., 2018, App. A.2). However, if we want to investigate the possibility of non-vanishing neutrino chemical potentials $\mu$ (Pitrou et al., 2018; Serpico and Raffelt, 2005; Iocco et al., 2009; Simha and Steigman, 2008), once must use instead

[TABLE]

IV.2 Isotropy of distributions

At low temperature, it is enough to assume that nucleons follow an isotropic Maxwellian distribution of velocities at the plasma temperature $T$. Hence the following integrals are obtained

[TABLE]

In particular contracting with $\delta_{ij}$ we recover the expression for the pressure of nucleons in the low temperature limit

[TABLE]

For electron or neutrino distributions, since we have assumed isotropy, we deduce the property

[TABLE]

where $\alpha$ and $\beta$ are some numbers. From isotropy we also find that

[TABLE]

Hence for all practical purposes, we can perform the replacements

[TABLE]

on all species, resulting in great simplifications.

IV.3 Expansion in the energy transfer

The integral in (IV.1) is $8$-dimensional when on removes the Dirac function. Due to the isotropy of all distributions, this can be reduced to a $5$-dimensional integral. This is the method followed by Lopez et al. (1997). Here we follow a much simpler route by performing a Fokker-Planck expansion, that is an expansion in the momentum transferred to the nucleons. It consists in expanding the energy difference between the nucleons, $E_{n}-E_{p}$ around the lowest order value

[TABLE]

As we shall see, this results in one-dimensional integrals which are much faster to evaluate.

We evaluate the rates by performing an expansion in powers of $\epsilon\equiv\sqrt{T/m_{N}}$. To evaluate the order of each term, we consider that the momentum or energies of neutrinos are of order $T\sim\Delta$, that is factors of the type $E_{e}/m_{N}$ or $E_{\nu}/m_{N}$ are of order $\epsilon^{2}$. Furthermore, from (177) a factor ${\bm{p}}_{n}/m_{n}$ is of order $\sqrt{T/M}\sim\sqrt{\Delta/M}$ and thus $\epsilon$. However since only even powers of the spatial momentum of nucleons must appear [see Eqs. (177)], we shall encounter terms of the type $|({\bm{p}}_{p}/m_{n})^{2}|$ which are of order $\epsilon^{2}$.

Keeping only the lowest corrections this expansion reads

[TABLE]

[TABLE]

where ${\bm{q}}\equiv{\bm{p}}_{p}-{\bm{p}}_{n}=\alpha_{\nu}{\bm{p}}_{\nu}+\alpha_{e}{\bm{p}}_{e}$ is the spatial momentum transfered. The first term in (183) is the lowest order, or Born approximation, that is the only appearing when considering the infinite nucleon mass approximation. The second term is an order $\epsilon$ correction, and the third term is an order $\epsilon^{2}$ correction. Finally the last term is of order $T\Delta/m_{N}$ so it is an order $\epsilon^{2}$ correction as well. It is the only corrective term for which it is crucial to take into account the difference of mass between neutrons and protons. Using Eq. (183), we expand the Dirac delta function on energies as

[TABLE]

where $\Sigma\equiv\Delta+\alpha_{e}E_{e}+\alpha_{\nu}E_{\nu}$.

We must then expand the matrix element and the energies appearing in Eq. (IV.1). It proves much easier to expand all these contributions together. Furthermore, whenever a term is already of order $\epsilon^{2}$, we know that it should multiply only the Born term of the expansion (185), so we can apply the simplification rule (181). With this method we find

[TABLE]

The second term in Eqs. (186) and (186) is of order $\epsilon$ and the last term in these equations is of order $\epsilon^{2}$. Hence the second term needs to be coupled with the order $\epsilon$ term in the Dirac delta expansion (185) which is $\delta^{\prime}(\Sigma)\delta Q_{1}$, and simplified with the rules (181).

There are four steps to complete this Fokker-Planck expansion.

First, using Eqs. (186) and (185) in the reaction rates (IV.1) we perform the integral on the initial neutron momentum with the rules (177). 2. 2.

Second, we can replace the differential elements for the integral on electron and neutrino momenta with ${\rm d}^{3}p\to 4\pi p^{2}{\rm d}p$ because we have already performed all angular averages. 3. 3.

We are left with a two dimensional integral on the electron and neutrino momentum magnitudes $p_{e}=|{\bm{p}}_{e}|$ and $p_{\nu}=|{\bm{p}}_{\nu}|$. Let us note $E_{\nu}=p_{\nu}$ in order to write the result in a easily readable form. Third, we perform the integral on $E_{\nu}$ using the Dirac delta and its derivatives. Whenever a Dirac delta derivative appears, it means that we have to perform integration by parts to convert it into a normal Dirac delta. This will introduce derivatives with respect to the $E_{\nu}$ applied on the neutrino distribution function or Pauli-blocking factor. Also for a given reaction it might appear that the value of $E_{\nu}$ constrained by the Dirac delta is not physical for that reaction if $\alpha_{\nu}=1$ and physical if $\alpha_{\nu}=-1$, or vice-versa. This is the reason why we consider the total reaction rate of the reactions (168) and (168). Once their rates are added, the Dirac delta automatically selects either the neutrino in the initial state, with the corresponding distribution function, or the antineutrino in the final state, with the associated Pauli-blocking factor. Eventually once the rates (168) and (168) are added, we might forget about $\alpha_{\nu}$, that is about the position of the neutrino. We need only to compute two rates, one where the electrons is in the initial state [reaction (168)], and one where it is a positron which is in the final state [the sum of reactions (168) and (168)]. 4. 4.

Finally, we need to determine the procedure to convert the rate with a neutron in the initial state into the reverse rate with a proton in the initial state. Even if the matrix element is the same for all reactions, the method to perform a finite mass expansion is not symmetric under the interchange $p\leftrightarrow n$. Indeed we chose to expand the momentum of the final nucleon around the initial one, and we remove the integral on the final nucleon momenta. It is apparent on Eqs. (172) that the electron (resp. neutrino) momentum is contracted with the neutron (resp. proton) in the $LL$ term but this is the opposite in the $RR$ term. Since the coupling factors of these terms are interchanged by the replacement $g_{A}\to-g_{A}$, we can deduce the rates with an initial proton from those with an initial neutron using the rule $g_{A}\to-g_{A}$. Obviously the argument of the Dirac delta contains now $E_{p}-E_{n}=-\Delta+\dots$ instead of $E_{n}-E_{p}=\Delta+\dots$ so we must also apply the rule $\Delta\to-\Delta$. Finally when considering a reverse reaction, the electron in the initial state turns into a positron in the final state so we must also apply the rule $E_{e}\to-E_{e}$, that is change the electron distribution function to a Pauli-blocking factor or vice-versa.

Having sketched the details of the procedure, we are in position to give the results. In the next section, we report the lowest order reaction rates in § IV.4, also called Born approximation rates. The first corrections, that we call finite nucleon mass corrections, are reported in appendix C.

IV.4 Lowest order $n\leftrightarrow p$ reaction rate

Let us note $g(E)$ the Fermi-Dirac distribution at temperature of electrons $T$ and $g_{\nu}(E)$ the Fermi-Dirac distribution at the neutrino temperature $T_{\nu}$, that is

[TABLE]

At lowest order in the Fokker-Planck expansion, the reaction rates take simple forms. First, the factors entering the matrix element reduce to

[TABLE]

as seen from Eqs. (186). The last equality is correct only if it is understood that an angular average either on electrons momentum or neutrino momentum is performed, that is using the rule (181). Hence from Eq. (IV.1), we find the Born rates Brown and Sawyer (2001); Lopez and Turner (1999); Weinberg (1972); Bernstein et al. (1989); Pitrou et al. (2018)

[TABLE]

with $E=\sqrt{p^{2}+m_{e}^{2}}$ and

[TABLE]

The first contribution in Eq. (189) corresponds to the $n\to p$ processes (168) and (168) added, that is for all processes where the electron is in the final state. It can be checked indeed that the electron distribution is evaluated as $g(-E)=1-g(E)$. Furthermore, if the neutrino is in the initial state (when $E>\Delta$) its energy is $E_{\nu}=E-\Delta$ and its distribution function appears as $g_{\nu}(E_{\nu})$, but if it is in the final state (when $E<\Delta$) its energy is $E_{\nu}=\Delta-E$ and the neutrino distribution function is evaluated as $g_{\nu}(E-\Delta)=1-g_{\nu}(\Delta-E)$.

The second term of Eq. (189) corresponds to the reaction (168), that is to the process where the positron is in the initial state. The energy of the positron is $E$ and its distribution function appears as an initial state [$g(E)$], whereas the neutrinos in the final state have energy $E_{\nu}=\Delta+E$ and their distribution function appear thus as Pauli-blocking factor $g_{\nu}(-E-\Delta)=1-g_{\nu}(E+\Delta)$.

The reaction rate for protons, that is $\Gamma_{p\to n}$, is obtained by the simple replacement $\Delta\to-\Delta$, which amounts to $\chi_{+}\to\chi_{-}$. We give it for completeness

[TABLE]

Similarly the second term corresponds to the reverse processes (168) and (168) added since the electron distribution function is always in an initial state [$g(E)$], and the neutrino is in the initial or final state depending on the sign of $E_{\nu}=-E+\Delta$. The first term corresponds to the reverse process (168) with the positron always in the final state [$g(-E)=1-g(E)$] and the neutrino always in the initial state [$g_{\nu}(E+\Delta)$].

Finally, note that using

[TABLE]

we get in the case of thermal equilibrium between neutrinos and the plasma (that is when $T_{\nu}=T$)

[TABLE]

This implies that if neutrinos have the same temperature as the plasma, the reaction rates satisfy the Born approximation detailed balance relation (Brown and Sawyer, 2001; Pitrou et al., 2018)

[TABLE]

Compton scattering

The Fokker-Planck expansion exposed in § IV.3 was inspired from a similar expansion often used for Compton scattering. We have already stressed the numerous similarities between the construction of distribution functions for fermions and bosons. We now turn to the computation of the collision term for photons associated with Compton scattering onto electrons. As we detail in the next section, the structure of the Compton collision term is nearly identical to the weak interaction collision term except that stimulated emission factors replace Pauli-blocking ones. The collision term obtained is exposed in § VI for isotropic distributions (along with a discussion on its cosmological implications) and in § VII for the general case of anisotropic distributions.

V Compton collision term

V.1 Extended Klein-Nishina formula

We consider the Compton reaction

[TABLE]

The initial photon and electron momenta are decomposed as

[TABLE]

with similar decompositions for the final particles. Throughout this part, the electron mass is noted $m$.

Even though the Hamiltonian of QED accounts for a vertex between the electronic current and a single photon, that is it is a three-leg vertex, it is more adapted to consider an effective QED Hamiltonian assuming that the electron propagates freely between two interactions with photons (Beneke and Fidler, 2010, §. II.D). This is essentially similar to our treatment of weak interactions in the Fermi theory of § III.1, except that we do not use that the propagator of the internal electron line is dominated by the electron mass. The effective interaction Hamiltonian takes the form (III.2.2) with $a,c=\gamma$ and $b,d=e^{-}$ and a matrix element

[TABLE]

[TABLE]

The two terms of $\cal O$ correspond to the two possible Feynman diagrams associated with the reaction (197). Then, the procedure to obtain a collision operator exactly follows our derivation in § III.2.2, that is we also obtain Eq. (143), with the only difference that the hatted notation on photons now refers to stimulated emission factors instead of Pauli-blocking factors, because we use Eq. (113) instead of Eq. (112) when ordering operators. Once expressed as a collision term for an operator by contraction with ${\bm{\epsilon}}^{\star}_{s}{\bm{\epsilon}_{s^{\prime}}}$ [as in Eq. (58)], it takes a form fully similar to Eq. (III.2.2). We prefer to report it with explicit indices for all operators. Furthermore, we assume that electrons are unpolarized and they are thus described by their distribution function (per helicity) $g(q)=I_{e^{-}}(q)/2$. We also assume that we can neglect the associated Pauli-blocking factors as the baryon-to-photon ration is very low. The equivalent of Eq. (III.2.2) for photons under Compton scattering finally reads

[TABLE]

where the dependence of ${\cal H}^{\mu}_{\,\,\nu}$ and $K^{\mu}_{\,\,\nu}$ on $p$ has been omitted. The operators involved are

[TABLE]

where as in Eq. (144) we defined

[TABLE]

and where the detailed form of ${\cal M}$ is reported below in Eq. (V.1). Note that for photons, the equivalent of the identity operator (145) is ${\cal H}^{\mu}_{\,\nu}=\sum_{ss^{\prime}}\delta^{K}_{ss^{\prime}}\epsilon^{\star\mu}_{s}(\epsilon_{s^{\prime}})_{\nu}$ as seen on Eq. (51), and the structure is totally similar to Eq. (III.2.2), but with the stimulated emission operators

[TABLE]

If $g(q^{\prime})=g(q)$ was satisfied, the general collision term (V.1) would take an even simpler form. It is thus convenient to define

[TABLE]

so as to recast it under the sum of two contributions as

[TABLE]

The first contribution is the gain minus loss contribution (Beneke and Fidler, 2010, § IV.B)

[TABLE]

where we omitted to write the $(p,q,p^{\prime},q^{\prime})$ dependence on $\cal M$ and with the dependence of ${\cal H}^{\mu}_{\,\,\nu}$ omitted when it is on $p$. This part of the collision term is linear in the distribution function and is not impacted by stimulated effects.

The second contribution is the pure gain term

[TABLE]

This gain term is however affected by stimulated effects, and it requires $\delta_{g}\neq 0$. Hence we already see that for very heavy electrons (in the sense that their rest mass energy is much larger than the typical energy of photons) corresponding to the Thomson limit of Compton scattering, we would have $\delta_{g}(q,q^{\prime})\simeq 0$ and no stimulated emission effects.

The squared amplitude of the QED process associated to Compton collisions is

[TABLE]

Using (199), and using the initial electron frame to define polarization vectors, we checked that [using xAct Martín-García (2004) to handle products of Dirac matrices and contractions of vectors]

[TABLE]

in agreement with Portsmouth and Bertschinger (2004, Eq. 187) or Stedman and Pooke (1982). Eq. (V.1) is the extended Klein-Nishina formula for Compton scattering. The Thomson cross-section ${}\sigma_{\rm T}=8\pi\alpha_{\rm FS}^{2}/(3m^{2})$, where $\alpha_{\rm FS}\equiv e^{2}/(4\pi)$, has been introduced and it is related to the electron charge by $6\pi{}\sigma_{\rm T}m^{2}=e^{4}$. This is covariantized by contraction with photon polarization vectors as

[TABLE]

and it then takes the form (Portsmouth and Bertschinger, 2004, Eq. 189)

[TABLE]

where we use the combinations of the screen projectors (note that these are screen projectors in the frame of the initial electron)

[TABLE]

V.2 Kinematics

Let us detail the kinematics of the reaction (197) enforced by the Dirac function on momenta. From $(p+q)^{2}=(p^{\prime}+q^{\prime})^{2}$, which leads to $p\cdot q=p^{\prime}\cdot q^{\prime}$, and expressing the energy and spatial momentum of the final electron using energy and momentum conservation, we find (Chluba et al., 2012, Eq. C7c)

[TABLE]

We define a fractional photon energy shift

[TABLE]

which is deduced from Eq. (211). Note that in the initial electron frame $\Gamma=1$, $\beta^{i}=0$ and

[TABLE]

In particular, when considering the initial electron frame, we find for the prefactor of the last term of Eq. (V.1)

[TABLE]

V.3 Electron velocity distribution

Given that we neglect the degeneracy of electrons, the electron distribution function is well described by a Boltzmann distribution, and Pauli-blocking effects can be ignored. In the frame comoving with the bulk velocity we have (Challinor et al., 2000)

[TABLE]

where $K_{2}$ is a modified Bessel function. The energy of a given electron with momentum $q^{\mu}$ in the bulk electron frame is given by

[TABLE]

where $u^{\mu}_{\rm b}=\gamma(1,{\bm{v}})$ (with $\gamma\equiv 1/\sqrt{1-v^{2}}$) is the electron bulk velocity. The distribution of the final electron $g(q^{\prime})$ is related to $g(q)$ using (Challinor and Lasenby, 1999; Chluba et al., 2012)

[TABLE]

However for our purpose it is enough to consider that electrons follow a Maxwell-Boltzmann distribution

[TABLE]

The first few moments of this distribution are simply

[TABLE]

As long as we restrict our expansion to second order contributions in the electron velocities (see next section), using the distribution (218) and the integrations (219) is valid. In the more general case, considered in e.g. Challinor (1998); Challinor and Lasenby (1999); Challinor et al. (2000), one must use the distribution (215) and Eq. (217) to expand $g(q^{\prime})$ around $g(q)$.

V.4 Fokker-Planck expansion

The Fokker-Planck expansion is an expansion in the momentum transferred to the heavy species. As in the case of weak interactions detailed in § IV.3, we consider that terms of order ${\bm{v}}$ or ${\bm{\beta}}$ are of order $\epsilon$, whereas terms of order $T_{e}/m$ or $E/m$ are of order $\epsilon^{2}$. We list all steps required to perform such expansion.

V.4.1 Energy shift

There are two equivalent methods to perform this expansion. The first method, as used by Chluba et al. (2012) [see also Peskin and Schroeder (1995, § 5.5)] consists in simplifying the integrations (202) using

[TABLE]

where Eq. (211) was used in the last step. The integration on $E^{\prime}$ is then trivially performed and it amounts to the replacement $E^{\prime}\to\phi(E,{\bm{n}},{\bm{n}}^{\prime},{\bm{\beta}})$. In practice all quantities which depend on $E^{\prime}$ (among which the distribution functions) are expanded in the small parameter $\phi-E$.

The second method is the long-standing method for the CMB computations (Hu et al., 1994), among which Dodelson and Jubas (1995); Bartolo et al. (2006); Pitrou (2009a); Beneke and Fidler (2010). After integration on the spatial momentum of the final electron ${\bm{q}}^{\prime}$, the remaining Dirac function ensuring conservation of energy is handled via a Taylor expansion. Eq. (211) cannot be used as it results precisely from the Dirac function, and the final photon energy $E^{\prime}$ must remain unknown at this stage. The Taylor expansion reads

[TABLE]

where the second line is of order $\epsilon$ and the last line of order $\epsilon^{2}$. The expression of the energy transfer that we have used to write this Taylor expansion is essentially similar to Eqs. (184), except that the equivalent of the term (184) does not exist since the final massive particle is the same as the initial massive particle (both are electrons).

It must be understood that integration by parts must be performed to remove derivatives on Dirac functions. Since ${\bm{p}}^{\prime}=E^{\prime}{\bm{n}}^{\prime}$, these integration by parts also act on the expansion (221) itself, in addition to acting on the distribution function. In both methods, it is assumed that the spectrum is smooth enough such that the Taylor expansion of the spectrum is meaningful. In the case of a spectrum containing narrow lines, the Fokker-Planck expansion might not be applicable (Sazonov and Sunyaev, 2000). Furthermore, for the Taylor expansion to remain meaningful, one must ensure that the various moments of the scattering Kernel remain small (Sarkar et al., 2019). This is equivalent to asking the conditions $E/m\ll 1$ and $p\ll 1$ or equivalently $T_{e}/m\ll 1$.

V.4.2 Electron distribution

Furthermore, in both methods, one also expands $\delta_{g}(q,q^{\prime})$ and ${\cal M}$ [Eq. (V.1)]. The first expansion is given by

[TABLE]

In the first method, it is understood that ${\bm{p}}^{\prime}=\phi{\bm{n}}^{\prime}$ and one must further expand in $\phi-E$. In the second method we use ${\bm{p}}^{\prime}=E^{\prime}{\bm{n}}^{\prime}$ so that ${\bm{p}}^{\prime}$ is acted upon by the integrations by parts as $\partial_{E^{\prime}}{\bm{p}}^{\prime}={\bm{n}}^{\prime}$.

V.4.3 Squared amplitude

The expansion of Eq. (V.1) is simplified when restricting to orders smaller or equal to $\epsilon^{3}$, as we need only consider the first term in its r.h.s given the property (214). Its expansion is deduced from its relation to the screen projector in the definitions (210). The expansion of the screen projector in the initial electron frame ${\cal H}_{\mu\nu}(q/m,p)$ is obtained from the transformation rule (53) considered for $\tilde{u}^{\mu}=q^{\mu}/m$, and using $\tilde{E}=\Gamma E(1-{\bm{n}}\cdot{\bm{\beta}})$. For instance restricting to order $\epsilon$, we find

[TABLE]

While the expansion of Eq. (V.1) up to order $\epsilon^{2}$ is rather sizable, it is instructive to consider the case in which the distribution function is unpolarized, restricting also to the intensity part of the collision term. One must thus consider

[TABLE]

Considering only the contribution of the first term in the r.h.s of Eq. (V.1), that is neglecting terms of order $\epsilon^{4}$, we find using property (54) that

[TABLE]

where the direction in the initial electron frame of a momentum $p$ is defined as in Eq. (50), that is $p^{\mu}=\tilde{E}(\tilde{u}^{\mu}+\tilde{n}^{\mu})$. This is the usual Thomson squared amplitude in the electron frame $\propto 1+\cos^{2}\theta$ where $\theta$ is the deflection angle. Eventually, we must expand $\tilde{{\bm{n}}}\cdot\tilde{{\bm{n}}}^{\prime}$ around of ${\bm{n}}\cdot{\bm{n}}^{\prime}$ in order to obtain its expansion. From the invariance of $p^{\prime}_{\mu}p^{\mu}$ we get

[TABLE]

This is the procedure followed by Dodelson and Jubas (1995) or Hu (1995, §2.2).

V.4.4 Integration measure

Finally, the expansion is completed by expanding the energies of electrons which appear in the relativistic integration elements contained in (202). Up to order $\epsilon^{2}$ we find

[TABLE]

Furthermore, when using the first method for handling the energy shift, the factor $1/\Gamma/(1-{\bm{\beta}}\cdot{\bm{n}})$ appearing in (220) must be expanded in powers of ${\bm{\beta}}$.

V.5 Structure of the expansion

We checked that both methods for the energy shifts lead to the same result777The computations were performed with xAct Martín-García (2004). and we present and analyze them in the subsequent sections. Let us briefly comment on the prominent features of the method.

In the integration (202), we had nine integrals (on ${\bm{q}}$, ${\bm{p}}^{\prime}$ and ${\bm{q}}^{\prime}$). Those on ${\bm{q}}^{\prime}$ were removed using the spatial part of the Dirac function. The integration on ${\bm{q}}$ is subsequently performed using the moments (219), once the expansion in the momentum transferred is written. Finally, the integration on $E^{\prime}$ is performed either as in (220) with the first method, or from the Taylor expansion (221) in the second method. Eventually we are left with an integration on ${\bm{n}}^{\prime}$, that is the direction of the final photon. As we shall detail, this is handled from the multipolar decomposition of the distribution function.

The gain minus loss terms are always linear in the distribution function. The loss terms are particularly simple to compute because they do not depend on the distribution of the final photon, hence the residual integration on ${\bm{n}}^{\prime}$ is simple.

The stimulated emission factors appear only as a result of $\delta_{g}\neq 0$. As a consequence, they do not appear in the Thomson limit of the collision term, even for a general non-isotropic distribution function.

The various contributions of the collision term can be classified as follows.

• •

Thomson terms are the lowest order ones (order $\epsilon^{0}$). They are entirely due to the structure of ${\cal P}^{\mu\nu\alpha\beta}(q,p,p^{\prime})$ (when the velocity of the initial electron is the same as the one of the observer) coupled with the structure of the gain minus loss term (V.1).

• •

Thermal terms are of order $\epsilon^{2}$, and they appear whenever we average over the electron distribution products of the type $q^{i}q^{j}$, by using Eq. (219).

• •

Kinetic terms are proportional to $v^{i}$ and arise from the order $\epsilon$ terms using Eq. (219). Non-linear kinetic terms arise in the same condition as the thermal terms, that is they are of order $\epsilon^{2}$, and are proportional to $v^{i}v^{j}$ or $v^{2}$. One can conveniently check thermal terms by replacing $v^{i}v^{j}\to\delta^{ij}T_{e}/m$.

• •

Recoil terms arise from the fact that the electron mass is not infinite. They are of order $\epsilon^{2}$ and are proportional to $E/m$. They originate from the energy shift ${\delta_{E}}$, and also from $\delta_{g}$ since it is also affected by ${\delta_{E}}$.

Finally, the factor $n_{e}{}\sigma_{\rm T}$ appears as a prefactor to all terms, since collisions are proportional to the number densities of electrons and to the Thomson cross section. It is customary to define the optical depth $\tau$ by

[TABLE]

Once divided by ${\rm d}\tau/{\rm d}t$, the collision term is reinterpreted as a rate of variation per unit of $\tau$ instead of per unit of time.

VI Evolution of isotropic distributions

In this section we restrict to the case of an isotropic distribution function. Given the absence of preferred directions, there is no linear polarization (${\cal P}_{ij}=0$), and we further assume that there is no circular polarization ($V=0$). Isotropy also implies that intensity is equal to its monopole $I=I_{\emptyset}$.

VI.1 Kompaneets equations

Under these symmetries, it is shown that there is no contribution from the Thomson terms, since scattering out and scattering in term exactly cancel. Furthermore let us consider the situation in the bulk frame of baryons (${\bm{v}}=0$). The only possible terms, when considering the expansion up to order $\epsilon^{2}$, are the thermal and recoil terms.

The collision term reduces then to the celebrated Kompaneets equation  Kompaneets (1957)

[TABLE]

We notice that if $f_{\gamma}\equiv I_{\emptyset}/2$ (as in § III.3, this is the distribution per helicity state) is a Planck distribution at temperature $T_{\gamma}$

[TABLE]

then it satisfies

[TABLE]

and the Kompaneets collision term is recast as

[TABLE]

The factor multiplying $T_{e}-T_{\gamma}$ (where $I_{\emptyset}/2$ is a Planck spectrum) is exactly the spectral shape of a $y$-type spectral distortion. If the photon temperature is equal to the electron temperature, it vanishes. Or said differently, the Planck spectrum at the electron temperature is a fixed point of the associated Boltzmann equation. It is customary to define a Compton optical depth by

[TABLE]

The Boltzmann equation associated with the Kompaneets collision term then takes the very simple form

[TABLE]

In fact, property (231) holds even if the distribution is a Bose-Einstein distribution, that is with a constant chemical potential $f_{\gamma}(E)={\cal B}[(E+\mu)/T_{\gamma}]$. Indeed, it is obvious under the form (229) that the Kompaneets term conserves the number of photons [as it should since Compton collisions of the type (197) do conserve in general photons]. Since a Planck spectrum at a given temperature $T_{\gamma}$ has a number density of photons uniquely determined by $T_{\gamma}$, a chemical potential must develop when photons thermalize with electrons, so as to ensure the conservation of photons, and the final spectrum is a Bose-Einstein spectrum with non-vanishing chemical potential, see Sunyaev and Zeldovich (1970), Burigana et al. (1991) and Hu (1995, §3.2.2).

Higher order effects in the Fokker-Planck expansion have been considered in the isotropic case, either without electron bulk velocity (Challinor, 1998; Itoh et al., 1998) or with a bulk velocity on top of thermal effects (Challinor and Lasenby, 1999; Nozawa et al., 1998; Itoh et al., 2000; Sazonov and Sunyaev, 1998). These corrections are particularly relevant for the Sunyaev-Zel’dovich effect (Sunyaev and Zeldovich, 1972) of hot galaxy clusters on the CMB, whose lowest order description is Eq. (232) restricted to the limit $T_{\gamma}\ll T_{e}$.

VI.2 Thermalization on electrons

Let us discuss briefly the thermalization process from the Kompaneets equations in two limiting cases. First we consider that electrons are completely dominating the energy content, and then the case where the photon energy density dominates over baryons. The first case applies in the late universe, whereas the second case applies in the radiation dominated era.

VI.2.1 Test distribution

In the case where the energy density of baryons dominate, their temperature is not affected by the back-reaction of the Kompaneets collision term on them. Electrons set the final value of the photons temperature. In Fig. (3) we plot the response of a Planck spectrum to a rapid increase of $0.01\%$ of the electron temperature.

We find that approximately for $\tau^{\rm C}\gtrsim 4$, the spectrum is reasonably converged to a Bose-Einstein distribution with the appropriate chemical potential. However, the case of dominant baryonic energy density applies to the late universe, for instance in the galaxy inter-cluster hot gas, for which the effect of the Kompaneets equation is known as the Sunyaev Zel’dovich (SZ) effect (Sunyaev and Zeldovich, 1972) and in such cases the Compton optical depth is much lower than unity. Hence the associated spectral distortions are expected to be of the $y$-type, that is with a distortion given by the Kernel (232). Note that distortions from all clusters should in principle contribute collectively to a global monopolar distortion of the SZ type including relativistic temperature corrections (Refregier et al., 2000; Hill et al., 2015).

VI.2.2 Dominant distribution

In the opposite case where the energy density of baryons is subdominant, the photons cannot gain energy from their interactions with electrons. Hence the total energy density transfer rate must vanish, that is the electron temperature must be such that $\int I_{C}^{\rm Kom}E^{3}{\rm d}E=0$. It is then given by

[TABLE]

The effect of the Kompaneets equation is to redistribute photons while conserving the total energy density so as to approach a Bose-Einstein spectrum. This case applies to the radiation dominated era. The photon spectrum is redistributed via Compton interactions with electrons, but there is no net creation of photons, nor energy gained or lost. Note that Comptonization of electrons is much faster (Iwamoto, 1983), precisely because their energy density is subdominant and they always possess a thermal spectrum, that is a well defined temperature.

One can still define the Compton optical depth via Eq. (233). In the early universe, any energy injected at a time corresponding to $\tau^{\rm C}\gtrsim 4$ results in a Bose-Einstein spectrum, that is a spectral distortion of the $\mu$-type, whereas energy injected later, and thus corresponding to a low value of $\tau^{\rm C}$ (that is $\tau^{\rm C}\lesssim 0.1$) are of the $y$-type. For the most recent cosmological parameters Ade et al. (2016), the Compton optical depths scales for redshifts $z$ belonging to the radiation era (that is for $z\gtrsim 3\times 10^{3}$) as888The scaling in $(1+z)^{2}$ is deduced from ${\rm d}\tau^{\rm C}/{\rm d}a=n_{e}T_{e}/(Ha)\propto 1/a^{3}$, where $H$ is the cosmic time Hubble function and $a$ the scale factor.

[TABLE]

as depicted on Fig. (4). Hence it is found that $\tau^{\rm C}\simeq 4$ corresponds to $z\simeq 3\times 10^{5}$, above which distortions are of the $\mu$-type, whereas $\tau^{\rm C}\simeq 0.1$ corresponds to $z\simeq 5\times 10^{4}$ below which they are of the $y$-type.

However this simplistic picture is altered by the processes which do not conserve the number of photons (double Compton scattering and Bremsstrahlung). These create photons of low energy and force the low energy part of the spectrum to stick to the Planck distribution at $T_{e}$. Hence, these photon creating processes also solve the issue of negative chemical potentials. Photons in this lower part of the spectrum are subsequently up-scattered by Compton interactions, that is by the Kompaneets collision term, and this tends to decrease the chemical potential. Eventually, the spectrum relaxes fully to a Planck spectrum (Hu, 1995, §3.4.1). For very large redshifts ($z\gtrsim z_{\rm max}\simeq 2\times 10^{6}$), the distortions of the $\mu$-type are no-longer visible (Chluba and Sunyaev, 2012) as they are erased by these photon creating processes. More details about spectral distortions generated during the radiation dominated era can be found in Hu (1995); Chluba et al. (2012); Chluba (2016), and numerical solutions are presented in Chluba and Sunyaev (2012); Chluba (2015).

VI.2.3 Compton cooling

The temperature of massive particles tends to decay like $1/a^{2}$ where $a$ is the scale factor, whenever their temperature is much lower than their mass, which is the case for $z\ll 10^{10}$. Hence, baryons (electrons and nuclei), tend to cool faster than photons. From the Kompaneets term (232) it is seen that colder baryons will tend to extract energy from photon, via the Compton interactions between photons and electrons. Let us evaluate briefly the energy extracted from the photons by the faster adiabatic cooling of baryons, using basic thermodynamics.

The pressure of baryons is very well approximated by $P_{B}=n_{B}T$. Hence the adiabatic evolution of baryons kinetic energy when the volume $V$ expands is

[TABLE]

where $\delta Q$ is the energy brought to the photons by the electrons, which is expected to be negative. If we use that for baryons (that is massive particles) $U_{B}\simeq 3/2n_{B}VT_{B}$, and if we assume that baryons are forced to follow nearly exactly the photon temperature because of Compton interactions (${\rm d}T_{B}/T_{B}={\rm d}T_{\gamma}/T_{\gamma}=-1/3{\rm d}V/V$), then we get

[TABLE]

This results in a total cooling of the photons

[TABLE]

Estimating that Compton cooling is efficient for $z\gtrsim z_{\rm min}\simeq 200$, we can estimate the total energy extracted in the CMB spectrum by integrating from $z_{\rm min}$ up to $z_{\rm max}$. In a Planck spectrum, the average energy of photons is about $2.70T_{\gamma}$, hence $\rho_{\gamma}\simeq 2.70n_{\gamma}T_{\gamma}$. We then use that the baryon-to-photon number density ratio is of order $\eta\simeq 6.1\times 10^{-10}$ Pitrou et al. (2018). However this is the ratio of the number of nucleons to the number of photons. Taking into account that $24.7\%$ in mass is in the form of Helium nuclei Pitrou et al. (2018), and adding also the contribution of electrons we estimate $n_{B}/n_{\gamma}\simeq 1.81\eta$. Hence the energy extracted from the CMB is estimated to be

[TABLE]

in very good agreement with the estimation performed in Khatri et al. (2012, §3), previously published in Chluba and Sunyaev (2012). In fact only the contribution above roughly $z=10^{5}$ has the time to generate a $\mu$-type distortion, as energy extracted at lower redshift is mainly in the form of $y$-type distortions. Hence the energy extracted and which results in a $\mu$-type distortion is reduced to $\delta\rho_{\gamma}/\rho_{\gamma}\simeq-1.8\times 10^{-9}$. This translates into

[TABLE]

where we used $\mu/T_{\gamma}\simeq 1.4\delta\rho_{\gamma}/\rho_{\gamma}$ (Hu, 1995, Eq. 3.48). It should not be surprising that energy extraction by Compton cooling results in a negative chemical potential, since it is also accompanied by a small decrease of the photon temperature on top of the cosmological redshifting $\propto 1/a$, and this chemical potential ensures the conservation of photons.

VII Anisotropic distribution functions

The Kompaneets equation (229) takes a simple form because of the high symmetries of an isotropic distribution. However the general form of the collision term at order $\epsilon^{2}$ in the Fokker-Planck expansion is much more involved as we detail in this section, since it involves the various angular moments of the distribution.

VII.1 Choice of frame

Given that we know how the distribution transforms from one frame (a tetrad) to another frame, we are free to choose the one which is the most adapted to the description of the photon spectrum. It is even possible to use the boost operator approach of Dai and Chluba (2014) to relate harmonic space multipoles of different frames at all orders in the boost velocity. It is customary in cosmology to consider the cosmological frame, which is defined by the fact that the time-like vector of the tetrad ${\bm{e}}_{0}$, is normal to constant coordinate time hyper-surfaces (Pitrou, 2009a).

However it appears that this choice is not optimal, and it is simpler and physically more transparent in many situations to work in the baryons frame, that is the one in which the bulk velocity of electrons vanishes. Of course this statement is arbitrary and one might prefer to work with the cosmological frame. Let us list the advantages of the baryon frame.

• •

The collision term remains at most quadratic in the distribution function, and even linear if we keep only the Thomson term. This would remain true even if we were to perform the Fokker-Planck expansion up to an arbitrary large power of $\epsilon$. However, if we were to transform the result to the cosmological frame, one would need to add terms with arbitrary high powers of the baryons velocity. One would thus needs another criteria to cut this expansion. This is discussed in appendix E.

• •

The variations of the electrons temperature is better analyzed in their frame. Since the energy transferred to electrons is deduced from the Compton collision term, it is easier to work directly in the baryons frame in order to have the energy transfer directly in the needed frame.

• •

The number density of electrons appears in Compton collisions, and it is much more natural to define this density in the frame of electrons, in order to avoid extra Lorentz factors.

• •

The effect of the baryons velocity cannot be fully removed by a change of frame. Once the bulk velocity is removed from the collision term, it reappears as a modification of the ${\rm d}\ln E/{\rm d}t$ in the Liouville term, that is it changes how free-streaming affects the energy of photons. However this can be seen as the effect of two boosts, one being at the last scattering event of the photon, and one at its reception. This is seen evidently in the line-of-sight reformulation (Zaldarriaga and Seljak, 1997; Hu and White, 1997) of the Boltzmann equation. Furthermore, a boost only aberrates the distribution, that is it changes directions, and shifts all $\ln E$ by the same quantity. Hence, if the spectrum is seen as a function of $\ln E$ instead of $E$, the effect of the boost related to the change of frame is simple. This idea is at the basis of our decomposition of spectra presented in § VIII, for which a shift in $\ln E$ is not a distortion, but simply a change of temperature.

• •

Eventually, working in the baryon frame allows for a neat separation between collisional effects, and special and general relativistic effects. What is unusual is not using the baryon frame, but instead using a different frame as it obfuscates this separation.

• •

All cosmological predictions are given up to a boost due to our peculiar velocity with respect to the cosmological frame, so we should always use the one which makes computations easier, or possibly the one which is more closely related to our peculiar velocity.

• •

Working in the baryons frame brings results which are valid up to order $\epsilon^{3}$ even though we performed an expansion up to order $\epsilon^{2}$, because there can be no odd powers. In the baryon frame one would need to keep terms which are cubic in the baryons velocity.

• •

More generally there is an inflation of terms when not working in the baryons frame, because they all appear as the effect of a boost on a simpler collision term computed in the baryon frame (Pitrou, 2009a; Chluba et al., 2012).

• •

Finally, it has also been argued in Chluba et al. (2012) that the optical depth (228) is better defined in the baryon frame, so as to simplify the interpretation of the SZ signal.

If we make such an important case for working in the baryon frame, it is because so far all CMB theoretical computations are formulated with the photon distributions considered in the cosmological frame. It appeared natural since numerical integrations are also associated to the cosmological frame. Describing the photon spectrum in the baryon frame, while computing its evolution following the cosmological time, would not seem a natural choice at first. This reformulation would not require much work for linear perturbation theory, as e.g. Ma and Bertschinger (1995); Hu and White (1997); Zaldarriaga and Seljak (1997). It would also be useless as our case of using the baryon frame is crucial only to describe the spectrum correctly, and at first order in cosmological perturbations no spectral distortions can arise, but only temperature variations. It is rather easy to see that once expressed with the line-of-sight method, the modification of the energy evolution ${\rm d}\ln E/\ln t\to{\rm d}\ln E/\ln t-{\rm d}(v^{i}n_{i})/{\rm d}t$ would give the same contribution as a collision term expressed in the cosmological frame, because perturbative effects at the observer’s position are always ignored.

However, using the baryon frame to describe the photon spectrum would require to revise more substantially the literature (and its associated numerical codes) on the second-order Boltzmann equation, among which Beneke and Fidler (2010); Huang and Vernizzi (2013, 2014); Pettinari et al. (2013); Su et al. (2012) but also Pitrou (2009a); Pitrou et al. (2010b). For completeness, we give in appendix E the collision term in the baryon frame, but only up to second order in cosmological perturbations.

VII.2 Thomson term

The lowest order terms in the Fokker-Planck expansion are those of order $\epsilon^{0}$. For an isotropic and unpolarized spectrum, they vanish. However, as soon as we allow for an angular structure in the spectrum, they lead to the Thomson contribution to the collision term. They are compactly written in the form

[TABLE]

The first term is the effect of scattering out events and is thus proportional to the distribution function itself. The second term accounts for the scattering in events. Its general expression is

[TABLE]

where the $\cal T$ indicates that the projected traceless part must be taken thanks to the projector (63). We have clearly separated the intensity, linear and polarization parts in a decomposition similar to Eq. (61) and we note them $I_{\cal R}[f]$, ${\cal P}_{\cal R}[f]_{ij}$ and $V_{\cal R}[f]$. With this notation, the various components of the Thomson collision term are

[TABLE]

VII.3 Thermal effects

The thermal effects are of order $\epsilon^{2}$ and arise from the integration (219). They also take a compact form with appropriate definitions. We find

[TABLE]

where we defined

[TABLE]

[TABLE]

Again, we have written these expressions such that their intensity ($I_{\cal S}[f]$ and $I_{\cal U}[f]$), circular polarization ($V_{\cal S}[f]$ and $V_{\cal U}[f]$), and linear polarization (${\cal P}_{\cal S}[f]_{ij}$ and ${\cal P}_{\cal U}[f]_{ij}$) components can be read directly. Hence the components of the collision term due to thermal effects are

[TABLE]

Eq. (248) matches exactly the thermal terms inside Eq. (C19) of Chluba et al. (2012) when ignoring linear and circular polarization.

VII.4 Recoil effects

The contribution of the electron recoil to the collision term also takes a very compact form and requires no further definition. It reads simply as

[TABLE]

The last term corresponds to a reduction of the scattering out contribution. Contrary to the other contributions, there are both linear and quadratic terms in the distribution functions due to the stimulated emission effects. The components of the recoil term are easily read. The intensity component is

[TABLE]

and it matches the recoil terms inside Eq. (C19) of Chluba et al. (2012), when ignoring linear and circular polarization. The polarization components of the recoil term are

[TABLE]

[TABLE]

We did not yet open explicitly the terms which are quadratic in the distribution function and which were on the second lines of the three previous equations. These purely quadratic terms are

[TABLE]

[TABLE]

[TABLE]

We remark that linear polarization does not affect circular polarization in the recoil term, nor does circular polarization affect linear polarization.

In practice, since the recoil terms are already small, since they are reduced by $E/m$, they might be linearized around an isotropic background distribution. In that case only the contributions from the first term of Eq. (VII.4) or the two first terms of Eq. (VII.4) survive, and the quadratic contributions reduce to

[TABLE]

[TABLE]

[TABLE]

In Pitrou (2009a, Eq. 6.18), a linearization of the effects coming from stimulated emission was used, hence matching only these linearized quadratic terms.

The thermal contributions of § VII.3 and the recoil terms of this section, once added, constitute the generalized Kompaneets equation, whose derivation is original in the case of polarized radiation, and extends Eq. (C19) of Chluba et al. (2012). Combined with the Thomson term of § VII.2, it rules the thermalization of an anisotropic and polarized distribution over an electron distribution. The total collision term, valid up to order $\epsilon^{2}$ is the sum of the Thomson contributions and the extended Kompaneets equations, that is

[TABLE]

Since we worked in the baryon frame, there are in fact no contribution at order $\epsilon^{3}$ so this collision term is only corrected by order $\epsilon^{4}$ contributions which have factors $ET_{e}/m^{2}$, $E^{2}/m^{2}$ or $T_{e}^{2}/m^{2}$. The decomposition of the angular dependence in spherical harmonics rather than in STF tensor is easily obtained, especially if the quadratic terms are linearized. Indeed, in that case one needs only the relations (D.2) and (340) to express the result with spherical harmonics.

The modifications of the Thomson contribution when considered in a general frame instead of the baryon frame are gathered for completeness in appendix E, even though as argued in § VII.1 it is preferable to work in the baryon frame. When decomposing the result in spherical harmonics, the procedure is much more involved and one must use various relations of § (D.4), yet another reason for not working in the baryon frame.

Finally let us comment that if circular polarization is initially vanishing, it is not generated by Compton collisions since $V_{\cal R}[f]$, $V_{\cal S}[f]$ and $V_{\cal U}[f]$ depend only on the circular polarization multipoles, and the quadratic terms in the recoil term (255) are linear in circular polarization. It is therefore customary to ignore circular polarization, unless in contexts where Faraday rotation sources it from linear polarization thanks to birefringence [see e.g. Montero-Camacho and Hirata (2018); Kamionkowski (2018)]. In such a case, the circular part of the collision term presented in this part should be used to describe properly its subsequent evolution under Compton scattering.

Spectral distortions

The angular dependence of the distribution function is expanded in moments, either with STF tensors or with spherical harmonics. However, it would be convenient to find an expansion of the dependence in the photon energies $E$ in suitable moments so as to reduce the number of degrees of freedom. In the next section we review the proposition of Stebbins (2007); Pitrou and Stebbins (2014) for such a decomposition, and we argue that when restricting to the Thomson terms, only the first few spectral moments are necessary to describe the spectrum. The dynamical evolution of the spectral moments deduced from the Thomson collision term is detailed in § IX.

VIII Spectrum parameterization

VIII.1 Distribution of Planck spectra

VIII.1.1 Temperature transform

Restricting first to unpolarized radiation, the distribution of photons is characterized only by its intensity $I$. It is a function of the position in space-time, the direction of propagation ${\bm{n}}$ and the energy $E$ of radiation, hence it is of the form $I(E,\dots)$, where dots indicate all the non-spectral dependence which we omit in most cases. In previous literature Zel’dovich et al. (1972); Chan and Jones (1975); Salas (1992); Chluba and Sunyaev (2004) the starting point for the description of the spectral dependence is to consider that $I$ is a superposition of Planck spectra with different temperatures, given by the distribution $p(T,\dots)$, such that

[TABLE]

with ${\cal B}(x)\equiv 1/(\exp(x)-1)$. If $\int_{0}^{\infty}p(T)dT\neq 1$ the distribution is said to be “gray”. Stebbins (2007) gives a full treatment of grayness and there it is shown that an initially non-gray distribution with only Compton-type interactions will remain non-gray. Henceforth we consider only non-gray distributions (see Ellis et al. (2013) for an example of a process inducing grayness). One can characterize the shape of the spectrum by the moments of the distribution $p(T)$. One thus defines

[TABLE]

Different authors have concentrated on following only specific moments. The most commonly used are the Rayleigh-Jeans temperature (which usually coincides with the electron temperature), $\bar{T}_{\rm RJ}\equiv\bar{T}_{(1)}$ Chluba and Sunyaev (2004); the number density temperature, $\bar{T}_{\rm n}\equiv\bar{T}_{(3)}$ Pitrou et al. (2010a); Naruko et al. (2013); Renaux-Petel et al. (2014); and the bolometric temperature $\bar{T}_{\rm b}\equiv\bar{T}_{(4)}$ Pitrou et al. (2010b); Creminelli et al. (2011); Huang and Vernizzi (2013); giving respectively the low frequency brightness, the number density of photons, and the energy density in photons. Indeed, using Eq. (260) we find

[TABLE]

for $p\geq 2$. Note that if the distribution function has a chemical potential, as in the case of a general Bose-Einstein distribution, the low energy limit is then a constant ($[\exp(\mu/T)-1]^{-1}$), and it is thus impossible to describe such distribution as a superposition of Planck spectra like in Eq. (260) whose low energy limit is $\propto\bar{T}_{\rm RJ}/E$.

However in Stebbins (2007) an alternative description of the spectrum based on different moments was proposed. At the basis of the formalism, is the use of the variable ${\cal T}\equiv\ln T$ (where a reference unit of temperature is implicit) whose distribution is $q({\cal T})\equiv Tp(T)$. The logarithmically averaged temperature is then simply defined by

[TABLE]

VIII.1.2 Spectral moments

The spectral distortions are characterized by the logarithmically averaged moments (LAM) of $q({\cal T})$: the moments about 0, $\{\eta_{p}\}$; the central moments, $\{u_{p}\}$; and the moments about a reference temperature, $\{d_{p}\}$, i.e

[TABLE]

where ${\cal T}_{0}\equiv\ln T_{0}$ and $T_{0}$ is an arbitrary reference temperature, usually chosen close to the mean. By construction, $u_{1}=0$, and since the spectrum is non-gray $\eta_{0}=d_{0}=u_{0}=1$.

Using ${\cal T}=({\cal T}-\bar{\cal T})+\bar{\cal T}$ and ${\cal T}=({\cal T}-{\cal T}_{0})+{\cal T}_{0}$, the moments (264) are related by Leibniz-type relations

[TABLE]

where

[TABLE]

The meaning of the moments is clear as one can reconstruct the spectrum by

[TABLE]

where

[TABLE]

Thus $\{d_{m}\}$ and $\{u_{m}\}$ are the coefficients of a generalized Fokker-Planck expansion around $T_{0}$ and $\bar{T}$, respectively. The $u_{p}$ are frame independent, but this is not the case for the other types of moments Stebbins (2007). The observed spectrum as a function of frequency and direction requires knowledge of the observer frame because of the Doppler effect and associated aberration, so one must also know $\bar{\cal T}$ and thus

[TABLE]

This is the first moment and it is directly related to the “temperature relative perturbation” which is $\exp(d_{1})-1$ since Eq. (270) is also

[TABLE]

The second moment gives the Compton $y$ distortion,

[TABLE]

These two moments are the ones most relevant for current observations.

VIII.1.3 Spectral moment of a polarized spectrum

Linear polarization will be generated by Compton scattering and the previous formalism can be extended to describe the polarization spectrum Stebbins (2007); Pitrou and Stebbins (2014), as we only need to consider the tensor-valued distribution function $f_{ij}(E,\dots)$ to describe both intensity and polarization. In practice, one would only consider the components in the two-dimensional sub-space, transverse to the photon direction and the observer velocity, that is one would use the $2\times 2$ matrix (64) noted $f_{ab}$, where the indices $a,b$ refer to a basis in this subspace. For simplicity we ignore circular polarization so $f_{ab}$ reduces to a trace ($I$) and a STF part $f_{\langle ab\rangle}$ for linear polarization.

A tensor-valued distribution of Planck spectra $q^{ab}({\cal T})$ is defined by

[TABLE]

and its matrix-valued moments $\{d^{ab}_{p}\},\{\eta^{ab}_{p}\},\{u^{ab}_{p}\}$ can be generalized from the $\{d_{p}\},\{\eta_{p}\},\{u_{p}\}$, for which a trace and a STF part can be defined. The relations (265,266) are then straightforwardly extended for linear polarization.

From the structure of the Compton collision term, it can be shown Stebbins (2007) that $d_{0}^{{\langle ab\rangle}}=\eta_{0}^{{\langle ab\rangle}}=0$ if initially so, but $u_{1}^{{\langle ab\rangle}}\neq 0$. The set of variables for the polarized part is thus simply the set of $\{u_{p}^{\langle ab\rangle}\}_{p\geq 1}$ as they are frame independent. Compared to the intensity, the main difference is that there is no temperature to be defined for polarization, but there is the non-vanishing moment $u_{1}^{\langle ab\rangle}$ which is the dominant one. A common misstatement or misunderstanding consists in treating this moment as a temperature perturbation, and to use the definition $\Theta^{\langle ab\rangle}\equiv u_{1}^{\langle ab\rangle}$, but strictly speaking, it is a pure spectral distortion, and as such frame independent. In Naruko et al. (2013), it is called the “temperature part” of the polarization, as opposed to the primary spectral distortion

[TABLE]

VIII.1.4 Discussion on the choice of a set of variables

It is clear that since the $\{u_{p}\}_{p\geq 2}$ are frame invariant they are good candidates to describe the spectral distortions. The use of $d_{1}$ for the temperature perturbation is then natural as it fits into this formalism. However, one might wonder if this is the only set of variables with such appealing properties. Starting from the moments defined in (261), we can relate these to the $\{d_{p}\}$ and $\{u_{p}\}$ by

[TABLE]

It appears clearly that, for a given $p$, the temperature $\bar{T}_{(p)}$ can be used to define a temperature perturbation and the moments

[TABLE]

To illustrate how these temperature perturbations are related to $d_{1}$, let us keep only the moments $m\leq 2$ to express the bolometric temperature perturbation $\Theta_{\rm b}$ ($p=4$), and the number density temperature $\Theta_{n}$ ($p=3$). We find they are related to $(d_{1},d_{2})$ or $(d_{1},u_{2})$ by

[TABLE]

and in particular

[TABLE]

The $\{M_{(p),m}\}_{m\geq 2}$ would be as good as the $\{u_{m}\}_{m\geq 2}$ to describe the spectral distortions, since they are obviously frame invariant as they involve only an (infinite) sum of products of the $\{u_{p}\}$. In the next section, we argue that to decide which set of variables should be used, one should examine the dynamical evolution, and choose the one which has the simplest structure, and for which numerical integration is simplified.

VIII.2 Spectral moments evolution

VIII.2.1 General form of the Boltzmann equation

The general form of the Boltzmann equation is (again we omit the dependence in $(E,\dots)$ for brevity)

[TABLE]

where the convective derivative ${\cal D}/{\cal D}t$ acts on all the dependence except the spectral dependence, and accounts for the effect of free streaming.

The collision term can also be described by its moments $\{\eta_{p}^{C,\,ab}\},\{u_{p}^{C,\,ab}\},\{d_{p}^{C,\,ab}\}$ which are related by relations similar to (265) and (266), that is

[TABLE]

In order to find the evolution of the $\{u_{p}^{ab}\}$, it proves simpler to first derive from (279) the evolution of the $\{d_{p}^{ab}\}$, and we get

[TABLE]

So for the temperature perturbation ($d_{1}$), the trace of $m=1$ gives

[TABLE]

If the spectrum is initially non-gray, and radiation is only subject to Compton scattering, it remains so and this property translates to $d_{0}^{C,\,ab}=0$. In that case, this implies

[TABLE]

and the first $u^{C,\,ab}_{m}$ are related to the first $d^{C,\,ab}_{m}$ by

[TABLE]

The moments $\{d_{p}^{C,ab}\}$ can be read off the collision term as long as we do not consider recoil terms. Indeed the Thomson and thermal terms involve only $\partial^{n}_{\ln E}$, which is exactly what is used in the expansion (268). In fact if thermal effects are ignored, and if we work in the baryon frame, the Thomson term is extremely simple, whereas one must use the general frame expression of appendix E if not working in the baryon frame. What is crucial is that the moments $\{d_{p}^{C,ab}\}$ are linear in the variables $\{d^{ab}_{p}\}$ which describe the radiation spectrum. However, they still couple non-linearly to the baryons bulk velocity if one insists in not working in the baryon frame.

From the relations (280), one infers that

[TABLE]

This system of equation is closed at any order $p$, since the equation-of-motion for $u^{ab}_{p}$ depends only on $u^{ab}_{p^{\prime}}$ for $p^{\prime}\leq p$. One can truncate this system of equations at any order, but one must bear in mind that we have neglected recoil terms.

It was crucial in these derivations that ${\rm d}\ln E/{\rm d}t$ does not depend on $E$ but only on the metric and the direction of propagation ${\bm{n}}$. Seen as a function of $\ln E$ instead of $E$, the spectrum is only shifted by free-streaming but the overall shape remains unchanged. Since the temperature transform (260) depends on $E/T$, that is on $\ln E-\ln T$, this property of global shifting is transferred to the distribution $q({\cal T})$ of superimposed Planck spectra. Centered moments are thus very well adapted since only the center ($d_{1}$) is affected by a global shift in Eq. (283), but not the centered moments $u_{n}^{ab}$ in Eq. (286). The structure is exactly similar for the effect of a boost, if we let aside the aberration which affects directions. Indeed it also shifts $\ln E$ by a constant quantity, and therefore $d_{1}$ is affected by a boost but not the centered moments $u^{ab}_{n}$.

VIII.2.2 Doppler, SZ effect and $y$-type distortion

At first order one needs only the temperature perturbation $d_{1}$ and $u_{1}^{\langle ab\rangle}=d_{1}^{\langle ab\rangle}$. At second order, one adds the spectral distortions $u_{2}$ and $u_{2}^{\langle ab\rangle}$, and this distortion, known in this context as the non-linear kinetic SZ effect Pitrou et al. (2010a); Renaux-Petel et al. (2014), is generated by the r.h.s. of (286) with $p=2$.

The distortion generated by the thermal SZ effect Zeldovich and Sunyaev (1969) is also captured by $u_{2}$ and the usual $y$ parameter associated with it is related by the relation (272). Note that it is apparent on Eq. (232) when compared to the expansion (268) that the thermal SZ effect also induces a shift in $d_{1}$, but we also check from Eq. (277) that it does not affect $\Theta_{\rm n}$ since Compton collisions conserve the number of photons.

A polarized $y$-type distortion can also be defined Sunyaev and Zeldovich (1980); Naruko et al. (2013); Renaux-Petel et al. (2014) and is related to the moments by

[TABLE]

VIII.2.3 Structure of the numerics

Eq. (286) shows that

spectral distortions are affected only by the collision term, as they remain unaffected by metric perturbations [see also Stebbins (2007); Pitrou et al. (2010a); Naruko et al. (2013)]; 2. 2.

metric perturbations, which enter through the redshifting term ${\rm d}\ln E/{\rm d}t$ affect only the evolution of the temperature perturbation $d_{1}$, and more importantly do not couple non-linearly with $d_{1}$ [Eq. (283)]; 3. 3.

the collision term for the evolution of $u_{p}^{ab}$ [the r.h.s of (286)], contains only terms of the form $d_{1}^{p-k}u_{k}^{ab}$ with $k\leq p$ (see Stebbins (2007) for more details) multiplied by powers of the baryons bulk velocity. Therefore it restricts the non-linearities to products of at most $p$ factors of spectral moments, when considering the evolution of the moment of order $p$. N.B. for $p=1$ the collision term ($d_{1}^{C,ab}$) is linear in the moments.

Any other parameterization of the distortion based on the $M_{(p),n}$ defined in (276) would conserve property (1). However, property (3) would be lost with the $M_{(p),n}$. The loss of this property is, in principle, not a serious problem for the numerical integration, since interactions are localized in time by the visibility function. However, this would lead to unnecessary complications when going to higher orders of perturbations and thus higher moments. Our first argument here is that the simplest is the best.

Our second argument is that property (2) is crucial for the numerical integration since redshifting effects are not localized in time. Indeed, by avoiding a non-linear coupling between the temperature perturbations and the metric perturbations, the numerical integration is made possible even at the non-linear level as it avoids coupling between the angular moments of the temperature perturbations with the metric perturbation Huang and Vernizzi (2013). Finding a form of the Boltzmann equation that satisfies this property, was the key to a successful numerical integration at second order Huang and Vernizzi (2013); Pettinari et al. (2013). With the present formalism, this property arises naturally for the variable $d_{1}$. Metric perturbations would also affect the geodesic and lead to time-delay and lensing effects, but these can be treated separately Hu and Cooray (2001); Huang and Vernizzi (2014). There would be of course other variables for which property (2) holds. For instance, defining

[TABLE]

one obtains from (275) that the variables

[TABLE]

obviously satisfy property (2) but not property (3). Up to second order in cosmological perturbations (neglecting $\{u_{p}\}_{p\geq 3}$) the definitions for the most common temperatures are related by

[TABLE]

This motivated the use of $\widetilde{\Theta}_{\rm b}$ instead of $\Theta_{\rm b}$ in the final output of Huang and Vernizzi (2013), since property (2) is satisfied for the former and not for the latter.

Similarly, for the fractional perturbation to the energy density, one finds up to second order in cosmological perturbations

[TABLE]

and using $\tilde{\Delta}\equiv\ln(1+\Delta)$, we find

[TABLE]

Again this motivated the use of $\tilde{\Delta}$ instead of $\Delta$ in the intermediate numerics of Huang and Vernizzi (2013), so as to keep property (2) satisfied. A final example can be made with the fractional energy density perturbation of linear polarization. One finds

[TABLE]

and the non-linear term $d_{1}^{\langle ab\rangle}d_{1}$ will induce a non-linear coupling of the type $d_{1}^{\langle ab\rangle}{\rm d}\ln E/{\rm d}t$ in the evolution equation of $\Delta^{\langle ab\rangle}$. However, using

[TABLE]

this non-linear coupling disappears Fidler et al. (2014) and property (2) is recovered. In all these three examples, property (2) can be restored with an ad-hoc change of variable, but property (3) is not satisfied, due to the term in $u_{2}$ for the first two examples, and due to the term $u_{2}^{\langle ab\rangle}$ for the last one. It implies in particular that the evolution equation for the lowest order moment in this description, i.e. their temperature perturbation, has a collision term which is not linear in the moments of radiation.

VIII.3 Summary and notation

The essential properties described above for the structure of dynamical equations are only met with the set of variables made of $d_{1}$, $\{u_{p}\}_{p\geq 2}$ and $\{u_{p}^{\langle ab\rangle}\}_{p\geq 1}$. Furthermore, the moments which characterize the spectral distortions are frame independent and thus do not depend on our local velocity. Only the angular dependence is affected by the choice of frame due to aberration effects. We strongly recommend that these moments should be used to parameterize the CMB spectrum when recoil effects are neglected. However we shall use names which are more reminiscent of temperature and for the next section we define

[TABLE]

We must remember that the temperature associated with intensity is recovered from Eq. (271), and that $\theta^{{\langle ab\rangle}}$ is strictly speaking not a temperature, but rather the lowest spectral distortion of polarization. Similarly we also do not work directly with $u_{2}$ nor $u_{2}^{{\langle ab\rangle}}$ for the spectral distortions, but rather with their halves, the $y$ and $Y^{{\langle ab\rangle}}$ variables defined in Eqs. (272) and (287). We also restore spatial tetrad indices $i,j\dots$ instead of indices $a,b$ referring two the screen projected space, that is we use $\theta_{ij}$ and $Y_{ij}$. We now restrict the expansion (268) of the spectrum to these moments, that is to second order effects only. Hence our spectrum parameterization is

[TABLE]

where the argument of the Planck spectrum is ${\cal B}(E/T_{0})$ and with the logarithmic derivative

[TABLE]

Note that the temperature defined in Pitrou et al. (2010a); Naruko et al. (2013) is exactly $\Theta_{n}$. The expansions (2.24) and (2.26) of Naruko et al. (2013) appear at first sight different, but using the relation (277) we check that they agree with Eqs. (296).

Finally, the moments of $\Theta$ and $y$ in an expansion of the type (68) are noted $\Theta_{L}$ and $y_{L}$. As for $\theta^{ij}$ and $Y^{ij}$, we decompose them as in Eq. (81), that is in $E$-modes (noted $E^{\theta}_{L}$ and $E^{Y}_{L}$) and $B$-modes (noted $B^{\theta}_{L}$ and $B^{Y}_{L}$).

IX Angular correlations of spectral distortions

IX.1 Collision term of spectral moments

We consider only the effect of the Thomson collision term and we work in the baryon frame. One should use the collision terms of appendix (E) if one wishes to rephrase these results in a general frame. However we argued that distortions are better analyzed in the baryon frame. This description with only the Thomson terms applies to the reionization epoch or around recombination when dealing with the dissipation of baryon acoustic oscillations, since as a first approximation the effect of the extended Kompaneets collision term can be ignored on anisotropies Chluba et al. (2012).

From the discussion in § VIII.2 the evolution of $\Theta$ and $\theta_{ij}$ is extremely simple. It is just read from the Thomson term exposed in § VII.2 where we use the replacement rules

[TABLE]

and the same rules for the associated STF multipoles. For completeness we repeat the result here which is

[TABLE]

There is no term quadratic in the spectrum entering the collision term for $\Theta$ nor $\theta_{{\langle ij\rangle}}$ as already stressed in § VIII.2.3! In Naruko et al. (2013) it is found that quadratic terms arise, but it is only because of the use of $\Theta_{n}$. Furthermore, red-shifting effects (encoded formally by ${\rm d}\ln E/{\rm d}t$) do not affect $\theta_{ij}$, as seen on Eq. (284). With our definitions, the red-shifting of energy only affects $\Theta$.

This shows once more the importance of choosing appropriate variables. Of course, one must not forget that non-linearities arising from metric perturbations enter the evolution of $\Theta$ through ${\rm d}\ln E/{\rm d}t$, or from the transport operator ${\cal D}/{\cal D}t$, and eventually the temperature is obtained by the non-linear relation (271).

We now turn to the collision terms of $y$ and $Y_{ij}$. These contain quadratic terms, even when restricting to the Thomson collision contribution in the baryon frame, and this can be seen from the relation  (285).

From results of §VIII.2, and using Eq. (285) we find

[TABLE]

where the notation $[\dots]_{\emptyset}$ means that we must extract the monopole when decomposing the angular dependence of the expression in square brackets in STF tensors. Similarly $[\dots]_{ij}$ indicates we must extract the quadrupole of the scalar quantity inside the brackets and $E[\dots]_{ij}$ that we must extract the electric type quadrupole of the tensorial quantity inside brackets. Physically, Thomson scattering remaps directions and thus mixes Planck spectra of different temperatures if the distribution is not isotropic, yielding a y-type distortion at lowest order.

The collision term associated with $Y_{ij}$ is

[TABLE]

Note that the dissipation of baryon acoustic oscillations originates when the temperature $\Theta$ ceases to be equal to its monopole $\Theta_{\emptyset}$, hence feeding the evolution of $y$ in Eq. (300) but also sourcing the distortion of polarization through Eq. (301). A tight-coupling expansion (Pitrou, 2011) of Eqs. (IX.1) and (299) allows to obtain $\Theta-\Theta_{0}$ but also the quadrupoles which can be used to estimate the effect (Chluba et al., 2012).

IX.2 Non-linear kSZ effect during reionization

After recombination, that is below $z\simeq 10^{3}$, baryons start to decouple from photons and the velocity difference between them and photons starts to grow. Since we performed computations in the baryons frame, this velocity difference is hidden in the dipole of the distribution. Hence we define

[TABLE]

and we assume that it is sufficient to characterize the angular dependence of temperature ($\Theta\simeq\Theta_{\emptyset}+{\cal V}_{i}n^{i}$). During the reionization era, ${\cal V}_{i}$ is growing because matter collapses whereas radiation free-streams. The collision terms of the previous section considerably simplify as they reduce to

[TABLE]

The velocity difference ${\cal V}_{i}$ sources the monopole and the quadrupole of the $y$ spectral distortion via quadratic terms. As for the distortion in polarization $Y_{ij}$, its quadrupolar electric type multipole is also sourced by quadratic terms in the velocity. This effect is nothing but the non-linear kinetic Sunyaev Zel’dovich (kSZ) effect. The angular correlations of $y$ but also of the $E$ and $B$ modes of $Y_{ij}$ generated during the recombination era due to the large scale velocity of baryons in the intergalactic medium has been computed in Pitrou et al. (2010a); Renaux-Petel et al. (2014) with the line-of-sight method, and we reproduce the figures. In Fig. 5 we plot the $E$-type multipoles (the $C_{\ell}$’s) associated with the temperature-like signal ($\theta_{ij}$) along with those from the distortion ($Y_{ij}$). Then in Fig. 6 we also compare the $B$-type multipoles of $\theta_{ij}$ and $Y_{ij}$, but also those arising from primordial gravitational waves with tensor-to-scalar ratio $r=0.001$. Finally $E$-type multipoles of distortions and its correlations with $y$-type distortions are plotted in Fig. 7. For the $y$-distortion the effect should be subdominant compared to the thermal $y$-distortions from all unresolved clusters which has been already detected Aghanim et al. (2016); Hill et al. (2014). However there is no thermal counterpart for the distortion in polarization $Y_{ij}$.

Conclusion

We have emphasized the similarities in the construction of distribution functions for fermions and bosons. While polarization for fermions is naturally described by a vector, we need a tensor to describe the polarization state of photons. In the case of massless fermions, there are even more similarities since linear polarization is described in a screen-projected space and circular polarization is defined separately, exactly as for a gas of photons. We can then anticipate the description of a gas of gravitons when considered as spin-$2$ particles. Indeed, the stochastic background of gravitational waves [see e.g. Cusin et al. (2017); Cusin et al. (2018c, b)] could be considered as a gas whose transfer could be deduced from a Liouville equation. Once covariantized with the polarization tensors $\epsilon_{\pm}^{i}\epsilon_{\pm}^{j}$, its polarization state built as in Eq. (62) would be a $4$-index symmetric transverse (to both $p^{\mu}$ and $u^{\mu}$) traceless tensor, whose multipolar decomposition is performed as in Eq. (78) but with spin-$4$ spherical harmonics (Cusin et al., 2018a).

The derivation of the collision terms are also extremely similar between the Fermi theory of weak interactions and the effective $4$-vertex description of QED. They only differ in the statistical factor of the final distributions, which are Pauli-blocking factors for fermions in weak interactions but stimulated emission factors for photons in Compton scattering. Apart from this change of sign in the statistics, the structure of the effect of the final state distribution [Eq. (143)], which cannot be guessed a priori, is completely similar.

We then emphasized a third similarity in the treatment of weak interactions exchanging neutrons and protons before BBN, and Compton scattering. In both cases there is a massive particle in the initial and final state and the collision term can be computed with one-dimensional integrals using a Fokker-Planck expansion, that is an expansion in the momentum transferred to the massive particle. For the neutron-proton conversions, this allowed to compute the so-called finite nucleon mass corrections with a method previously introduced in Pitrou et al. (2018). For Compton scattering, it allows to obtain the thermal and recoil corrections to the lowest order approximation known as Thomson scattering. These corrections when considered for an anisotropic photon distribution are only consistent when polarization is included (since the quadrupole of the distribution generates linear polarization) and lead to the extended Kompaneets equation presented in § VII. It extends the results of Chluba et al. (2012, Eq. C19) which were derived in the anisotropic but unpolarized case. Furthermore, we argued that spectral distortions should be computed in the baryon frame even though any frame is in principle suitable since it is always possible to boost the distribution functions.

Finally, we discussed that by remapping directions, the Thomson collision term generates spectral distortions through the mixing of Planck spectra when the distribution is not isotropic. We argued that a parameterization based on logarithmic and centered moments of the distribution of Planck spectra should be preferred as it leads to a simple separation of the collision term into thermal and spectral contributions. Furthermore, apart from the variable describing temperature fluctuations, all spectral moments are frame invariant in the sense that only the directions are aberrated by a Lorentz transformation. With this parameterization, no quadratic term arises in the collision term governing the evolution of the quantity which characterizes temperature fluctuations. We summarized the equations governing the dissipation of baryon acoustic oscillations and the non-linear kinetic SZ effect, and we stressed that distortions exist also in polarization.

Acknowledgements.

It is a pleasure to thank Jean-Philippe Uzan, Francis Bernardeau, Thiago Pereira, Pierre Fleury, Sébastien Renaux-Petel, Guillaume Faye, Giulia Cusin, Christian Fidler, Atsushi Naruko and Julien Froustey for collaborations and continuous discussions on the topic. I thank particularly Guilhem Lavaux and Lucie Gastard for their long standing encouragements to complete this work.

Appendices

Appendix A Spinor valued operators

A spinor valued operator has $16$ degrees of freedom and we thus need a $16$-dimensional basis to decompose operators in spinor space. The set (33) is a complete basis for the space of operators in spinor space, where we defined the matrices

[TABLE]

The operators are orthogonal and we find for any two different operators $X$ and $Y$(${X_{\mathfrak{a}}}^{\mathfrak{b}}$ and ${Y_{\mathfrak{a}}}^{\mathfrak{b}}$ in spinor components) in the set ${\cal O}$ that ${\rm Tr}[X\cdot Y]={X_{\mathfrak{a}}}^{\mathfrak{b}}{Y_{\mathfrak{b}}}^{\mathfrak{a}}=0$. Using this property, any operator can be decomposed onto this basis with the help of the Fierz identity Nishi (2005) (see also Eq. G.1.99 of Dreiner et al. (2010) taking into account a factor $2$ difference in the definition of $\Sigma^{\mu\nu}$)

[TABLE]

The last equality defines the coefficients $c_{X}$ of the expansion which, by construction, satisfy $c_{X}=1/{\rm Tr}[X\cdot X]$. Note also that this identity can be used with the matrices $\widetilde{\Sigma}^{\mu\nu}$ instead of $\Sigma^{\mu\nu}$ by employing $(\Sigma^{\mu\nu})_{\mathfrak{a}}^{\,\,\mathfrak{d}}(\Sigma_{\mu\nu})_{\mathfrak{c}}^{\,\,\mathfrak{b}}=-(\widetilde{\Sigma}^{\mu\nu})_{\mathfrak{a}}^{\,\,\mathfrak{d}}(\widetilde{\Sigma}_{\mu\nu})_{\mathfrak{c}}^{\,\,\mathfrak{b}}$.

Any operator $O$ in spinor space is decomposed on the basis ${\cal O}$ thanks to Eq. (A) as

[TABLE]

In particular, any bilinear tensor product of the form $u_{r}(p)\bar{u}_{s}(p)$ or $v_{r}(p)\bar{v}_{s}(p)$ [with the standard notation $\bar{u}_{s}=u_{s}^{\dagger}\gamma^{0}$], is a spinor-space operator and can be decomposed as

[TABLE]

Hence, we only need to compute the $\bar{u}_{s}(p)Xu_{r}(p)$ and $\bar{v}_{s}(p)Xv_{r}(p)$ [these expressions are computed in Fidler and Pitrou (2017, App. D)] for all $X\in{\cal O}$ to decompose the operator ${F_{\mathfrak{a}}}^{\mathfrak{b}}$, defined in Eq. (31).

Appendix B Construction of the classical collision term

We define the operator $N_{rs}^{(0)}$ characterising the ingoing states prior to the collision. To first order in the interaction we obtain

[TABLE]

The interpretation of this equation is that at the time $t=0$ the system is starting to interact, but as the background does not yet contain any correlations between the interacting species we can still evaluate the collisions using the zeroth order number operator. This first order solution describes forward scatterings and we need to go to second order to find the first non-forward interactions.

We insert the first order solution (310) into Eq. (106) and find to second order

[TABLE]

The second order contribution describes an interaction which is active between the time $t^{\prime}$ and $t$. We identify this timescale with our microscopic timescale $t_{\rm mic}=t-t^{\prime}$, quantifying the timescale of individual particle interactions. The averaged fluid however does not change significantly on this timescale and evolves on the much larger mesoscopic time-scale $t_{\rm mes}\gg t_{\rm mic}$. Since we compute the derivative of the number operator $N$ with respect to the time $t$ we may identify the mesoscopic time $t_{\rm mes}=t$. Expressed in these parameters we obtain

[TABLE]

Ideally we would like to evaluate this equation at the initial time and set $t_{\rm mes}=0$ to compute the change of our initial states under the considered interactions. This choice however is mathematically inconsistent as we are mixing mesoscopic and microscopic timescales in the integration. Instead we average the resulting time-derivative, considering the time-reversal symmetry, over a box that is centered on the initial time and has a length of $2\epsilon$ which is chosen to be small compared to the scale of macroscopic evolution.

[TABLE]

We may split this integration into three regions. First, the central region $\epsilon\approx t_{\rm Compton}$, where $t_{\rm Compton}$ the typical Compton time scale of particles, is highly non-trivial, but this region is negligible compared to our entire integration volume. In the remaining positive and negative regions the integrand is constant in time. The reason is that the integral over the microscopic time already has sufficient support and is converged. The remaining time-dependence based on the mesoscopic time $t_{\rm mes}$ is not relevant as we have chosen the box small compared to the mesoscopic evolution and we may now set $t_{\rm mes}=0$ yielding

[TABLE]

Finally we may extend the integration limit $\epsilon$ to infinity compared to the microscopic evolution using a separation of scales.

We note that the interaction Hamiltonian appearing in this equation may always be evaluated based on the non-interacting field value as we only utilise times which are small compared to the mesoscopic time. Our expression is equivalent to those used in Sigl and Raffelt (1993); Kosowsky (1996); Beneke and Fidler (2010).

We finally deduce from Eq. (106) that the classical evolution of the distribution function is dictated by

[TABLE]

where the first term on the rhs is the forward scattering term. It is responsible for refractive effects or flavor oscillations in matter (see Lesgourgues and Pastor (2006, 2012) for neutrino oscillations in cosmology) such as the MSW effect Wolfenstein (1978); Mikheev and Smirnov (1986); Marciano and Parsa (2003); Sigl and Raffelt (1993); Volpe (2015)). The second term is the collision term and we define

[TABLE]

such that the Boltzmann equation (315) (when neglecting forward scattering and restoring the notation of the momentum dependence) is written

[TABLE]

Appendix C Finite nucleon mass corrections

It is not fully correct to consider that nucleons have an infinite mass. Indeed, the typical energy transfer in weak interactions to electrons and neutrinos is of the order of the mass gap $\Delta\simeq 1.29\,{\rm MeV}$, which is $1.4\times 10^{-3}$ smaller than the nucleon mass. It corresponds to a temperature of $1.5\times 10^{10}{\rm K}$ which is not much larger than the freeze-out temperature. In the infinite nucleon mass approximation, we have thus neglected factors of the type $E_{\nu}/m_{N}$, $E_{e}/m_{N}$ or $\Delta/m_{N}$ (where $m_{N}$ is the average nucleon mass $m_{N}\equiv(m_{p}+m_{n})/2$) which represent order $10^{-3}$ corrections with respect to the leading one around $10^{10}{\rm K}$ and even larger corrections at higher temperature. Our method consists in expanding the full reaction rate in power of a small parameter $\epsilon$ related to the momentum transfer. Given the relation between kinetic energy and momenta, $T/m_{N}$ is of order $\epsilon^{2}$. Terms of the type $E_{\nu}/m_{N}$ or $E_{e}/m_{N}$ are also of order $\epsilon^{2}$ and terms of the type $\Delta/m_{N}$ are also treated as being of order $\epsilon^{2}$. Our implementation of the finite mass corrections consists in including all the terms up to order $\epsilon^{2}$, but neglecting terms of higher order. This means that we neglect terms whose importance is of order $10^{-6}$.

If we ignore radiative corrections at null temperature, these corrections take the form

[TABLE]

and the functions $\chi^{\rm FM}_{\pm}$ are

[TABLE]

where $p=\sqrt{E^{2}-m_{e}^{2}}$, $E_{\nu}^{\mp}=E\mp\Delta$. We defined the reduced couplings

[TABLE]

and the functions [with the notation (187)]

[TABLE]

However, the finite nucleon mass corrections must be coupled with radiative corrections, and one must also account for the weak magnetism in the neutron/proton current as it is also of the same order as finite nucleon mass corrections. The full set of corrections is reported in Pitrou et al. (2018).

Appendix D Symmetric trace-free (STF) tensors

D.1 Notation

We introduce the multi-index notation

[TABLE]

and when no ambiguity can arise we use $L$ instead of $I_{\ell}$. When $\ell=0$ we use the notation $I_{\emptyset}$.

The symmetric trace-free part of a set of indices is noted $\langle\dots\rangle$ and it can be used with multi-index notation, e.g. $n^{\langle I_{\ell}\rangle}$. General formula for extracting symmetric and then traceless parts can be found in e.g. Thorne (1980); Blanchet and Damour (1986).

The angular integration of a product of direction vectors is

[TABLE]

D.2 Relation to spherical harmonics

Let us define for functions $A({\bm{n}})$ and $B({\bm{n}})$

[TABLE]

It is possible to obtain the orthogonality relations

[TABLE]

where

[TABLE]

Eq. (325) is a particular case of Faye et al. (2015, Eq. C2). Defining

[TABLE]

we can expand the directional dependence either on spherical harmonics or $n^{\langle L\rangle}$ using

[TABLE]

The inverse relation is

[TABLE]

From the closure relation

[TABLE]

we get the closure relation

[TABLE]

Explicitly the ${\cal Y}_{\ell m}^{L}$ are given by

[TABLE]

where

[TABLE]

Since we use a Cartesian or triad basis we also define ${\cal Y}_{\ell m}^{I_{\ell}}={\cal Y}^{\ell m}_{I_{\ell}}$ and we have property ${\cal Y}_{\ell m}^{\star\,I_{\ell}}=(-1)^{m}{\cal Y}_{\ell\,-m}^{I_{\ell}}$. The ${\cal Y}_{\ell m}^{L}$ satisfy the orthogonality property

[TABLE]

D.3 Relation to spin-weighted spherical harmonics

The ${\cal Y}_{\ell m}^{L}$ are also related to spin-weighted spherical harmonics. To that purpose, we use the polarization basis

[TABLE]

Let us define (for $s>0$) the compact notation (Pitrou and Pereira, 2019)

[TABLE]

which generalizes the products (322). For $s>0$ the spin-weighted spherical harmonics are related by

[TABLE]

where

[TABLE]

The relations (338) are inverted as

[TABLE]

D.4 Products and contractions of STF tensor

The symmetrized products of the ${\cal Y}_{\ell m}^{L}$ are directly related to the products of spherical harmonics. Indeed, it can be shown that 999This formula, though appearing first in Pitrou (2009a), has been derived by Guillaume Faye. It is obtained by contracting the l.h.s with $n^{I_{\ell_{1}+\ell_{2}}}$, and using Eq. (328) to recognize the products of spherical harmonics whose expressions in terms of Clebsch-Gordan coefficients is known. The formula is then recovered by taking $\partial_{I_{\ell_{1}+\ell_{2}}}$ with some algebraic manipulations.

[TABLE]

where $m_{3}=m_{1}+m_{2}$ and the sum runs only over even $\ell_{1}+\ell_{2}-\ell_{3}$. The Clebsch-Gordan coefficients are related to Wigner-3j symbols by

[TABLE]

with $\ell_{3}=\ell_{1}+\ell_{2}$ and $m_{3}=m_{1}+m_{2}$.

In particular, we deduce a relation similar to the Gaunt integral of three spherical harmonics, which is

[TABLE]

As a first application, it allows to obtain the symmetrized and trace-free product (still for $\ell_{3}=\ell_{1}+\ell_{2}$)

[TABLE]

When considering functions $A({\bm{n}})$ and $B({\bm{n}})$, which are expanded in spherical harmonics multipoles $A_{\ell m}$ and $B_{\ell m}$ as in Eq. (68) or in STF tensor $A_{L}$ and $B_{L}$ as in Eq. (69), this relation is the key to extract the multipoles (in terms of spherical harmonics) of a product of the type $A_{\langle I_{\ell_{1}}}B_{J_{\ell_{2}}\rangle}$ in terms of the $A_{\ell m}$ and $B_{\ell m}$. It allows to work entirely with STF tensors, and only convert to spherical harmonics multipoles at the very end if desired.

The second application of Eq. (343) are the contractions

[TABLE]

Again this allows to get the spherical harmonics multipoles of a product of the type $A_{I_{\ell_{1}}J_{\ell_{2}}}B^{J_{\ell_{2}}}$ in terms of the $A_{\ell m}$ and $B_{\ell m}$.

When the Levi-Civita is involved we must use

[TABLE]

which is deduced from (Thorne, 1980, Eqs. 2.26c-e)

[TABLE]

where

[TABLE]

Their explicit expression is (Pitrou, 2009a, Eq. 7.33)

[TABLE]

Finally let us report the useful relation (Blanchet and Damour, 1986, Eq. A.22a)

[TABLE]

We extend it for products of the type (337) and we find (Pitrou and Pereira, 2019, App. C.2)

[TABLE]

D.5 The ${\cal Y}_{\ell m}^{L}$ in the literature

The ${\cal Y}_{\ell m}^{L}$ are often present in the literature even though not under the present notation nor with the same normalizations. For instance, the $Q_{i}^{(m)}$ and $Q_{ij}^{(m)}$ of Hu and White (1997) and Pitrou (2009a) or the $\xi^{i}_{m}$ and $\chi_{2m}^{ij}$ of Beneke and Fidler (2010) are proportional to the ${\cal Y}^{1m}_{i}$ and ${\cal Y}^{2m}_{ij}$ when evaluated at vanishing radial distance. In particular, the relations (66) of Beneke and Fidler (2010) are particular cases of Eqs. (340).

Appendix E Collision term in a different frame

In this part we report the collision term when expressed in a general frame (that is not in the baryon frame) such that baryons possess a spatial bulk velocity ${\bm{v}}$. We consider that this bulk velocity is a factor of order $\epsilon$ in the Fokker-Planck expansion, even though this is not really the case as it was initially an expansion in the momentum transferred to the electron. But this choice allows a first set of simplifications, since we ignore the coupling of recoil and thermal terms with the baryon bulk velocity as these would be of order $T_{e}/m|{\bm{v}}|$ or $E/m|{\bm{v}}|$ and thus of order $\epsilon^{3}$. That is restricting to second order in $\epsilon^{2}$, one would use the thermal and recoil terms of § VII.3 and VII.4. Hence the modifications introduced by the bulk baryon velocity arise only from the Thomson term.

Furthermore, as expressions can still be very sizable, we also perform a secondary expansion in which all multipoles (except the monopole of intensity $I_{\emptyset}$) and the baryon bulk velocity, are considered as first order quantities. In the cosmological context this amounts to an expansion in cosmological perturbations. We choose to restrict to second order in this expansion.

The intensity, linear and circular polarization parts of the collision term arising from bulk baryon velocity in the Thomson term are summarized in the next sections, in which we use the short-hand notation (297) and we omit to write the dependence of the multipoles on $E$.

E.1 Intensity

The intensity part of the Thomson contribution in a general frame is

[TABLE]

The non-vanishing multipoles in the right hand side are

[TABLE]

[TABLE]

The full multipolar decomposition is then obtained by decomposing the first term in Eq. (369), using

[TABLE]

where in the left hand side is meant the STF components of $Iv_{i}n^{i}$ in an expansion of the type (69). This relation is easily shown using Eq. (367). Once Eq. (369) is integrated over $E^{3}{\rm d}E$ so as to get a collision term for the brightness only, we can check that we recover Eq. (6.24) of Pitrou (2009a).

E.2 Linear polarization

The linear polarization part of the Thomson contribution in a general frame is

[TABLE]

The non vanishing multipoles in the right hand side are

[TABLE]

[TABLE]

The full multipolar decomposition of the form (81) is then obtained by decomposing the first term of Eq. (E.2) using

[TABLE]

[TABLE]

where in the left hand sides are meant the STF components of $E$ and $B$ type of ${\cal P}_{ij}v_{k}n^{k}$ in an expansion of the type (81). The proof of these identities follows from the use of Eq. (368) with Eq. (82). Once Eq. (E.2) integrated over $E^{3}{\rm d}E$ so as to get a collision term for the brightness only, we can check that we recover Eqs. (6.25-6.26) of Pitrou (2009a).

E.3 Circular polarization

The circular polarization part of the Thomson contribution in a general frame is

[TABLE]

The non-vanishing multipoles in the right hand side are

[TABLE]

and one should use a relation of the form (370) to obtain the decomposition of the first term of Eq. (374) in STF tensors.