Early universe and dark matter phenomenon

S R Jayaganeshan

Indian Institute of Science Education and Research, Maruthamala Post Office, Vithura, Thiruvananthapuram, Kerala 695551

Dr. Biplob Bhattacherjee

Indian Institute of Science, C V Raman Avenue, Bengaluru, Karnataka 560012

Abstract

According to the Cosmological principle, the Universe we observe today is almost homogeneous and isotropic. It is composed of matter-Luminous and Dark matter, Radiation- Relativistic particles like photons , etc. and Vaccum Energy. But everything started from the Big Bang where the energy densities were infinite. Later when the universe started expanding with scale factor a (which is dependent on the density distribution of its constituents), it started cooling down. Many processes like pair production, Big bang Nucleosynthesis, Recombination, and other remarkable events took place at different time periods based on the temperature of the universe. Hence, the temperature history of the universe gives away a lot of facts about the early universe. So, the first part of the project is to study about the fundamentals of the early universe like the evolution of the universe, BBN, Observable Geometries of the space which is dependent on the density parameter Ω_today which is the ratio of the present density of the Universe to the Critical Density of the universe. The Matter we see only constitutes 15% of the total mass and the rest is Dark matter. WIMP (Weakly Interacting Massive Particles) which is thermally produced according to the BigBang Cosmology and comes under the nonrelativistic or Cold Dark matter. When the temperature of the universe was sufficient , the particles were in thermal equilibrium with each other. Hence, dark matter particle and its corresponding antimatter particle would annihilate and give particle-antiparticle pair. As the temperature dropped down and also due to the expansion, the density of the dark matter particles got fixed as the energy for such annihilations was not enough and the density of dark matter remained constant. The second part of the project is to study the basics of Dark Matter, evidence of dark matter and the derivation of the basic equation using Boltzmann's equation and show how the abundance of WIMPs freezes out.

Keywords: Scalefactor, Bigbang Nucleosynthesis, Density parameter, Boltzmann Equation, WIMP, Relic Density.

Abbreviations

INTRODUCTION

Background/Rationale ​

The first part of this project focusses on the early universe, various eras and how the scale factor of universe varies with time by using the Friedmann's Equations which is derived using the basic Newtonian mechanics. Other concepts related to the Early Universe which is needed for the understanding of the second part of the project like Hubble's law, Geometry of the Universe, Density Parameter Ω and the Temperature history of the Universe.

One of the open questions in physics is to identify the nature of DM. All evidences in favour of the particle DM comes from observations of its gravitational effects on baryonic matter. So, the second part of the report provides the fundamentals of DM like its properties, evidences supporting it's evidences and most importantly we focus on the Relic Density Calculation using the Boltzmann Equation. Also here, we have narrowed down to a Thermal DM Candidate called WIMPs. Freezeout of these WIMPs took place in a Radiation dominated Universe, so several parameters related to this era are briefly discussed for the Relic Density calculation.

Objectives of the Research

✽ To study the fundamentals of Early Universe.

✽ Calculation of the Relic Density of Dark matter.

LITERATURE REVIEW

Fundamentals of Early Universe

Introduction

One of models which is used for study of Universe is Friedmann-Robertson-Walker Model(FRW Model). This model is so successful that it is otherwise called as the "Standard Cosmology". The main consideration in this FRW model is that universe is isotropic and it is homogeneous on the large scale.

For an observer in the universe, universe looks same in all the direction and this is known as being isotropic. Whereas if the universe looks the same from all the points i.e., if all the points are same, then it is called homogeneous. Universe needs to be homogeneous,thus Universe is isotropic to every observer from the evidences such as what we see from the cosmic microwave background, which is uniformly distributed over the universe.This means that there is no center for the Universe. It also means that the universe is expanding homogeneously.

Co moving Coordinates and Hubble Law

Universe is homogeneous and it is expanding in all the directions equally as the time progresses. So in order to account for the expansion of the universe we define a new coordinate system with coordinates called co moving coordinates. These coordinates unlike Euclidean coordinates keeps changing with time such that the coordinate grid which is taken for study remains static. To put it simply during the expansion of universe even the grid moves along with the universe so that the galaxies remain at their same positions with respect to the frame we study in.

Mathematically, the coordinate system of consideration expands with a scale parameter a, which is dependent on time. So in co moving coordinates the length of two moving galaxies Δx remain the same with time. But the actual distance between them is D=a(t)Δx.

So the Actual Relative velocity would be derivative of D, V=a(t)Δx. Taking the ratio of V and D we get

$\frac VD=\frac{\displaystyle\dot a}a$<math xmlns="http://www.w3.org/1998/Math/MathML"/>\frac VD=\frac{\displaystyle\dot a}a

where, $\dot{a\;}=\;\frac{\operatorname da}{\operatorname dt}$<math xmlns="http://www.w3.org/1998/Math/MathML"/>\dot{a\;}=\;\frac{\operatorname da}{\operatorname dt} and $\ddot a\;=\;\frac{\operatorname d\dot a}{\operatorname dt}$<math xmlns="http://www.w3.org/1998/Math/MathML"/>\ddot a\;=\;\frac{\operatorname d\dot a}{\operatorname dt}

This ratio is the same for any two bodies taken into account, hence this ratio is called as Hubble Parameter and it is denoted by H. From this we know that V/D = H and V = H × D that is relative Velocity of recession between two bodies will be proportional to the distance between them and this law is known as the "Hubble's Law". Hubble's Law works the best when the two bodies are in a large separation. Nevertheless,gravitational force between the bodies gets them closer and the Hubble's Law seems to be wrong but we take the averaged result into account then Hubble's law will be correct i.e, the universe is expanding.
And the density of the observable matter ρ is proportional to a^-3.

Newtonian Universe

We are now interested to find the dynamics of some galaxy or an object which is at a distance d from the origin i.e., from the point we take even though we know that universe is same everywhere. We later show that the result is same for all the points in the space. Universe is expanding hence we need to go back to our co moving coordinates where everything remains static. Let the distance between them in our coordinate
grid is d, hence the actual distance is D=a(t)d ,velocity and acceleration of both could be found out by simply differentiating and double differentiating D,

$\begin{array}{l}v=\dot ad\;\;\\A=\ddot ad\end{array}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>v</mi><mo>=</mo><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mi>d</mi><mo>&#xA0;</mo><mo>&#xA0;</mo><mspace linebreak="newline"/><mi>A</mi><mo>=</mo><mover><mi>a</mi><mo>&#xA8;</mo></mover><mi>d</mi></math>\begin{array}{l}v=\dot ad\;\;\\A=\ddot ad\end{array}

Hence to find the gravitational force on the object which is at the distance of d from the center, we need to imagine a sphere centred at origin and of radius d thus passing through the object and also due to the isotropic nature of the origin, we imagine the entire mass in the sphere to be concentrated at the centre of the sphere. So Gravitational force between the object and the Mass at the origin is given by the expression

$mA=\frac{-GMm}{D^2}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>m</mi><mi>A</mi><mo>=</mo><mfrac><mrow><mo>-</mo><mi>G</mi><mi>M</mi><mi>m</mi></mrow><msup><mi>D</mi><mn>2</mn></msup></mfrac></math>mA=\frac{-GMm}{D^2}

Using the above results, we get,

$\frac1D\frac{d^2D}{dt^2}=\frac{-GM}{D^3}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mn>1</mn><mi>D</mi></mfrac><mfrac><mrow><msup><mi>d</mi><mn>2</mn></msup><mi>D</mi></mrow><mrow><mi>d</mi><msup><mi>t</mi><mn>2</mn></msup></mrow></mfrac><mo>=</mo><mfrac><mrow><mo>-</mo><mi>G</mi><mi>M</mi></mrow><msup><mi>D</mi><mn>3</mn></msup></mfrac></math>\frac1D\frac{d^2D}{dt^2}=\frac{-GM}{D^3}

Modifying the equations to have a density term in it

$\frac1D\frac{d^2D}{dt^2}=\frac{-4\mathrm\pi G\rho}3$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mn>1</mn><mi>D</mi></mfrac><mfrac><mrow><msup><mi>d</mi><mn>2</mn></msup><mi>D</mi></mrow><mrow><mi>d</mi><msup><mi>t</mi><mn>2</mn></msup></mrow></mfrac><mo>=</mo><mfrac><mrow><mo>-</mo><mn>4</mn><mi mathvariant="normal">&#x3C0;</mi><mi>G</mi><mi>&#x3C1;</mi></mrow><mn>3</mn></mfrac></math>\frac1D\frac{d^2D}{dt^2}=\frac{-4\mathrm\pi G\rho}3

Since it is now a function of ρ, it is the same for any two places in the universe. So our previous claim that all the points are equivalent is proven here. As the final equation is negative on the RHS it implies that the universe is decelerating yet it is still expanding which is in contradiction to the actual case where the universe is actually expanding. And a very important implication that one needs to understand is that if suppose we assume that the universe is static that means the scale parameter doesnot vary with time so LHS in the equation $\frac1D\frac{d^2D}{dt^2}=\frac{-4\mathrm\pi G\rho}3$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mn>1</mn><mi>D</mi></mfrac><mfrac><mrow><msup><mi>d</mi><mn>2</mn></msup><mi>D</mi></mrow><mrow><mi>d</mi><msup><mi>t</mi><mn>2</mn></msup></mrow></mfrac><mo>=</mo><mfrac><mrow><mo>-</mo><mn>4</mn><mi mathvariant="normal">&#x3C0;</mi><mi>G</mi><mi>&#x3C1;</mi></mrow><mn>3</mn></mfrac></math>\frac1D\frac{d^2D}{dt^2}=\frac{-4\mathrm\pi G\rho}3 would turn zero which means the density on the RHS would become zero which is true only when the universe is empty and it is clearly a contradiction hence universe is expanding.

Friedmann's Equations(For Matter Dominated Universe)

Friedmann's approach is similar to that of Newton's except for a fact that it uses Energy conservation for deriving the equation. We take the similar consideration of the total mass of the sphere to be concentrated at the origin and the object to be at a distance of d from it. Then to study the expansion we take the escape velocity into play i.e.,when the velocity with which universe expands exceeds or is equal to the escape velocity then the universe expands. Else it would collapse after some time.

Nevertheless, we are now going to apply the Energy conservation at the critical case that is when the velocity is equal to it's escape velocity. So the Energy that a body will have would be zero, From the law of conservation of energy, one could tell that it is going to be zero at any place and time, So,

$\frac{-GMm}D+\frac12mv^2=0$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mrow><mo>-</mo><mi>G</mi><mi>M</mi><mi>m</mi></mrow><mi>D</mi></mfrac><mo>+</mo><mfrac><mn>1</mn><mn>2</mn></mfrac><mi>m</mi><msup><mi>v</mi><mn>2</mn></msup><mo>=</mo><mn>0</mn></math>\frac{-GMm}D+\frac12mv^2=0

Solving the equation by making use of D = a(t)d and $v=\dot ad$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>v</mi><mo>=</mo><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mi>d</mi></math>v=\dot ad

But, $\frac{\dot a}a=H$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mi>a</mi></mfrac><mo>=</mo><mi>H</mi></math>\frac{\dot a}a=H

So,
$H^{2\;}=\;\frac{8\mathrm{πGρ}}3$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>H</mi><mrow><mn>2</mn><mo>&#xA0;</mo></mrow></msup><mo>=</mo><mo>&#xA0;</mo><mfrac><mrow><mn>8</mn><mi>&#x3C0;G&#x3C1;</mi></mrow><mn>3</mn></mfrac></math>H^{2\;}=\;\frac{8\mathrm{πGρ}}3

we know that ρ =ka^-3 ,where k is a constant. Hence the last equation becomes,

$H^{2\;}=\;\frac{8\mathrm{πG}k}{3a^3}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>H</mi><mrow><mn>2</mn><mo>&#xA0;</mo></mrow></msup><mo>=</mo><mo>&#xA0;</mo><mfrac><mrow><mn>8</mn><mi>&#x3C0;G</mi><mi>k</mi></mrow><mrow><mn>3</mn><msup><mi>a</mi><mn>3</mn></msup></mrow></mfrac></math>H^{2\;}=\;\frac{8\mathrm{πG}k}{3a^3}

Upon solving this equation, we get that scale factor a during this period when the matter dominated the energy density of the universe as,

$a(t)\;=\;ct^\frac23$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>a</mi><mo>(</mo><mi>t</mi><mo>)</mo><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mi>c</mi><msup><mi>t</mi><mfrac><mn>2</mn><mn>3</mn></mfrac></msup></math>a(t)\;=\;ct^\frac23

Hence the double derivative of the above scale factor will give a negative value indicating that the expansion would slow down as the time progresses. We came across Friedmann's Equations in a case when the total energy is Zero. But in real universe that is not really needed.

Universe with non zero energy

It is all similar to that of last case except for a fact that we will be equating the total energy to a constant. Taking the Cosmological Principle into account we get the total energy of a body which is at a distance of R in the comoving coordinates will be,

$\frac{-GMm}D+\frac12mv^2\;=\;E$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mrow><mo>-</mo><mi>G</mi><mi>M</mi><mi>m</mi></mrow><mi>D</mi></mfrac><mo>+</mo><mfrac><mn>1</mn><mn>2</mn></mfrac><mi>m</mi><msup><mi>v</mi><mn>2</mn></msup><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mi>E</mi></math>\frac{-GMm}D+\frac12mv^2\;=\;E

The above expression for total energy(E) is constant for a particular body under consideration. Therefore,

$v^2-\frac{2GM}D=\frac{2E}m$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>v</mi><mn>2</mn></msup><mo>-</mo><mfrac><mrow><mn>2</mn><mi>G</mi><mi>M</mi></mrow><mi>D</mi></mfrac><mo>=</mo><mfrac><mrow><mn>2</mn><mi>E</mi></mrow><mi>m</mi></mfrac></math>v^2-\frac{2GM}D=\frac{2E}m

$\left(\dot a\right)^2\;-\frac{2GM}a=\;k$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mfenced><mover><mi>a</mi><mo>&#x2D9;</mo></mover></mfenced><mn>2</mn></msup><mo>&#xA0;</mo><mo>-</mo><mfrac><mrow><mn>2</mn><mi>G</mi><mi>M</mi></mrow><mi>a</mi></mfrac><mo>=</mo><mo>&#xA0;</mo><mi>k</mi></math>\left(\dot a\right)^2\;-\frac{2GM}a=\;k

,where k is the constant.

Dividing both the sides with a² we get,

$H^2=\frac{8\mathrm{πGρ}}3+\frac C{a^2}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>H</mi><mn>2</mn></msup><mo>=</mo><mfrac><mrow><mn>8</mn><mi>&#x3C0;G&#x3C1;</mi></mrow><mn>3</mn></mfrac><mo>+</mo><mfrac><mi>C</mi><msup><mi>a</mi><mn>2</mn></msup></mfrac></math>H^2=\frac{8\mathrm{πGρ}}3+\frac C{a^2}

where Constant C and the Energy have the same signs.

Variation in the scale factor with the signs of C

i) If C > 0, then Velocity of Recession V > V_escape

Therefore the universe expands continuously and the scale factor is an increasing function of time and it expands more rapidly than the case where the total energy is zero.

Hence,

$H^2=\frac{8\mathrm{πGρ}}3+\frac C{a^2}\;=\;\frac{8\mathrm{πG}k}{3a^3}+\frac C{a^2}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>H</mi><mn>2</mn></msup><mo>=</mo><mfrac><mrow><mn>8</mn><mi>&#x3C0;G&#x3C1;</mi></mrow><mn>3</mn></mfrac><mo>+</mo><mfrac><mi>C</mi><msup><mi>a</mi><mn>2</mn></msup></mfrac><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><mrow><mn>8</mn><mi>&#x3C0;G</mi><mi>k</mi></mrow><mrow><mn>3</mn><msup><mi>a</mi><mn>3</mn></msup></mrow></mfrac><mo>+</mo><mfrac><mi>C</mi><msup><mi>a</mi><mn>2</mn></msup></mfrac></math>H^2=\frac{8\mathrm{πGρ}}3+\frac C{a^2}\;=\;\frac{8\mathrm{πG}k}{3a^3}+\frac C{a^2}

In the early stages of matter domiated era, as the scalefactor was small in magnitude, clearly term having 1/a³ dominates and hence the equation again reduces to Equation 10. Therefore it could be inferred that the early universe expansion followed the zero energy case i.e., a(t) =ct^2/3. In the Later stage of the era, exactly opposite happens . Hence on solving the differential eqn we get a simple eqn of scalefactor as a linear function of time. This means that after a long time the objects get seperated to large distances such that they move with constant velocity.

ii) If C < 0

Even in this case, 1/a³ tends to dominate the first part of expansion, hence the expansion here too follows the same relation with time i.e., a(t) =ct^2/3. After a certain time, the 1/a² term seems to dominate. But since the equation obtained looks like H² = -C/a² , the expansion seizes and it collapses.

The above 2 Subsections were discussed considering a fact that the universe is matter dominated.

Let us see what would be the case when the universe is radiation dominated.

Radiation Domiated Universe

When the Universe is radiation dominated, we mean that the universe is filled with radiation i.e., photons, and other relativistic particles which are massless. So in order to know the trend of expansion during the radiation dominated epoch we need to know how radiation is affected due to the expansion.

We know that the FRW metric is conformal with the Minkowski metric i.e.,the line element of FRW metric could be written as the "product of Conformal time and the line element of minkowski space", where conformal time is dη = dt/ a(t) times the minkowski line element. So the metric of FRW model could be
written as

$ds^2=R^2(\eta)\;\left(d\eta^2-\left[\frac{dr^2}{(1-Kr^2)}+r^2(d\theta^2+(\sin\theta)^2d\phi^2\right]\right)$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>d</mi><msup><mi>s</mi><mn>2</mn></msup><mo>=</mo><msup><mi>R</mi><mn>2</mn></msup><mo>(</mo><mi>&#x3B7;</mi><mo>)</mo><mo>&#xA0;</mo><mfenced><mrow><mi>d</mi><msup><mi>&#x3B7;</mi><mn>2</mn></msup><mo>-</mo><mfenced open="[" close="]"><mrow><mfrac><mrow><mi>d</mi><msup><mi>r</mi><mn>2</mn></msup></mrow><mrow><mo>(</mo><mn>1</mn><mo>-</mo><mi>K</mi><msup><mi>r</mi><mn>2</mn></msup><mo>)</mo></mrow></mfrac><mo>+</mo><msup><mi>r</mi><mn>2</mn></msup><mo>(</mo><mi>d</mi><msup><mi>&#x3B8;</mi><mn>2</mn></msup><mo>+</mo><mo>(</mo><mi>sin</mi><mi>&#x3B8;</mi><msup><mo>)</mo><mn>2</mn></msup><mi>d</mi><msup><mi>&#x3D5;</mi><mn>2</mn></msup></mrow></mfenced></mrow></mfenced></math>ds^2=R^2(\eta)\;\left(d\eta^2-\left[\frac{dr^2}{(1-Kr^2)}+r^2(d\theta^2+(\sin\theta)^2d\phi^2\right]\right)

We are interested to nd the time dependence of Wavelength of radiation. So lets take our co moving coordinates where we are observing a source to be at a co moving radial distance of r = r1 from the origin which is at r = 0. Since we are on the same radial line, dθ and dφ are zero.

Let the source at r₁ emit a wave at time t₁ and the detector receives it only at time t₀. When we look at the geodesic of light we know that ds² = 0, Hence the FRW metric modifies as,

$\int_{t1}^{t_0}\frac{dt}{a(t)}=\int_{r1}^0\frac{dr^2}{(1-Kr^2)}=f(r_1)$<math xmlns="http://www.w3.org/1998/Math/MathML"><msubsup><mo>&#x222B;</mo><mrow><mi>t</mi><mn>1</mn></mrow><msub><mi>t</mi><mn>0</mn></msub></msubsup><mfrac><mrow><mi>d</mi><mi>t</mi></mrow><mrow><mi>a</mi><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mfrac><mo>=</mo><msubsup><mo>&#x222B;</mo><mrow><mi>r</mi><mn>1</mn></mrow><mn>0</mn></msubsup><mfrac><mrow><mi>d</mi><msup><mi>r</mi><mn>2</mn></msup></mrow><mrow><mo>(</mo><mn>1</mn><mo>-</mo><mi>K</mi><msup><mi>r</mi><mn>2</mn></msup><mo>)</mo></mrow></mfrac><mo>=</mo><mi>f</mi><mo>(</mo><msub><mi>r</mi><mn>1</mn></msub><mo>)</mo></math>\int_{t1}^{t_0}\frac{dt}{a(t)}=\int_{r1}^0\frac{dr^2}{(1-Kr^2)}=f(r_1)

Lets suppose a wave crest is produced at time t₁ + δt₁ and the detector detects the signal at time t₀ + δt₀, so

$\int_{t_1+\delta t_1}^{t_0+\delta t_0}\frac{dt}{a(t)}=\int_{r1}^0\frac{dr^2}{(1-Kr^2)}=f(r_1)$<math xmlns="http://www.w3.org/1998/Math/MathML"><msubsup><mo>&#x222B;</mo><mrow><msub><mi>t</mi><mn>1</mn></msub><mo>+</mo><mi>&#x3B4;</mi><msub><mi>t</mi><mn>1</mn></msub></mrow><mrow><msub><mi>t</mi><mn>0</mn></msub><mo>+</mo><mi>&#x3B4;</mi><msub><mi>t</mi><mn>0</mn></msub></mrow></msubsup><mfrac><mrow><mi>d</mi><mi>t</mi></mrow><mrow><mi>a</mi><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mfrac><mo>=</mo><msubsup><mo>&#x222B;</mo><mrow><mi>r</mi><mn>1</mn></mrow><mn>0</mn></msubsup><mfrac><mrow><mi>d</mi><msup><mi>r</mi><mn>2</mn></msup></mrow><mrow><mo>(</mo><mn>1</mn><mo>-</mo><mi>K</mi><msup><mi>r</mi><mn>2</mn></msup><mo>)</mo></mrow></mfrac><mo>=</mo><mi>f</mi><mo>(</mo><msub><mi>r</mi><mn>1</mn></msub><mo>)</mo></math>\int_{t_1+\delta t_1}^{t_0+\delta t_0}\frac{dt}{a(t)}=\int_{r1}^0\frac{dr^2}{(1-Kr^2)}=f(r_1)

Equating Equations we get,

$\int_{t1}^{t_0}dt/a(t)=\int_{t_1+\delta t_1}^{t_0+\delta t_0}dt/a(t)$<math xmlns="http://www.w3.org/1998/Math/MathML"><msubsup><mo>&#x222B;</mo><mrow><mi>t</mi><mn>1</mn></mrow><msub><mi>t</mi><mn>0</mn></msub></msubsup><mi>d</mi><mi>t</mi><mo>/</mo><mi>a</mi><mo>(</mo><mi>t</mi><mo>)</mo><mo>=</mo><msubsup><mo>&#x222B;</mo><mrow><msub><mi>t</mi><mn>1</mn></msub><mo>+</mo><mi>&#x3B4;</mi><msub><mi>t</mi><mn>1</mn></msub></mrow><mrow><msub><mi>t</mi><mn>0</mn></msub><mo>+</mo><mi>&#x3B4;</mi><msub><mi>t</mi><mn>0</mn></msub></mrow></msubsup><mi>d</mi><mi>t</mi><mo>/</mo><mi>a</mi><mo>(</mo><mi>t</mi><mo>)</mo></math>\int_{t1}^{t_0}dt/a(t)=\int_{t_1+\delta t_1}^{t_0+\delta t_0}dt/a(t)

On arranging the limits we will be able to get a result as

$\int_{t1}^{t_1+\delta t_1}dt/a(t)=\int_{t_0}^{t_0+\delta t_0}dt/a(t)$<math xmlns="http://www.w3.org/1998/Math/MathML"><msubsup><mo>&#x222B;</mo><mrow><mi>t</mi><mn>1</mn></mrow><mrow><msub><mi>t</mi><mn>1</mn></msub><mo>+</mo><mi>&#x3B4;</mi><msub><mi>t</mi><mn>1</mn></msub></mrow></msubsup><mi>d</mi><mi>t</mi><mo>/</mo><mi>a</mi><mo>(</mo><mi>t</mi><mo>)</mo><mo>=</mo><msubsup><mo>&#x222B;</mo><msub><mi>t</mi><mn>0</mn></msub><mrow><msub><mi>t</mi><mn>0</mn></msub><mo>+</mo><mi>&#x3B4;</mi><msub><mi>t</mi><mn>0</mn></msub></mrow></msubsup><mi>d</mi><mi>t</mi><mo>/</mo><mi>a</mi><mo>(</mo><mi>t</mi><mo>)</mo></math>\int_{t1}^{t_1+\delta t_1}dt/a(t)=\int_{t_0}^{t_0+\delta t_0}dt/a(t)

If we take suffciently small values for the time difference then we can take that as a constant. Hence,

​$\delta t_0/a(t_0)=\delta t_1/a(t_1)$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>&#x3B4;</mi><msub><mi>t</mi><mn>0</mn></msub><mo>/</mo><mi>a</mi><mo>(</mo><msub><mi>t</mi><mn>0</mn></msub><mo>)</mo><mo>=</mo><mi>&#x3B4;</mi><msub><mi>t</mi><mn>1</mn></msub><mo>/</mo><mi>a</mi><mo>(</mo><msub><mi>t</mi><mn>1</mn></msub><mo>)</mo></math>\delta t_0/a(t_0)=\delta t_1/a(t_1)​

and we know that wavelength of a radiation is inversely proportional to the time, hence

$\lambda_0/a(t_0)=\lambda_1/a(t_1)$<math xmlns="http://www.w3.org/1998/Math/MathML"><msub><mi>&#x3BB;</mi><mn>0</mn></msub><mo>/</mo><mi>a</mi><mo>(</mo><msub><mi>t</mi><mn>0</mn></msub><mo>)</mo><mo>=</mo><msub><mi>&#x3BB;</mi><mn>1</mn></msub><mo>/</mo><mi>a</mi><mo>(</mo><msub><mi>t</mi><mn>1</mn></msub><mo>)</mo></math>\lambda_0/a(t_0)=\lambda_1/a(t_1)

​Therefore the wavelength of a radiation is directly proportional to scale factor,a.We know that the energy and the energy density of the radiation varies inversely with wavelength, therefore the energy density of the radiation per comoving volume varies as 1/a⁴.

So, Using the above result in the Friedmann's Equation, we get

$(\dot a/a)^2=8\pi Gk/3a^4+C/a^2\;$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>(</mo><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mo>/</mo><mi>a</mi><msup><mo>)</mo><mn>2</mn></msup><mo>=</mo><mn>8</mn><mi>&#x3C0;</mi><mi>G</mi><mi>k</mi><mo>/</mo><mn>3</mn><msup><mi>a</mi><mn>4</mn></msup><mo>+</mo><mi>C</mi><mo>/</mo><msup><mi>a</mi><mn>2</mn></msup><mo>&#xA0;</mo></math>(\dot a/a)^2=8\pi Gk/3a^4+C/a^2\;

If C = 0, then $\left(\frac{\dot a}a\right)^2=\frac{8\mathrm{πGk}}{3a^3}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mfenced><mfrac><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mi>a</mi></mfrac></mfenced><mn>2</mn></msup><mo>=</mo><mfrac><mrow><mn>8</mn><mi>&#x3C0;Gk</mi></mrow><mrow><mn>3</mn><msup><mi>a</mi><mn>3</mn></msup></mrow></mfrac></math>\left(\frac{\dot a}a\right)^2=\frac{8\mathrm{πGk}}{3a^3}, on solving the differential equation we get the scale factor as a(t) = ct^1/2. If both the matter and radiation are present in a at universe then the equation looks like$\left(\frac{\dot a}a\right)^2=\frac{C_{matter}}{a^3}+\frac{C_{Radiation}}{a^4}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mfenced><mfrac><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mi>a</mi></mfrac></mfenced><mn>2</mn></msup><mo>=</mo><mfrac><msub><mi>C</mi><mrow><mi>m</mi><mi>a</mi><mi>t</mi><mi>t</mi><mi>e</mi><mi>r</mi></mrow></msub><msup><mi>a</mi><mn>3</mn></msup></mfrac><mo>+</mo><mfrac><msub><mi>C</mi><mrow><mi>R</mi><mi>a</mi><mi>d</mi><mi>i</mi><mi>a</mi><mi>t</mi><mi>i</mi><mi>o</mi><mi>n</mi></mrow></msub><msup><mi>a</mi><mn>4</mn></msup></mfrac></math>\left(\frac{\dot a}a\right)^2=\frac{C_{matter}}{a^3}+\frac{C_{Radiation}}{a^4} So in initial stages of expansion, the 1/a⁴ term dominates, hence one could say that the universe was radiation dominated and in the later phases 1/a³ term becomes more dominant and thus became matter dominated later.

Fig 1 Variation of Scale Factor with time(Blue - Matter dominated & Orange - Radiation dominated) [Plotted using Mathematica]

Geometries of space

We know that the Universe we live in is homogeneous and isotropic on a large scale. Hence on deciding the geometry of the space we need to choose only such geometries which are consistent with the above fact. So one could possibly come up with 3 geometries. They are Flat,Spherical
and Hyperbolic Geometry. Before starting we need to know about the metric of each space. Metric is the distance between any two points in the space. So metric for a Flat 2-surface will be $ds^2\;=\;dx^2\;+\;dy^2\;+\;dz^2$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>d</mi><msup><mi>s</mi><mn>2</mn></msup><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mi>d</mi><msup><mi>x</mi><mn>2</mn></msup><mo>&#xA0;</mo><mo>+</mo><mo>&#xA0;</mo><mi>d</mi><msup><mi>y</mi><mn>2</mn></msup><mo>&#xA0;</mo><mo>+</mo><mo>&#xA0;</mo><mi>d</mi><msup><mi>z</mi><mn>2</mn></msup></math>ds^2\;=\;dx^2\;+\;dy^2\;+\;dz^2.

Similarly, metric for a 1-sphere(One dimensional,also represented by dΩ₁) of
unit radius will be simply dθ² . Metric for a 2-sphere(A sphere which is embedded in 3 dimensions) and its
metric could be found out by a simple transformation of coordinates in the at 2D Surface's metric by using x = sinθcosφ; y = cosθcosφ; z = cosθ:, where θ is the Azimuthal angle and φ is the polar angle.

We get, $ds^2=d\theta^2+(sin\theta)^2d\phi^2=d\theta^2+(sin\theta)^2d\Omega_1$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>d</mi><msup><mi>s</mi><mn>2</mn></msup><mo>=</mo><mi>d</mi><msup><mi>&#x3B8;</mi><mn>2</mn></msup><mo>+</mo><mo>(</mo><mi>s</mi><mi>i</mi><mi>n</mi><mi>&#x3B8;</mi><msup><mo>)</mo><mn>2</mn></msup><mi>d</mi><msup><mi>&#x3D5;</mi><mn>2</mn></msup><mo>=</mo><mi>d</mi><msup><mi>&#x3B8;</mi><mn>2</mn></msup><mo>+</mo><mo>(</mo><mi>s</mi><mi>i</mi><mi>n</mi><mi>&#x3B8;</mi><msup><mo>)</mo><mn>2</mn></msup><mi>d</mi><msub><mi>&#x3A9;</mi><mn>1</mn></msub></math>ds^2=d\theta^2+(sin\theta)^2d\phi^2=d\theta^2+(sin\theta)^2d\Omega_1

We can also extrapolate and derive the metric for 3 sphere which needs to be embedded in a 4-Dimensions,i.e., dθ² + sin²θdΩ², where dΩ² is the metric for a 2-Sphere. We can also derive an another metric for a 2 sphere and generalise it for all the possible space geometries,

We know that the metric,

$ds^2=-dt^2+a^2(t)\left[\begin{array}{cc}d\chi^2+\begin{pmatrix}\sin^2\chi\\\chi^2\\Sinh^2\chi\end{pmatrix}&amp;(d\theta^2+\sin^2\theta d\varphi^2\end{array})\right]$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>d</mi><msup><mi>s</mi><mn>2</mn></msup><mo>=</mo><mo>-</mo><mi>d</mi><msup><mi>t</mi><mn>2</mn></msup><mo>+</mo><msup><mi>a</mi><mn>2</mn></msup><mo>(</mo><mi>t</mi><mo>)</mo><mfenced open="[" close="]"><mrow><mtable><mtr><mtd><mi>d</mi><msup><mi>&#x3C7;</mi><mn>2</mn></msup><mo>+</mo><mfenced><mtable><mtr><mtd><msup><mi>sin</mi><mn>2</mn></msup><mi>&#x3C7;</mi></mtd></mtr><mtr><mtd><msup><mi>&#x3C7;</mi><mn>2</mn></msup></mtd></mtr><mtr><mtd><mi>S</mi><mi>i</mi><mi>n</mi><msup><mi>h</mi><mn>2</mn></msup><mi>&#x3C7;</mi></mtd></mtr></mtable></mfenced></mtd><mtd><mo>(</mo><mi>d</mi><msup><mi>&#x3B8;</mi><mn>2</mn></msup><mo>+</mo><msup><mi>sin</mi><mn>2</mn></msup><mi>&#x3B8;</mi><mi>d</mi><msup><mi>&#x3C6;</mi><mn>2</mn></msup></mtd></mtr></mtable><mo>)</mo></mrow></mfenced></math>ds^2=-dt^2+a^2(t)\left[\begin{array}{cc}d\chi^2+\begin{pmatrix}\sin^2\chi\\\chi^2\\Sinh^2\chi\end{pmatrix}&(d\theta^2+\sin^2\theta d\varphi^2\end{array})\right]

Now generalising this equation and writing in terms of the K( constant related to the curvature), Take the comoving radial distance r, and if

So, the generalised metric in FRW model,

$ds^2=-dt^2+a(t)^2\left(\frac{dr^2}{1-Kr^2}+r^2(d\theta^2+sin^2\theta d\phi^2)\right)$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>d</mi><msup><mi>s</mi><mn>2</mn></msup><mo>=</mo><mo>-</mo><mi>d</mi><msup><mi>t</mi><mn>2</mn></msup><mo>+</mo><mi>a</mi><mo>(</mo><mi>t</mi><msup><mo>)</mo><mn>2</mn></msup><mfenced><mrow><mfrac><mrow><mi>d</mi><msup><mi>r</mi><mn>2</mn></msup></mrow><mrow><mn>1</mn><mo>-</mo><mi>K</mi><msup><mi>r</mi><mn>2</mn></msup></mrow></mfrac><mo>+</mo><msup><mi>r</mi><mn>2</mn></msup><mo>(</mo><mi>d</mi><msup><mi>&#x3B8;</mi><mn>2</mn></msup><mo>+</mo><mi>s</mi><mi>i</mi><msup><mi>n</mi><mn>2</mn></msup><mi>&#x3B8;</mi><mi>d</mi><msup><mi>&#x3D5;</mi><mn>2</mn></msup><mo>)</mo></mrow></mfenced></math>ds^2=-dt^2+a(t)^2\left(\frac{dr^2}{1-Kr^2}+r^2(d\theta^2+sin^2\theta d\phi^2)\right)

From the Friedmann's Equations, we also know that$\left(\frac{\dot a}a\right)^2=\frac{8\mathrm{πGρ}}3-\frac K{a^2}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mfenced><mfrac><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mi>a</mi></mfrac></mfenced><mn>2</mn></msup><mo>=</mo><mfrac><mrow><mn>8</mn><mi>&#x3C0;G&#x3C1;</mi></mrow><mn>3</mn></mfrac><mo>-</mo><mfrac><mi>K</mi><msup><mi>a</mi><mn>2</mn></msup></mfrac></math>\left(\frac{\dot a}a\right)^2=\frac{8\mathrm{πGρ}}3-\frac K{a^2} , by modifying it we can write this as $\frac K{H^2R^2}=\frac\rho{\frac{3H^2}{8\pi G}}=\Omega-1$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mi>K</mi><mrow><msup><mi>H</mi><mn>2</mn></msup><msup><mi>R</mi><mn>2</mn></msup></mrow></mfrac><mo>=</mo><mfrac><mi>&#x3C1;</mi><mfrac><mrow><mn>3</mn><msup><mi>H</mi><mn>2</mn></msup></mrow><mrow><mn>8</mn><mi>&#x3C0;</mi><mi>G</mi></mrow></mfrac></mfrac><mo>=</mo><mi>&#x3A9;</mi><mo>-</mo><mn>1</mn></math>\frac K{H^2R^2}=\frac\rho{\frac{3H^2}{8\pi G}}=\Omega-1, where $\frac{3H^2}{8\pi G}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mrow><mn>3</mn><msup><mi>H</mi><mn>2</mn></msup></mrow><mrow><mn>8</mn><mi>&#x3C0;</mi><mi>G</mi></mrow></mfrac></math>\frac{3H^2}{8\pi G} is the Critical Density of the Universe.

If K is equal to +1, then the space is positively curved and it is closed. So from the Equation $\frac K{H^2R^2}=\frac\rho{\frac{3H^2}{8\pi G}}=\Omega-1$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mi>K</mi><mrow><msup><mi>H</mi><mn>2</mn></msup><msup><mi>R</mi><mn>2</mn></msup></mrow></mfrac><mo>=</mo><mfrac><mi>&#x3C1;</mi><mfrac><mrow><mn>3</mn><msup><mi>H</mi><mn>2</mn></msup></mrow><mrow><mn>8</mn><mi>&#x3C0;</mi><mi>G</mi></mrow></mfrac></mfrac><mo>=</mo><mi>&#x3A9;</mi><mo>-</mo><mn>1</mn></math>\frac K{H^2R^2}=\frac\rho{\frac{3H^2}{8\pi G}}=\Omega-1,so if K is positive then Ω > 1. If K is equal to 0, then the space is infinitely at and it is open. So from the Equation $\frac K{H^2R^2}=\frac\rho{\frac{3H^2}{8\pi G}}=\Omega-1$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mi>K</mi><mrow><msup><mi>H</mi><mn>2</mn></msup><msup><mi>R</mi><mn>2</mn></msup></mrow></mfrac><mo>=</mo><mfrac><mi>&#x3C1;</mi><mfrac><mrow><mn>3</mn><msup><mi>H</mi><mn>2</mn></msup></mrow><mrow><mn>8</mn><mi>&#x3C0;</mi><mi>G</mi></mrow></mfrac></mfrac><mo>=</mo><mi>&#x3A9;</mi><mo>-</mo><mn>1</mn></math>\frac K{H^2R^2}=\frac\rho{\frac{3H^2}{8\pi G}}=\Omega-1, so if K is zero then Ω=1. If K is equal to -1, then the space is negatively curved and it is open. So from the Equation $\frac K{H^2R^2}=\frac\rho{\frac{3H^2}{8\pi G}}=\Omega-1$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mi>K</mi><mrow><msup><mi>H</mi><mn>2</mn></msup><msup><mi>R</mi><mn>2</mn></msup></mrow></mfrac><mo>=</mo><mfrac><mi>&#x3C1;</mi><mfrac><mrow><mn>3</mn><msup><mi>H</mi><mn>2</mn></msup></mrow><mrow><mn>8</mn><mi>&#x3C0;</mi><mi>G</mi></mrow></mfrac></mfrac><mo>=</mo><mi>&#x3A9;</mi><mo>-</mo><mn>1</mn></math>\frac K{H^2R^2}=\frac\rho{\frac{3H^2}{8\pi G}}=\Omega-1 if K is negative then Ω <1.

As discussed above we could actually represent the metric of a space in many forms, let's express the metric of a space in hyperspherical coordinates where the metric is $ds^2\;=\;d\theta^2\;+\;(S_K(\theta))^2(d\Omega)^2$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>d</mi><msup><mi>s</mi><mn>2</mn></msup><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mi>d</mi><msup><mi>&#x3B8;</mi><mn>2</mn></msup><mo>&#xA0;</mo><mo>+</mo><mo>&#xA0;</mo><mo>(</mo><msub><mi>S</mi><mi>K</mi></msub><mo>(</mo><mi>&#x3B8;</mi><mo>)</mo><msup><mo>)</mo><mn>2</mn></msup><mo>(</mo><mi>d</mi><mi>&#x3A9;</mi><msup><mo>)</mo><mn>2</mn></msup></math>ds^2\;=\;d\theta^2\;+\;(S_K(\theta))^2(d\Omega)^2.

where S_K(θ) = sinθ, when K = +1

S_K(θ) = θ, when K = 0

S_K(θ) = sinhθ, when K = -1

So in different geometries the angular diameter differs i.e, as suppose if we want to measure the angular diameter of a distant galaxy it would be at the same radial distance so dr = 0, then dθ = D/S_K, and D is the diameter of the galaxy. If we take Flat universe,the angular diameter dθ varies as D/r , which is expected
as the radial distance increases the angle subtended decreases. If we take 3-Sphere,the angular diameter dθ varies as D/sinr , which means that the angle subtended decreases and then increases due to the behaviour of sine function. If we take Hyperbolic universe,the angular diameter dθ varies as D/sinhr , which means the angle subtended by the observer increases exponentially due to the behaviour of sinh r . But all the above geometries are spatially homogeneous and isotropic on a large scale.

Fig 2 Diagrams of the possible homogeneous and isotropic Geometries of Space in FRW Model.(Closed, Flat and Open)[ Hartle, et al ,2014 ]

Fig 3 Evoution of a(t) in closed, flat and open FRW models[ E. Kolb, et al, 1994 ]

First Law of thermodynamics and Density Parameters

First law of thermodynamics states that dU = dQ + dW. This is nothing but the energy conservation. Let's assume a cube of volume v in a comoving grid. Hence it's physical volume would be V_physical = a³(t)v. And also let the energy density at that time(or temperature) is ρ. Therefore the first law of thermodynamics reduces to, d(ρa³V ) = dQ + dW.

Let the pressure inside this box does work on increasing the volume of this box to V to V+dV, so the internal energy of the material present inside the box decreases.And also as the standard model says that the universe is isotropic. Flow of heat is also isotropic hence So, dQ could be taken as zero. Hence the equation reduces to

d(ρV ) = −PdV. We also substitute P = wρ, where w is called the equation of state. On solving the equation,

​$\rho dV+Vd\rho=-PdV$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>&#x3C1;</mi><mi>d</mi><mi>V</mi><mo>+</mo><mi>V</mi><mi>d</mi><mi>&#x3C1;</mi><mo>=</mo><mo>-</mo><mi>P</mi><mi>d</mi><mi>V</mi></math>\rho dV+Vd\rho=-PdV​

$(P+\rho )dV=-Vd\rho$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>(</mo><mi>P</mi><mo>+</mo><mi>&#x3C1;</mi><mo>)</mo><mi>d</mi><mi>V</mi><mo>=</mo><mo>-</mo><mi>V</mi><mi>d</mi><mi>&#x3C1;</mi></math>$$(P + \rho)dV = -Vd \rho$$

So we get,

${\int }_{}^{}\frac{da}{a}=-{\int }_{}^{}\frac{d\rho }{{\rho }^3(1+w)}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msubsup><mo>&#x222B;</mo><mrow/><mrow/></msubsup><mfrac><mrow><mi>d</mi><mi>a</mi></mrow><mi>a</mi></mfrac><mo>=</mo><mo>-</mo><msubsup><mo>&#x222B;</mo><mrow/><mrow/></msubsup><mfrac><mrow><mi>d</mi><mi>&#x3C1;</mi></mrow><mrow><msup><mi>&#x3C1;</mi><mn>3</mn></msup><mo>(</mo><mn>1</mn><mo>+</mo><mi>w</mi><mo>)</mo></mrow></mfrac></math>\int_{}^{}\frac{da}a=-\int_{}^{}\frac{d\rho}{\rho^3(1+w)} and

$\rho \propto \frac{1}{a^3(1+w)}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>&#x3C1;</mi><mo>&#x221D;</mo><mfrac><mn>1</mn><mrow><msup><mi>a</mi><mn>3</mn></msup><mo>(</mo><mn>1</mn><mo>+</mo><mi>w</mi><mo>)</mo></mrow></mfrac></math>\rho\propto\frac1{a^3(1+w)}​

​​For relativistic particles like neutrinos and the photons $\rho\propto\frac1{a^4}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>&#x3C1;</mi><mo>&#x221D;</mo><mfrac><mn>1</mn><msup><mi>a</mi><mn>4</mn></msup></mfrac></math>\rho\propto\frac1{a^4} and w = 1/3 . For cold and non relativistic particles, as the pressure due to them is negligible on the walls of the imaginary comoving volume we consider Pressure P = 0. Hence, the value of their w =0. Therefore, $\rho\propto\frac1{a^3}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>&#x3C1;</mi><mo>&#x221D;</mo><mfrac><mn>1</mn><msup><mi>a</mi><mn>3</mn></msup></mfrac></math>\rho\propto\frac1{a^3}.

Vaccum Energy

We also assume an another component in the total density of the universe and that is vacuum energy. The peculiar fact about is that it's energy density ρ₀ is constant through out the history of the universe. Hence, the expansion of universe doesnot play any role on it's energy density. So in the equation $\rho\propto\frac1{a^{3(1+w)}}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>&#x3C1;</mi><mo>&#x221D;</mo><mfrac><mn>1</mn><msup><mi>a</mi><mrow><mn>3</mn><mo>(</mo><mn>1</mn><mo>+</mo><mi>w</mi><mo>)</mo></mrow></msup></mfrac></math>\rho\propto\frac1{a^{3(1+w)}}, for the ρ to become constant w = -1. Which means that the pressure applied due to vacuum energy is negative and also from the rst law of thermodynamics dU = PdV, which means that the internal energy increases with the expansion of the universe. The overall Friedmann equation thus reduces to,

$H^2\;=\;(\frac{\dot a}a)^2\;=\;\frac{C_m}{a^3}\;+\;\frac{C_R}{a^4}\;+\;\Lambda\;-\;\frac K{a^2}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>H</mi><mn>2</mn></msup><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mo>(</mo><mfrac><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mi>a</mi></mfrac><msup><mo>)</mo><mn>2</mn></msup><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><msub><mi>C</mi><mi>m</mi></msub><msup><mi>a</mi><mn>3</mn></msup></mfrac><mo>&#xA0;</mo><mo>+</mo><mo>&#xA0;</mo><mfrac><msub><mi>C</mi><mi>R</mi></msub><msup><mi>a</mi><mn>4</mn></msup></mfrac><mo>&#xA0;</mo><mo>+</mo><mo>&#xA0;</mo><mi>&#x39B;</mi><mo>&#xA0;</mo><mo>-</mo><mo>&#xA0;</mo><mfrac><mi>K</mi><msup><mi>a</mi><mn>2</mn></msup></mfrac></math>H^2\;=\;(\frac{\dot a}a)^2\;=\;\frac{C_m}{a^3}\;+\;\frac{C_R}{a^4}\;+\;\Lambda\;-\;\frac K{a^2}

,where Λ is the Cosmological constant. $\Lambda=\frac{8\pi G\rho_0}3$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>&#x39B;</mi><mo>=</mo><mfrac><mrow><mn>8</mn><mi>&#x3C0;</mi><mi>G</mi><msub><mi>&#x3C1;</mi><mn>0</mn></msub></mrow><mn>3</mn></mfrac></math>\Lambda=\frac{8\pi G\rho_0}3 and ρ₀ is the energy density of vacuum energy.

Universe started off by being radiation dominated and then slowly it moved onto a matter dominated universe. Now the major contribution to the energy density comes from vacuum energy. From the observations, density parameter Ω was observed to be close to one and the density parameter for matter, radiation and vacuum energy was observed to be,

We determine the geometry of the universe by counting the number of galaxies in a particular distance and we relate it to the red shift and the luminosity. Then we analyse the cosmological parameters which are pivotal in understanding the evolution of universe i.e., Hubble Parameter H, Density Parameter of matter Ω_m , Density Parameter of Radiation Ω_r and Density Parameter of Vacuum Energy Ω_Λ, where Density parameter is equal to $\frac{Present\;Density}{Critical\;Density}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mrow><mi>P</mi><mi>r</mi><mi>e</mi><mi>s</mi><mi>e</mi><mi>n</mi><mi>t</mi><mo>&#xA0;</mo><mi>D</mi><mi>e</mi><mi>n</mi><mi>s</mi><mi>i</mi><mi>t</mi><mi>y</mi></mrow><mrow><mi>C</mi><mi>r</mi><mi>i</mi><mi>t</mi><mi>i</mi><mi>c</mi><mi>a</mi><mi>l</mi><mo>&#xA0;</mo><mi>D</mi><mi>e</mi><mi>n</mi><mi>s</mi><mi>i</mi><mi>t</mi><mi>y</mi></mrow></mfrac></math>\frac{Present\;Density}{Critical\;Density} .

We know that the critical density of the universe is when the total energy of the universe becomes zero or simply the universe becomes at i.e, K= 0. From the Fried man's Equations, $H^2=(\frac{\dot a}a)^2=\frac{8\rho_{critical}\pi G}3=0$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>H</mi><mn>2</mn></msup><mo>=</mo><mo>(</mo><mfrac><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mi>a</mi></mfrac><msup><mo>)</mo><mn>2</mn></msup><mo>=</mo><mfrac><mrow><mn>8</mn><msub><mi>&#x3C1;</mi><mrow><mi>c</mi><mi>r</mi><mi>i</mi><mi>t</mi><mi>i</mi><mi>c</mi><mi>a</mi><mi>l</mi></mrow></msub><mi>&#x3C0;</mi><mi>G</mi></mrow><mn>3</mn></mfrac><mo>=</mo><mn>0</mn></math>H^2=(\frac{\dot a}a)^2=\frac{8\rho_{critical}\pi G}3=0 and $\rho_{critical}=\frac{3H^2}{8\pi G}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msub><mi>&#x3C1;</mi><mrow><mi>c</mi><mi>r</mi><mi>i</mi><mi>t</mi><mi>i</mi><mi>c</mi><mi>a</mi><mi>l</mi></mrow></msub><mo>=</mo><mfrac><mrow><mn>3</mn><msup><mi>H</mi><mn>2</mn></msup></mrow><mrow><mn>8</mn><mi>&#x3C0;</mi><mi>G</mi></mrow></mfrac></math>\rho_{critical}=\frac{3H^2}{8\pi G}.

Thermal Equilibrium in Early Universe

Early universe was very hot and the radiation(photons), protons,and electrons(even before that quarks were in a thermal equilibrium i.e.,before QCD phase transition at T ≈214MeV were in thermal equilibrium with each other, i.e., these particles were scattering(Charged particles electrons and protons being very good scatterers of photons) and were interacting with each other such that the interaction rate was much greater than the expansion rate of the universe. Due to these interactions, photons of various wavelengths were produced and the radiation was of course a black body radiation. So there was a common temperature for the entire universe during this time(Apart from the previously decoupled particles which had their own temperature independent of the thermal bath).

As the universe expanded, and the internal energy of these constituents were utilised for the expansion, they got cooled down to a certain temperature. At this stage approximately about 3 minutes after the big bang(which is an expansion that took place all over the space isotropically ), many light nuclei and their isotopes were created and the stable lightest nuclei formed was Helium-4. Now the amount of hydrogen and other lighter nuclei are fixed and this particular stage in the universe is referred to as "Big Bang Nucleosynthesis". This happened at around 10 - 0.1 MeV.

After this, the temperature fell even low to such a value that the photons which were no more in thermal equilibrium and thus stopped interacting(scattering) due to the formation of Hydrogen atoms and other atoms, and were not having energy to ionise the hydrogen atom and other atoms. So we know that atoms are neutral. We also know that the neutral bodies are poor scatterers of photons and hence interactions between photons and atoms now became far lesser than the initial stage and the photons now got decoupled. Due to the photon decoupling, the universe became transparent. The temperature where this decoupling took place was around 4000K-5000K. Now the maximum intensity of the relic radiation which was decoupled from the matter during the recombination era exists in the microwave region, due to the cosmological red shift and the temperature dependence on scale factor as 1/a .

The relic radiation which is now observed has its maximum contribution of its energy density from the microwave region and thus this is also referred to as the cosmic microwave background.It is also noted that this relic radiation is a good approximation for a radiation coming from a black body at 2.97K.

We know that the wavelength of the light observed at a time is directly proportional to the scale factor at that time and the temperature is inversely proportional.

So,

$\frac{a_{emitted}}{a_{today}}=\frac{\lambda_{emitted}}{\lambda_{today}}=\frac{T_{today}}{T_{emitted}}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><msub><mi>a</mi><mrow><mi>e</mi><mi>m</mi><mi>i</mi><mi>t</mi><mi>t</mi><mi>e</mi><mi>d</mi></mrow></msub><msub><mi>a</mi><mrow><mi>t</mi><mi>o</mi><mi>d</mi><mi>a</mi><mi>y</mi></mrow></msub></mfrac><mo>=</mo><mfrac><msub><mi>&#x3BB;</mi><mrow><mi>e</mi><mi>m</mi><mi>i</mi><mi>t</mi><mi>t</mi><mi>e</mi><mi>d</mi></mrow></msub><msub><mi>&#x3BB;</mi><mrow><mi>t</mi><mi>o</mi><mi>d</mi><mi>a</mi><mi>y</mi></mrow></msub></mfrac><mo>=</mo><mfrac><msub><mi>T</mi><mrow><mi>t</mi><mi>o</mi><mi>d</mi><mi>a</mi><mi>y</mi></mrow></msub><msub><mi>T</mi><mrow><mi>e</mi><mi>m</mi><mi>i</mi><mi>t</mi><mi>t</mi><mi>e</mi><mi>d</mi></mrow></msub></mfrac></math>\frac{a_{emitted}}{a_{today}}=\frac{\lambda_{emitted}}{\lambda_{today}}=\frac{T_{today}}{T_{emitted}}

, the time when the photons were emitted corresponds to that period of time when recombination took place or it could be rephrased as the era when the temperature was such that photons had enough energy to ionise hydrogen atoms. The ionisation energy of hydrogen atom is 13.4 eV. Hence, we consider the Boltzmann Distribution where probability that a photon has enough energy to ionise the hydrogen atoms is equal to $e^{\left(\frac{-E_{ionisation}}{K_BT}\right)}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>e</mi><mfenced><mfrac><mrow><mo>-</mo><msub><mi>E</mi><mrow><mi>i</mi><mi>o</mi><mi>n</mi><mi>i</mi><mi>s</mi><mi>a</mi><mi>t</mi><mi>i</mi><mi>o</mi><mi>n</mi></mrow></msub></mrow><mrow><msub><mi>K</mi><mi>B</mi></msub><mi>T</mi></mrow></mfrac></mfenced></msup></math>e^{\left(\frac{-E_{ionisation}}{K_BT}\right)} . Incidentally, the number of photons to electrons ratio was 10⁸, So number of photons which have energy to ionise is $e^{\left(\frac{-E_{ionisation}}{K_BT}\right)}\times10^8$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>e</mi><mfenced><mfrac><mrow><mo>-</mo><msub><mi>E</mi><mrow><mi>i</mi><mi>o</mi><mi>n</mi><mi>i</mi><mi>s</mi><mi>a</mi><mi>t</mi><mi>i</mi><mi>o</mi><mi>n</mi></mrow></msub></mrow><mrow><msub><mi>K</mi><mi>B</mi></msub><mi>T</mi></mrow></mfrac></mfenced></msup><mo>&#xD7;</mo><msup><mn>10</mn><mn>8</mn></msup></math>e^{\left(\frac{-E_{ionisation}}{K_BT}\right)}\times10^8​

This number is taken to be critical that is assuming that there is one photon for every atom, so we equate the preceding number to one. A good approximation is to take 10⁸ as e²⁰, thus $e^{\left(\frac{-E_{ionisation}}{K_BT}\right)}\times e^{20}\;=\;1$<math xmlns="http://www.w3.org/1998/Math/MathML"><msup><mi>e</mi><mfenced><mfrac><mrow><mo>-</mo><msub><mi>E</mi><mrow><mi>i</mi><mi>o</mi><mi>n</mi><mi>i</mi><mi>s</mi><mi>a</mi><mi>t</mi><mi>i</mi><mi>o</mi><mi>n</mi></mrow></msub></mrow><mrow><msub><mi>K</mi><mi>B</mi></msub><mi>T</mi></mrow></mfrac></mfenced></msup><mo>&#xD7;</mo><msup><mi>e</mi><mn>20</mn></msup><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mn>1</mn></math>e^{\left(\frac{-E_{ionisation}}{K_BT}\right)}\times e^{20}\;=\;1i.e., $K_BT\;=\;\frac{E_{ionisation}}{20}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msub><mi>K</mi><mi>B</mi></msub><mi>T</mi><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><msub><mi>E</mi><mrow><mi>i</mi><mi>o</mi><mi>n</mi><mi>i</mi><mi>s</mi><mi>a</mi><mi>t</mi><mi>i</mi><mi>o</mi><mi>n</mi></mrow></msub><mn>20</mn></mfrac></math>K_BT\;=\;\frac{E_{ionisation}}{20},which corresponds to a statement that the one in twentieth part of temperature will be enough for the ionisation this is due to the excess of photons.​

From the previous Temperature and Scale factor relation, we know that$\frac{T_{decoupling}}{T_{today}}=\frac{\lambda_{today}}{\lambda_{decoupling}}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><msub><mi>T</mi><mrow><mi>d</mi><mi>e</mi><mi>c</mi><mi>o</mi><mi>u</mi><mi>p</mi><mi>l</mi><mi>i</mi><mi>n</mi><mi>g</mi></mrow></msub><msub><mi>T</mi><mrow><mi>t</mi><mi>o</mi><mi>d</mi><mi>a</mi><mi>y</mi></mrow></msub></mfrac><mo>=</mo><mfrac><msub><mi>&#x3BB;</mi><mrow><mi>t</mi><mi>o</mi><mi>d</mi><mi>a</mi><mi>y</mi></mrow></msub><msub><mi>&#x3BB;</mi><mrow><mi>d</mi><mi>e</mi><mi>c</mi><mi>o</mi><mi>u</mi><mi>p</mi><mi>l</mi><mi>i</mi><mi>n</mi><mi>g</mi></mrow></msub></mfrac></math>\frac{T_{decoupling}}{T_{today}}=\frac{\lambda_{today}}{\lambda_{decoupling}} , and we know that the temperature during the decoupling was around 5000K and now it is around 3K, therefore the ratio is approximately equal to 1000 and this also means that the recombination actually took place when the universe was 1000 times smaller than what is today.

Matter Radiation Crossover

In the early stages of universe, it was dominated by the radiation or the relativistic particles like photons and neutrinos. Later a the time progress the universe slowly became matter dominated. We are now interested to nd that part of time where this crossover between radiation and matter actually took place.

Let us see the value of the ratio of energy densities of matter and radiation.

We know that the energy density of matter varies as 1/a³ and the energy density of radiation varies as 1/a⁴ . So this ratio of matter over radiation gives the dependence over the scale factor a.Whose value today is approximately 10⁶ and the ratio would have taken a value equals to one at the time of crossover, hence the
scale factor, a_crossover = 10^-6 × a_today.

To generalise everything which has been discussed above, when today's scale factor is taken as a_today, then the scale factor at the time of crossover will be 10^-6 ×a_today and the temperature during the crossover would be a million times the today's temperature( Today's temperature we mean by the temperature
of the relic radiation left).

When the temperature was even more higher, at T ≈ 10¹⁴K(10 raised to 14), now the photons were very energetic and the photons which belonged tothe gamma region which had an energy E ≈ 0.5 MeV underwent annihilation and resulted in an electron-positron pair. During this stage, the particles, antiparticles and the photons were in thermal equilibrium with each other.

Dark Matter Phenomenon

Introduction

The ΛCDM model of Cosmology is a paradigm in which the Universe consists of three major components which are Vaccum or Dark Energy(also referred to as Cosmological Constant, Λ),Dark matter and the ordinary visible matter.

The Measurements by WMAP and Planck(2018 results) suggest that dark matter which is a major contributor to the matter in this universe(i.e., 85 % of the total matter). The Evidence of Dark matter mostly comes from their gravitational interactions with the ordinary baryonic matter.

Evidences

[​Cerdeno, 2018​]There are a few evidences that actually give us a belief in this enigma. Some of the astrophysical evidences concerning the presence of DM are:-

1) Rotation Curves of Spiral Galaxies:- ​​

With the help of Newtonian mechanics, we can tell that the velocity of body which is rotating under an influence of Gravitational force of a massive body decreases as the distance from it to the massive body increases. This could be seen easily by the relation,

,but on observing the rotation curves(which is a plot of velocity against the radial distance from the center of the spiral galaxy) the velocity as we move away from the center was almost constant and a flat rotation curve was obtained. This was first approached with MOND and this hypothesis didn’t seem to work for other phenomena.

So this was successfully explained by taking velocity to be constant in the above formula, which gives that

M(r) ∝ r. So the mass of the matter present increases radially. But this did not seem to make sense as the visible matter was more concentrated at the center of the galaxy. So in order to account for this observation, extra mass was considered to be present and was also found that it is a spherical, isothermal halo of DM which accounts for 80 percent of the total galaxy.

Fig 4 Rotation Curve of a spiral galaxy[​Cerdeno, 2018​]​

2) Bullet Clusters:-

The Bullet Cluster is considered as a typical example of the effect of dark matter in dynamical systems. Here, two galaxy clusters undergo collision. From a characteristic shock wave, the visible components of the cluster observed by the Chandra Xray satellite. But it was revealed that most of the mass of the system is displaced from the visible components from the weak-lensing analyses, data from the Hubble Space Telescope . The accepted interpretation is that the dark matter components of the clusters have crossed without interacting significantly .

The Bullet Cluster is considered one of the best arguments against MOND theories (since the gravitational effects occur where there is no visible matter) and also sets an upper bound on the self-interaction strength of dark matter particles.

Fig 5 Galaxy - Cluster collisions[​Teresa Marrodán Undagoitia​]

Properties

They are called dark because they do not emit light. With a variety of other astrophysical measurements,and make some general and quantitative statements about properties of Dark Matter:

1) Neutral:-

It is generally argued that Dark Matter particles must be electrically neutral or should not have any charge. Otherwise they would scatter light and thus not be dark.

2) Non Relativistic:-

Numerical Simulations of structure formation in the Early universe helps us understand that the dark matter has to be cold(or non relativistic). The major contribution of the dark matter comes from the non relativistic part.

3) Non Baryonic:-

Predictions from BBN about the total baryonic content and CMB results, suggest that only 4 to 5 percent of the total energy density of the universe comes from ordinary (baryonic) matter. Therefore, we must conclude that DM is mostly non-baryonic.

4) Stable:-

DM is a long-lived (if not stable) particle. CMB anisotropies show us the footprint of DM. Its presence is important for structure formation and feel its gravitational effects in clusters of galaxies and galaxies today . Also, DM particles can decay, if their lifetime is longer than the age of the universe.

Dark Matter Candidates

1) MACHOs

Actually astronomical objects which are too dark to be seen like Jupiters,brown dwarfs etc. are considered as the Dark matter candidates. But their contribution to the total dark matter present is less as we know that the DM is mostly non-baryonic. Generally, one refers to such candidates which are made of ordinary matter as MACHOs( Massive Compact Halo Objects). We usually infer the amount of MACHOs in our galactic halo by microlensing, where a massive body bends the light of the source just behind the MACHO[​Hitoshi Murayama, 2007​].

2) Neutrinos

Neutrinos being non baryonic also contributes to a portion of $\Omega_\nu h^2=\frac{\sum_i^{}{m_\nu}_{i_{}}}{94eV}$<math xmlns="http://www.w3.org/1998/Math/MathML"><msub><mi>&#x3A9;</mi><mi>&#x3BD;</mi></msub><msup><mi>h</mi><mn>2</mn></msup><mo>=</mo><mfrac><mrow><msubsup><mo>&#x2211;</mo><mi>i</mi><mrow/></msubsup><msub><msub><mi>m</mi><mi>&#x3BD;</mi></msub><msub><mi>i</mi><mrow/></msub></msub></mrow><mrow><mn>94</mn><mi>e</mi><mi>V</mi></mrow></mfrac></math>\Omega_\nu h^2=\frac{\sum_i^{}{m_\nu}_{i_{}}}{94eV}[​Hitoshi Murayama, 2007​], which gives their contribution to be only 0.2 percent of the matter density.​ And more over neutrinos are relativistic(Travelling at the speed of light) even today and thus contributing to the Hot DM. In anyways, neutrinos cannot saturate the amount of non baryonic DM estimated.

3) Thermal DM Candidates:-

Thermal DM candidates refer to those DM which were in thermal equilibrium with the SM thermal bath. Furthermore, the assumption of thermal equilibrium implies some level of interaction between DM and the SM, leading to a variety of interesting and testable signatures. Thermal candidates can be as light as keV and can be massive than 100 TeV. Therefore, the mass ranges under consideration for Thermal DM Candidates are :- Light DM which is usually in keV – 10 GeV and Weakly interacting massive particles (WIMPs) fall under the mass range of 10 GeV – 10 TeV:. WIMPs are also sometimes used to mean thermal DM candidates as a whole. WIMPs are cold DM.[ ​Tongyan Lin​ ]

Thermal candidates here are typically examples of dark sector models, and will be the emphasis of these lectures.

The relic abundance of cold DM obtained from the Planck 2018 analysis of CMB observations is

Ω_ch²= 0.120 ±0.001

,where h =H₀ in units of 100 km/s/Mpc and H0 is the rate of expansion of the universe today. Using h = 0.68, then Ω_c = 0.259 and DM is 26 % of the energy density today.

Freeze out of species from the thermal bath

We know that in the early universe all the particles were in chemical and thermal equilibrium. That means that the interaction between all the particles were taking place at an appreciably higher rate compared to the Hubble Expansion Factor. The temperature of the thermal bath having all the particles(SM particles and DM(Thermal candidates)) was so high that the massive particles were created at this period and of course annihilation of these massive particles with it’s antiparticle resulted in an another SM particle and antiparticle. Apart from the number changing processes(Creation and Annihilation, which are inelastic), elastic scattering processes took place and this enabled that the particles interacting were having the same temperature is maintained between them.

Now, we talk about the time when universe started cooling down and the lighter particles in the thermal bath no longer had energy enough to create the massive particle and this led to the cease in the production of massive particles. But these massive particles continued to annihilate and their numbers gradually dropped down. But they didn’t vanish away as the expansion rate dominated their annihilation rate at the time of freezeout(When we refer to freezeout here, it means that the interactions that took place between two massive particles or self annihilation is no longer effective). This is represented by Γ_A > H.

Let’s assume the massive particle which we were talking about in the previous paragraph is the Cold DM candidate i.e., WIMP(whose mass is of ≈ GeV).

.Then the abundance of the WIMP particles should have frozen out. Such Abundances could be calculated by using the Boltzmann Equation.

Summary

The main focus in the project is to do the Relic Density calculations for the WIMP Particles which are thermal candidates that have frozen out during the Radiation Dominated epoch as a result of their interaction rate becoming very negligible compared to the expansion rate, H .

METHODOLOGY

Calculations using the Boltzmann Equation

The evolution of the number density operator can be computed by applying the covariant form of Liouvilles operator to the corresponding phase space distribution function f.

$\widehat{L_{}}\left[f\right]=C\left[f\right]$<math xmlns="http://www.w3.org/1998/Math/MathML"><mover><msub><mi>L</mi><mrow/></msub><mo>^</mo></mover><mfenced open="[" close="]"><mi>f</mi></mfenced><mo>=</mo><mi>C</mi><mfenced open="[" close="]"><mi>f</mi></mfenced></math>\widehat{L_{}}\left[f\right]=C\left[f\right]

,where $\widehat L$<math xmlns="http://www.w3.org/1998/Math/MathML"><mover><mi>L</mi><mo>^</mo></mover></math>\widehat Lis the Liouville operator and C is the Collision Operator.

​$\widehat{L_{{{}}_{}}}\;=p^\mu\frac{\partial{}}{\partial x^\mu}\;-\;{\Gamma_{\sigma\rho}^\mu}_{}p^\sigma p^\rho\frac{\partial{}}{\partial p^\mu}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mover><msub><mi>L</mi><msub><maction actiontype="argument"><mrow/></maction><mrow/></msub></msub><mo>^</mo></mover><mo>&#xA0;</mo><mo>=</mo><msup><mi>p</mi><mi>&#x3BC;</mi></msup><mfrac><mrow><mo>&#x2202;</mo><maction actiontype="argument"><mrow/></maction></mrow><mrow><mo>&#x2202;</mo><msup><mi>x</mi><mi>&#x3BC;</mi></msup></mrow></mfrac><mo>&#xA0;</mo><mo>-</mo><mo>&#xA0;</mo><msub><msup><msub><mi>&#x393;</mi><mrow><mi>&#x3C3;</mi><mi>&#x3C1;</mi></mrow></msub><mi>&#x3BC;</mi></msup><mrow/></msub><msup><mi>p</mi><mi>&#x3C3;</mi></msup><msup><mi>p</mi><mi>&#x3C1;</mi></msup><mfrac><mrow><mo>&#x2202;</mo><maction actiontype="argument"><mrow/></maction></mrow><mrow><mo>&#x2202;</mo><msup><mi>p</mi><mi>&#x3BC;</mi></msup></mrow></mfrac></math>\widehat{L_{{{}}_{}}}\;=p^\mu\frac{\partial{}}{\partial x^\mu}\;-\;{\Gamma_{\sigma\rho}^\mu}_{}p^\sigma p^\rho\frac{\partial{}}{\partial p^\mu}​

In FRW Universe, as the universe is homogeneous and isotropic, $f\left(x^\mu,\;p^\mu\right)=\;f\left(t,E\right)=\;f(t,\boldsymbol p)$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>f</mi><mfenced><mrow><msup><mi>x</mi><mi>&#x3BC;</mi></msup><mo>,</mo><mo>&#xA0;</mo><msup><mi>p</mi><mi>&#x3BC;</mi></msup></mrow></mfenced><mo>=</mo><mo>&#xA0;</mo><mi>f</mi><mfenced><mrow><mi>t</mi><mo>,</mo><mi>E</mi></mrow></mfenced><mo>=</mo><mo>&#xA0;</mo><mi>f</mi><mo>(</mo><mi>t</mi><mo>,</mo><mi mathvariant="bold-italic">p</mi><mo>)</mo></math>f\left(x^\mu,\;p^\mu\right)=\;f\left(t,E\right)=\;f(t,\boldsymbol p)​

So, the Liouville operator modifies to

​$\widehat{L_{{}_{}}}\;=E\frac\partial{\partial t}\;-\;H\left|{}_{}p\right|^2\frac\partial{\partial E}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mover><msub><mi>L</mi><msub><mrow/><mrow/></msub></msub><mo>^</mo></mover><mo>&#xA0;</mo><mo>=</mo><mi>E</mi><mfrac><mo>&#x2202;</mo><mrow><mo>&#x2202;</mo><mi>t</mi></mrow></mfrac><mo>&#xA0;</mo><mo>-</mo><mo>&#xA0;</mo><mi>H</mi><msup><mfenced open="|" close="|"><mmultiscripts><mi>p</mi><mprescripts/><mrow/><none/></mmultiscripts></mfenced><mn>2</mn></msup><mfrac><mo>&#x2202;</mo><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac></math>\widehat{L_{{}_{}}}\;=E\frac\partial{\partial t}\;-\;H\left|{}_{}p\right|^2\frac\partial{\partial E}​

Multiplying $\frac g{\left(2\mathrm\pi\right)^3}\frac{d^3p}E$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mi>g</mi><msup><mfenced><mrow><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi></mrow></mfenced><mn>3</mn></msup></mfrac><mfrac><mrow><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></mrow><mi>E</mi></mfrac></math>\frac g{\left(2\mathrm\pi\right)^3}\frac{d^3p}Eon both sides of the Equation 29 and integrating it , where g is the intrinsic degrees of freedom of the particle. We define $f$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>f</mi></math>fas the phase space distribtution function ,where $f\left(\boldsymbol p,t\right)d^3\boldsymbol p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>f</mi><mfenced><mrow><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi></mrow></mfenced><msup><mi>d</mi><mn>3</mn></msup><mi mathvariant="bold-italic">p</mi></math>f\left(\boldsymbol p,t\right)d^3\boldsymbol p gives the probability of particles occupying the state with momentum $\mathit{p}\mathbf{ }$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi mathvariant="bold-italic">p</mi><mo mathvariant="bold">&#xA0;</mo></math>\boldsymbol p\boldsymbol\;to $\mathit{p}\mathbf{ }\mathbf{+}\mathbf{ }d\mathit{p}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi mathvariant="bold-italic">p</mi><mo mathvariant="bold">&#xA0;</mo><mo mathvariant="bold">+</mo><mo mathvariant="bold">&#xA0;</mo><mi>d</mi><mi mathvariant="bold-italic">p</mi></math>\boldsymbol p\boldsymbol\;\boldsymbol+\boldsymbol\;d\boldsymbol p.

So,

Taking the LHS, $\frac g{\left(2\mathrm\pi\right)^3}\int\overset{}{\frac{\widehat L\left[f\right]}Ed^3p}\;=\;\frac g{(2\mathrm\pi)^3}\int\frac1E\left(E\frac{\partial f}{\partial t}\;-\;H\left|{}_{}p\right|^2\frac{\partial f}{\partial E}\right)d^3p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mi>g</mi><msup><mfenced><mrow><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi></mrow></mfenced><mn>3</mn></msup></mfrac><mo>&#x222B;</mo><mover><mrow><mfrac><mrow><mover><mi>L</mi><mo>^</mo></mover><mfenced open="[" close="]"><mi>f</mi></mfenced></mrow><mi>E</mi></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></mrow><mrow/></mover><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><mi>g</mi><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mo>&#x222B;</mo><mfrac><mn>1</mn><mi>E</mi></mfrac><mfenced><mrow><mi>E</mi><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi></mrow><mrow><mo>&#x2202;</mo><mi>t</mi></mrow></mfrac><mo>&#xA0;</mo><mo>-</mo><mo>&#xA0;</mo><mi>H</mi><msup><mfenced open="|" close="|"><mmultiscripts><mi>p</mi><mprescripts/><mrow/><none/></mmultiscripts></mfenced><mn>2</mn></msup><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac></mrow></mfenced><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></math>\frac g{\left(2\mathrm\pi\right)^3}\int\overset{}{\frac{\widehat L\left[f\right]}Ed^3p}\;=\;\frac g{(2\mathrm\pi)^3}\int\frac1E\left(E\frac{\partial f}{\partial t}\;-\;H\left|{}_{}p\right|^2\frac{\partial f}{\partial E}\right)d^3p

$\;=\;\frac g{(2\mathrm\pi)^3}\int\frac{\partial f(\boldsymbol p,t)}{\partial t}\;-\;\frac HE\left|{}_{}p\right|^2\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><mi>g</mi><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mo>&#x222B;</mo><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>t</mi></mrow></mfrac><mo>&#xA0;</mo><mo>-</mo><mo>&#xA0;</mo><mfrac><mi>H</mi><mi>E</mi></mfrac><msup><mfenced open="|" close="|"><mmultiscripts><mi>p</mi><mprescripts/><mrow/><none/></mmultiscripts></mfenced><mn>2</mn></msup><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></math>\;=\;\frac g{(2\mathrm\pi)^3}\int\frac{\partial f(\boldsymbol p,t)}{\partial t}\;-\;\frac HE\left|{}_{}p\right|^2\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p

$\;=\;\frac g{(2\mathrm\pi)^3}\int\frac{\partial f(\boldsymbol p,t)}{\partial t}d^3p\;-\;\frac g{(2\mathrm\pi)^3}\int\frac HE\left|{}_{}p\right|^2\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><mi>g</mi><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mo>&#x222B;</mo><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>t</mi></mrow></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi><mo>&#xA0;</mo><mo>-</mo><mo>&#xA0;</mo><mfrac><mi>g</mi><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mo>&#x222B;</mo><mfrac><mi>H</mi><mi>E</mi></mfrac><msup><mfenced open="|" close="|"><mmultiscripts><mi>p</mi><mprescripts/><mrow/><none/></mmultiscripts></mfenced><mn>2</mn></msup><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></math>\;=\;\frac g{(2\mathrm\pi)^3}\int\frac{\partial f(\boldsymbol p,t)}{\partial t}d^3p\;-\;\frac g{(2\mathrm\pi)^3}\int\frac HE\left|{}_{}p\right|^2\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p

We know that Number density, $n==\frac{g}{(2\pi {)}^3}\int f(\mathbf{p},t)d^{3}p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>n</mi><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><mi>g</mi><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mo>&#x222B;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo mathvariant="bold">,</mo><mi>t</mi><mo>)</mo><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi><mo>&#xA0;</mo></math>n\;=\;\;=\;\frac g{(2\mathrm\pi)^3}\int f(\boldsymbol p\boldsymbol,t)d^3p\; and differentiating it with respect to

time we get ,$\frac{dn}{dt}=\dot{n}==\frac{g}{(2\pi {)}^3}\int \frac{\partial f(\mathbf{p},t)}{\partial t}{d}^{3}p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mrow><mo>d</mo><mi>n</mi></mrow><mrow><mo>d</mo><mi>t</mi></mrow></mfrac><mo>=</mo><mo>&#xA0;</mo><mover><mi>n</mi><mo>&#x2D9;</mo></mover><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><mstyle displaystyle="true"><mi>g</mi></mstyle><mstyle displaystyle="true"><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mstyle></mfrac><mo>&#x222B;</mo><mfrac><mstyle displaystyle="true"><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mstyle><mstyle displaystyle="true"><mo>&#x2202;</mo><mi>t</mi></mstyle></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi><mo>&#xA0;</mo></math>\frac{\operatorname dn}{\operatorname dt}=\;\dot n\;=\;=\;\frac{\displaystyle g}{\displaystyle(2\mathrm\pi)^3}\int\frac{\displaystyle\partial f(\boldsymbol p,t)}{\displaystyle\partial t}d^3p\;. Also,

$\;\int\frac{\left|p\right|^2}E\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p\;=\;\int\frac{\left|p\right|^2}E\frac{\partial f(\boldsymbol p,t)}{\partial E}(4\mathrm{πp}^2\mathrm{dp})$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>&#xA0;</mo><mo>&#x222B;</mo><mfrac><msup><mfenced open="|" close="|"><mi>p</mi></mfenced><mn>2</mn></msup><mi>E</mi></mfrac><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mo>&#x222B;</mo><mfrac><msup><mfenced open="|" close="|"><mi>p</mi></mfenced><mn>2</mn></msup><mi>E</mi></mfrac><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac><mo>(</mo><mn>4</mn><msup><mi>&#x3C0;p</mi><mn>2</mn></msup><mi>dp</mi><mo>)</mo></math>\;\int\frac{\left|p\right|^2}E\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p\;=\;\int\frac{\left|p\right|^2}E\frac{\partial f(\boldsymbol p,t)}{\partial E}(4\mathrm{πp}^2\mathrm{dp}), we have substituted Volume element d³p with 4 πp²dp.

= $4\mathrm\pi\int\frac{\left|\mathrm p\right|^4}E\frac{\partial f(\boldsymbol p,t)}{\partial E}\;=\;4\mathrm\pi\int_{}^{}\frac{\mathrm p^4}{\mathrm E}\;dp\frac{\partial f}{\partial p}\frac{\partial p}{\partial E}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mn>4</mn><mi mathvariant="normal">&#x3C0;</mi><mo>&#x222B;</mo><mfrac><msup><mfenced open="|" close="|"><mi mathvariant="normal">p</mi></mfenced><mn>4</mn></msup><mi>E</mi></mfrac><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mn>4</mn><mi mathvariant="normal">&#x3C0;</mi><msubsup><mo>&#x222B;</mo><mrow/><mrow/></msubsup><mfrac><msup><mi mathvariant="normal">p</mi><mn>4</mn></msup><mi mathvariant="normal">E</mi></mfrac><mo>&#xA0;</mo><mi>d</mi><mi>p</mi><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi></mrow><mrow><mo>&#x2202;</mo><mi>p</mi></mrow></mfrac><mfrac><mrow><mo>&#x2202;</mo><mi>p</mi></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac></math>4\mathrm\pi\int\frac{\left|\mathrm p\right|^4}E\frac{\partial f(\boldsymbol p,t)}{\partial E}\;=\;4\mathrm\pi\int_{}^{}\frac{\mathrm p^4}{\mathrm E}\;dp\frac{\partial f}{\partial p}\frac{\partial p}{\partial E}​

In the Relativistic Limit, E² = m²c⁴ + p²c² and taking c =1, we differentiate it inorder to get, EdE =pdp. So,

$\frac{\partial p}{\partial E}=\frac{E}{p}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mrow><mo>&#x2202;</mo><mi>p</mi></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac><mo>=</mo><mfrac><mi>E</mi><mi>p</mi></mfrac></math>\frac{\partial p}{\partial E}=\frac Ep .

Hence the Equation becomes,

$4\pi \int \frac{{|p|}^4}{E}\frac{\partial f(\mathbf{p},t)}{\partial p}\frac{\partial p}{\partial f(\mathbf{p},t)}=4\pi \int {\left|\mathrm{p}\right|}^3\frac{\partial f(\mathbf{p},t)}{\partial p}dp$<math xmlns="http://www.w3.org/1998/Math/MathML"><mn>4</mn><mi mathvariant="normal">&#x3C0;</mi><mo>&#x222B;</mo><mfrac><msup><mfenced open="|" close="|"><mi>p</mi></mfenced><mn>4</mn></msup><mi>E</mi></mfrac><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>p</mi></mrow></mfrac><mo>&#xA0;</mo><mfrac><mrow><mo>&#x2202;</mo><mi>p</mi></mrow><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow></mfrac><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mn>4</mn><mi mathvariant="normal">&#x3C0;</mi><mo>&#x222B;</mo><msup><mfenced open="|" close="|"><mi mathvariant="normal">p</mi></mfenced><mn>3</mn></msup><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>p</mi></mrow></mfrac><mo>&#xA0;</mo><mi>d</mi><mi>p</mi></math>4\mathrm\pi\int\frac{\left|p\right|^4}E\frac{\partial f(\boldsymbol p,t)}{\partial p}\;\frac{\partial p}{\partial f(\boldsymbol p,t)}\;=\;4\mathrm\pi\int\left|\mathrm p\right|^3\frac{\partial f(\boldsymbol p,t)}{\partial p}\;dp

We use multiplication Rule for differentiating it and we get,

We didnot consider the first term in 31c , because it reduces to zero on integration.

Using Equation 31d, in $\;\;-\;\frac g{(2\mathrm\pi)^3}\int\frac HE\left|{}_{}p\right|^2\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>&#xA0;</mo><mo>&#xA0;</mo><mo>-</mo><mo>&#xA0;</mo><mfrac><mi>g</mi><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mo>&#x222B;</mo><mfrac><mi>H</mi><mi>E</mi></mfrac><msup><mfenced open="|" close="|"><mmultiscripts><mi>p</mi><mprescripts/><mrow/><none/></mmultiscripts></mfenced><mn>2</mn></msup><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></math>\;\;-\;\frac g{(2\mathrm\pi)^3}\int\frac HE\left|{}_{}p\right|^2\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p, we get

​$\;\;-\;\frac g{(2\mathrm\pi)^3}\int\frac HE\left|{}_{}p\right|^2\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>&#xA0;</mo><mo>&#xA0;</mo><mo>-</mo><mo>&#xA0;</mo><mfrac><mi>g</mi><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mo>&#x222B;</mo><mfrac><mi>H</mi><mi>E</mi></mfrac><msup><mfenced open="|" close="|"><mmultiscripts><mi>p</mi><mprescripts/><mrow/><none/></mmultiscripts></mfenced><mn>2</mn></msup><mfrac><mrow><mo>&#x2202;</mo><mi>f</mi><mo>(</mo><mi mathvariant="bold-italic">p</mi><mo>,</mo><mi>t</mi><mo>)</mo></mrow><mrow><mo>&#x2202;</mo><mi>E</mi></mrow></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></math>\;\;-\;\frac g{(2\mathrm\pi)^3}\int\frac HE\left|{}_{}p\right|^2\frac{\partial f(\boldsymbol p,t)}{\partial E}d^3p= ​ $\;\;\;\frac{3gH}{(2\mathrm\pi)^3}\int4\mathrm{πp}^2\mathrm f(\mathbf p,\mathrm t)\mathrm{dp}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>&#xA0;</mo><mo>&#xA0;</mo><mo>&#xA0;</mo><mfrac><mrow><mn>3</mn><mi>g</mi><mi>H</mi></mrow><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mo>&#x222B;</mo><mn>4</mn><msup><mi>&#x3C0;p</mi><mn>2</mn></msup><mi mathvariant="normal">f</mi><mo>(</mo><mi mathvariant="bold">p</mi><mo>,</mo><mi mathvariant="normal">t</mi><mo>)</mo><mi>dp</mi></math>\;\;\;\frac{3gH}{(2\mathrm\pi)^3}\int4\mathrm{πp}^2\mathrm f(\mathbf p,\mathrm t)\mathrm{dp}

= $\;\;\;3H\int\frac g{(2\mathrm\pi)^3}^{}\mathrm f(\mathbf p,\mathrm t)d^3p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>&#xA0;</mo><mo>&#xA0;</mo><mo>&#xA0;</mo><mn>3</mn><mi>H</mi><mo>&#x222B;</mo><msup><mfrac><mi>g</mi><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mrow/></msup><mi mathvariant="normal">f</mi><mo>(</mo><mi mathvariant="bold">p</mi><mo>,</mo><mi mathvariant="normal">t</mi><mo>)</mo><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></math>\;\;\;3H\int\frac g{(2\mathrm\pi)^3}^{}\mathrm f(\mathbf p,\mathrm t)d^3p

= $3Hn$<math xmlns="http://www.w3.org/1998/Math/MathML"><mn>3</mn><mi>H</mi><mi>n</mi></math>3Hn

The above equations are for phase space distribution f_χ lets denote the WIMP particle with χ.

Therefore the LHS of the Equation 30 becomes $\overset{}{\frac g{\left(2\mathrm\pi\right)^3}\int{\overset{}{\frac{\widehat L\left[f_\chi\right]}Ed^3p_\chi}}_{}}\;=\;\frac{\operatorname dn_\chi}{\operatorname dt}+\;3Hn_\chi$<math xmlns="http://www.w3.org/1998/Math/MathML"><mover><mrow><mfrac><mi>g</mi><msup><mfenced><mrow><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi></mrow></mfenced><mn>3</mn></msup></mfrac><mo>&#x222B;</mo><msub><mover><mrow><mfrac><mrow><mover><mi>L</mi><mo>^</mo></mover><mfenced open="[" close="]"><msub><mi>f</mi><mi>&#x3C7;</mi></msub></mfenced></mrow><mi>E</mi></mfrac><msup><mi>d</mi><mn>3</mn></msup><msub><mi>p</mi><mi>&#x3C7;</mi></msub></mrow><mrow/></mover><mrow/></msub></mrow><mrow/></mover><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><mrow><mo>d</mo><msub><mi>n</mi><mi>&#x3C7;</mi></msub></mrow><mrow><mo>d</mo><mi>t</mi></mrow></mfrac><mo>+</mo><mo>&#xA0;</mo><mn>3</mn><mi>H</mi><msub><mi>n</mi><mi>&#x3C7;</mi></msub></math>\overset{}{\frac g{\left(2\mathrm\pi\right)^3}\int{\overset{}{\frac{\widehat L\left[f_\chi\right]}Ed^3p_\chi}}_{}}\;=\;\frac{\operatorname dn_\chi}{\operatorname dt}+\;3Hn_\chi​

Now, let's take the right side of the Equation 30, i.e., $\frac g{\left(2\mathrm\pi\right)^3}\;\int_{}^{}\frac{C\left[f\right]}Ed^3p$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mi>g</mi><msup><mfenced><mrow><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi></mrow></mfenced><mn>3</mn></msup></mfrac><mo>&#xA0;</mo><msubsup><mo>&#x222B;</mo><mrow/><mrow/></msubsup><mfrac><mrow><mi>C</mi><mfenced open="[" close="]"><mi>f</mi></mfenced></mrow><mi>E</mi></mfrac><msup><mi>d</mi><mn>3</mn></msup><mi>p</mi></math>\frac g{\left(2\mathrm\pi\right)^3}\;\int_{}^{}\frac{C\left[f\right]}Ed^3p .But before that let's think of a simple annihilation χ χ* ↔ψ ψ* .Now, the term we considered is equal to ,

∫dΠ_ψ dΠ_ψ* dΠ_χ dΠ_χ* (2π)⁴ δ⁴( P_χ+ P_χ* - P_ψ -P_ψ*) [(M²_{χ χ* ↔ψ ψ*} ) f_χ f_χ* (1±f_ψ )(1±f_ψ*) - (Μ²_{ψ ψ* ↔χ χ*})f_ψ f_ψ* (1±f_χ )(1±f_χ* )]

where dΠi = $\frac{g_i}{(2\mathrm\pi)^3}\frac{d^3p_i}{2E_i}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><msub><mi>g</mi><mi>i</mi></msub><mrow><mo>(</mo><mn>2</mn><mi mathvariant="normal">&#x3C0;</mi><msup><mo>)</mo><mn>3</mn></msup></mrow></mfrac><mfrac><mrow><msup><mi>d</mi><mn>3</mn></msup><msub><mi>p</mi><mi>i</mi></msub></mrow><mrow><mn>2</mn><msub><mi>E</mi><mi>i</mi></msub></mrow></mfrac></math>\frac{g_i}{(2\mathrm\pi)^3}\frac{d^3p_i}{2E_i} and the delta function ensures that the momentum and the energy in the scattering processes are conserved. M is the amplitude or also called as the transition matrix which gives the probability of the particular process to happen.

And also the Phase space distribution function $f(p)=\frac1{exp({\displaystyle\frac{E\;-\mu}T})\;\pm\;\;1}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>f</mi><mo>(</mo><mi>p</mi><mo>)</mo><mo>=</mo><mfrac><mn>1</mn><mrow><mi>e</mi><mi>x</mi><mi>p</mi><mo>(</mo><mstyle displaystyle="true"><mfrac><mrow><mi>E</mi><mo>&#xA0;</mo><mo>-</mo><mi>&#x3BC;</mi></mrow><mi>T</mi></mfrac></mstyle><mo>)</mo><mo>&#xA0;</mo><mo>&#xB1;</mo><mo>&#xA0;</mo><mo>&#xA0;</mo><mn>1</mn></mrow></mfrac></math>f(p)=\frac1{exp({\displaystyle\frac{E\;-\mu}T})\;\pm\;\;1},

The ± sign is for the type of particles we are dealing with i.e., + for bosons and - for fermions.

Now in the range of temperature we are interested in say T_i is very less compared to the difference between E and μ so at this condition (1± f) ≈ 1.

We also make following assumptions like :-

i) No CP Variance in the DM sector. Hence the Μ²_{ψ ψ* ↔χ χ*} = M²_{χ χ* ↔ψ ψ*} = M².

ii) Distribution of particles and antiparticles are almoslt same. Hence, μ = 0. So, f_ψ= f_ψ* and f_χ = f_χ*.

iii)We also assume that the particles which are being created upon the annihilation of massive WIMP particle are lighter and their interactions between the other particles(SM Particles) in the thermal bath is more compared to that of WIMPs. Hence, these particles stay in thermal equilibrium for a long time and also chemical potential is zero. Hence, the f_ψ = exp(-E_ψ/T) and f_ψ*= exp(-E_ψ*/T). This observation tells that the particles started to follow the Boltzmann Statistics.By using the Energy Conservation , we get

So,$\Rightarrow\;f_\psi f_{\psi\ast\;}=\;{f^{EQ}}_\chi{f^{EQ}}_{\chi\ast}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mo>&#x21D2;</mo><mo>&#xA0;</mo><msub><mi>f</mi><mi>&#x3C8;</mi></msub><msub><mi>f</mi><mrow><mi>&#x3C8;</mi><mo>*</mo><mo>&#xA0;</mo></mrow></msub><mo>=</mo><mo>&#xA0;</mo><msub><msup><mi>f</mi><mrow><mi>E</mi><mi>Q</mi></mrow></msup><mi>&#x3C7;</mi></msub><msub><msup><mi>f</mi><mrow><mi>E</mi><mi>Q</mi></mrow></msup><mrow><mi>&#x3C7;</mi><mo>*</mo></mrow></msub></math>\Rightarrow\;f_\psi f_{\psi\ast\;}=\;{f^{EQ}}_\chi{f^{EQ}}_{\chi\ast}

iv) We make use of the relation that $f\propto n$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>f</mi><mo>&#x221D;</mo><mi>n</mi></math>f\propto n.

Considering these assumptions the equation modifies to ∫dΠ_ψ dΠ_ψ* dΠ_χ dΠ_χ* (2π)⁴ δ⁴( P_χ+ P_χ* - P_ψ -P_ψ*) [(M²) f_χ f_χ* - (Μ²)f_ψ f_ψ* ],

= ∫dΠ_ψ dΠ_ψ* dΠ_χ dΠ_χ* (2π)⁴ δ⁴( P_χ+ P_χ* - P_ψ -P_ψ*) [(M²) f_χ f_χ* - (Μ²)f_χ ^EQf_χ*^EQ]

= ∫dΠ_ψ dΠ_ψ* dΠ_χ dΠ_χ* (2π)⁴ δ⁴( P_χ+ P_χ* - P_ψ -P_ψ*) M²[ f_χ f_χ* - f_χ^EQf_χ*^EQ]

= ∫dΠ_ψ dΠ_ψ* dΠ_χ dΠ_χ* (2π)⁴ δ⁴( P_χ+ P_χ* - P_ψ -P_ψ*) M²f_χ^EQf_χ*^EQ[ f_χ f_χ* /f_χ^EQf_χ*^EQ - 1]

= ∫dΠ_ψ dΠ_ψ* dΠ_χ dΠ_χ* (2π)⁴ δ⁴( P_χ+ P_χ* - P_ψ -P_ψ*) M²f_χ^EQf_χ*^EQ[ n_χ n_χ* /n_χ^EQn_χ*^EQ - 1] → Equation a

The thermally averaged annihilation cross section is

<σ_{χ χ* →ψ ψ*} v_Mol> = (n_χ^EQn_χ*^EQ)^-1 ∫dΠ_ψ dΠ_ψ* dΠ_χ dΠ_χ* (2π)⁴ δ⁴( P_χ+ P_χ* - P_ψ -P_ψ*) M²f_χ^EQf_χ*^EQ

If we include the above definition of cross section and consider the assumptions . The Equation a becomes

∫dΠ_ψ dΠ_ψ* dΠ_χ dΠ_χ* (2π)⁴ δ⁴( P_χ+ P_χ* - P_ψ -P_ψ*) M²f_χ^EQf_χ*^EQ[ n_χ n_χ* /n_χ^EQn_χ*^EQ - 1]

= <σ_{χ χ* →ψ ψ*} v_Mol> [(n_χ)² - (n_χ^EQ)²]

But it is not mandatory that this is the only possible annihilation. Indeed many annihilation channels( χ χ* ↔F, where F is the final stateand it need not be a 2 body state ) are possible and thus we take the averaged crosssection <σ_A v_Mol>.

Therefore, we have modified the Boltzmann Equation into

where v = v_Mol = √(v₁ - v₂) - (v₁ × v₂)

We know that the Abundance or the number of particles per comoving volume is given by Y ≡ n_χ/s, where n is the actual number density of the particle and s is the entropy density. Hence, Y is the actual abundance and Y_EQ is the abundance if the expansion is not considered so it becomes Y_EQ ≡ n_χ^EQ/s .

Boltzmann Equation in terms of Abundance Y

In the previous section we have mentioned that the abundance Y is the number of particles present in a comoving volume. And also we make a statement that the total entropy per comoving volume is a constant, for which the proof will be given in the appendix section.

But we also know that

So, substituting equation 35 in 34, we get,

We now define a Dimensionless parameter x = m/T , where T is the temperature and m is the mass of the particle.

Now, $\begin{array}{l}a\;\propto\;\lambda\;,\;\;\lambda\;\propto\;\frac1T\Rightarrow\;\;\;a\propto\frac1T\;\\\end{array}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>a</mi><mo>&#xA0;</mo><mo>&#x221D;</mo><mo>&#xA0;</mo><mi>&#x3BB;</mi><mo>&#xA0;</mo><mo>,</mo><mo>&#xA0;</mo><mo>&#xA0;</mo><mi>&#x3BB;</mi><mo>&#xA0;</mo><mo>&#x221D;</mo><mo>&#xA0;</mo><mfrac><mn>1</mn><mi>T</mi></mfrac><mo>&#x21D2;</mo><mo>&#xA0;</mo><mo>&#xA0;</mo><mo>&#xA0;</mo><mi>a</mi><mo>&#x221D;</mo><mfrac><mn>1</mn><mi>T</mi></mfrac><mo>&#xA0;</mo><mspace linebreak="newline"/></math>\begin{array}{l}a\;\propto\;\lambda\;,\;\;\lambda\;\propto\;\frac1T\Rightarrow\;\;\;a\propto\frac1T\;\\\end{array}

So,Entropy Density s ∝ a^-3 ∝ T³ . Therefore,

​$m\frac d{dt}\left(\frac ax\right)\;=\frac{\displaystyle d}{\displaystyle dt}\left(\frac{\displaystyle a}{\displaystyle x}\right)\;=\;0$<math xmlns="http://www.w3.org/1998/Math/MathML"><mi>m</mi><mfrac><mi>d</mi><mrow><mi>d</mi><mi>t</mi></mrow></mfrac><mfenced><mfrac><mi>a</mi><mi>x</mi></mfrac></mfenced><mo>&#xA0;</mo><mo>=</mo><mfrac><mstyle displaystyle="true"><mi>d</mi></mstyle><mstyle displaystyle="true"><mi>d</mi><mi>t</mi></mstyle></mfrac><mfenced><mfrac><mstyle displaystyle="true"><mi>a</mi></mstyle><mstyle displaystyle="true"><mi>x</mi></mstyle></mfrac></mfenced><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mn>0</mn></math>m\frac d{dt}\left(\frac ax\right)\;=\frac{\displaystyle d}{\displaystyle dt}\left(\frac{\displaystyle a}{\displaystyle x}\right)\;=\;0 and $\begin{array}{l}\frac{x\dot a-a\dot x}{x^2}=0\;and\;finally,\;Hx\;=\;\frac{dx}{dt}\\\end{array}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mfrac><mrow><mi>x</mi><mover><mi>a</mi><mo>&#x2D9;</mo></mover><mo>-</mo><mi>a</mi><mover><mi>x</mi><mo>&#x2D9;</mo></mover></mrow><msup><mi>x</mi><mn>2</mn></msup></mfrac><mo>=</mo><mn>0</mn><mo>&#xA0;</mo><mi>a</mi><mi>n</mi><mi>d</mi><mo>&#xA0;</mo><mi>f</mi><mi>i</mi><mi>n</mi><mi>a</mi><mi>l</mi><mi>l</mi><mi>y</mi><mo>,</mo><mo>&#xA0;</mo><mi>H</mi><mi>x</mi><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mfrac><mrow><mi>d</mi><mi>x</mi></mrow><mrow><mi>d</mi><mi>t</mi></mrow></mfrac><mspace linebreak="newline"/></math>\begin{array}{l}\frac{x\dot a-a\dot x}{x^2}=0\;and\;finally,\;Hx\;=\;\frac{dx}{dt}\\\end{array}

Here, H is the Hubble parameter as a function of temperature . As the freezeout takes place in a Radiation Dominated Universe(Much Before BBN), H(T) ≈ $1.66g_\ast^{1/2}\frac{T^2}{M_{pl}},\;where\;M_{pl}\;=\;1.22\;\times10^{19}\;GeV$<math xmlns="http://www.w3.org/1998/Math/MathML"><mn>1</mn><mo>.</mo><mn>66</mn><msup><msub><mi>g</mi><mo>*</mo></msub><mrow><mn>1</mn><mo>/</mo><mn>2</mn></mrow></msup><mfrac><msup><mi>T</mi><mn>2</mn></msup><msub><mi>M</mi><mrow><mi>p</mi><mi>l</mi></mrow></msub></mfrac><mo>,</mo><mo>&#xA0;</mo><mi>w</mi><mi>h</mi><mi>e</mi><mi>r</mi><mi>e</mi><mo>&#xA0;</mo><msub><mi>M</mi><mrow><mi>p</mi><mi>l</mi></mrow></msub><mo>&#xA0;</mo><mo>=</mo><mo>&#xA0;</mo><mn>1</mn><mo>.</mo><mn>22</mn><mo>&#xA0;</mo><mo>&#xD7;</mo><msup><mn>10</mn><mn>19</mn></msup><mo>&#xA0;</mo><mi>G</mi><mi>e</mi><mi>V</mi></math>1.66g_\ast^{1/2}\frac{T^2}{M_{pl}},\;where\;M_{pl}\;=\;1.22\;\times10^{19}\;GeV

Let H(m) = $1.66g_\ast^{1/2}\frac{m^2}{M_{pl}}$<math xmlns="http://www.w3.org/1998/Math/MathML"><mn>1</mn><mo>.</mo><mn>66</mn><msubsup><mi>g</mi><mo>*</mo><mrow><mn>1</mn><mo>/</mo><mn>2</mn></mrow></msubsup><mfrac><msup><mi>m</mi><mn>2</mn></msup><msub><mi>M</mi><mrow><mi>p</mi><mi>l</mi></mrow></msub></mfrac></math>1.66g_\ast^{1/2}\frac{m^2}{M_{pl}}, therefore H(m) =H(T)x²