So I wanted know how to theoretically calculate this phase space invariant (equation $3.31a$ )$R$ in our universe (FLRW metric) during the cosmological nucleosysthesis:

$$R = \int_{p} \frac{\mathcal{N}(x,p,t)}{\mathcal{E}} dp_x dp_y dp_z$$

where $\mathcal{N}(x,p,t)$ is the phase space distribution of matter in our universe, $\mathcal{E}= p_0$ is the energy and $p_x$, $p_y$ and $p_z$ are the spatial momentum.

Feel free to assume extra assumptions (but assume the formation of the molecules has taken place)?


2 Answers 2


I'm going to start by using some relations that appear explicitly in the text you linked. Specifically:

\begin{equation} \mathcal{N} = \frac{g_{s}}{h^{3}} \eta \end{equation}

where, for the era of nucleosynthesis, the Universe is matter-dominated ($\rho \propto a^{-3} \, , \, P=0$ with $a$ being the time-dependent scale factor), and thus the particle statistics is dominated by non-relativistic protons, electrons and Helium nuclei ($g_{s} = 2s +1$ with $s$ being typically $1/2$). Since matter is non-relativistic and at that energy level behave roughly classically, the particles follow Boltzmann statistics as such:

\begin{equation} \eta \approx e^{-\beta (E - \mu)} = e^{\beta \mu} \, e^{-\beta p^{2}/2m} \end{equation}

where we defined $\beta \equiv (k_{B}T)^{-1}$ and approximated the kinetic energy to the Newtonian form. The integral for the field $R$ now reads:

\begin{equation} R = \frac{g_{s}}{h^{3}} e^{\beta \mu} \int \frac{e^{-\beta p^{2}/2m}}{\frac{p^{2}}{2m} + m} \, d\mathcal{V}_{p} = \frac{4\pi g_{s}}{h^{3}} e^{\beta \mu} \int _{0} ^{\Lambda} \frac{p^{2} e^{-\beta p^{2}/2m}}{\frac{p^{2}}{2m} + m} \, dp \end{equation}

The parameter $\Lambda$ is a cutoff we introduce since the momentum of each particle has a form of cap due to the assumption they behave non-relativistically. Next we can use the trick:

\begin{equation} p^{2} e^{-\beta p^{2}/2m} = -2m \frac{\partial }{\partial \beta} e^{-\beta p^{2}/2m} \end{equation}

so the integral now reads:

\begin{equation} R = -4m^{2} \frac{4\pi g_{s}}{h^{3}} e^{\beta \mu} \frac{\partial}{\partial \beta} \int _{0} ^{\Lambda} \frac{e^{-\beta p^{2}/2m}}{p^{2} + 2m^{2}} \, dp \end{equation}

The integral now can be calculated using the Cauchy integral formula by choosing an appropriate contour that circumvents the two imaginary poles $p = \pm i\sqrt{2} m$:

\begin{equation} \oint \frac{e^{-\beta p^{2}/2m}}{(p+i\sqrt{2}m)(p-i\sqrt{2}m)} \, dp = 2\pi i \left [ \frac{e^{\beta m}}{i2\sqrt{2}m} - \frac{e^{\beta m}}{-i2\sqrt{2}m} \right ] = \sqrt{2} \pi \frac{e^{\beta m}}{m} \end{equation}

The effect of the derivative with respect to $\beta$ only acts on the exponential which yields $m \, e^{\beta m}$, so the final result is:

\begin{equation} R = -\frac{16\pi g_{s} m^{2}}{\sqrt{2}h^{3}} e^{\beta \tilde{\mu}} \end{equation}

where as defined in the text $\tilde{\mu} = \mu + m$. The result comes out negative, but since the quantity $R$ is a scalar field and not a positive-definite quantity like energy or the number of particles, I expect the result should be allowed.


As an addendum and in response to your comment about not taking the non-relativistic limit, let me provide some additional calculations. Let's assume that we are at some point during the radiation-dominated era of the Universe where most particles behave relativistically. This time they can be either photons/other massless bosons or hard fermions. The coupling between them is very weak due to being of high energy, so $\beta \mu \ll 1$ and the particle statistics can be both Bose-Einstein or Fermi-Dirac:

\begin{equation} \eta = \frac{1}{e^{\beta(E - \mu)} \pm 1} \approx \frac{1}{e^{\beta p} \pm 1} \end{equation}

where we further assumed $p^{2} \gg m^{2} \Rightarrow E \approx \mathcal{E} \approx p$. The scalar field $R$ will now be:

\begin{equation} R = \frac{g_{s}}{h^{3}} \int \frac{1}{e^{\beta(E - \mu)} \pm 1} \frac{1}{p} \, d\mathcal{V}_{p} = \frac{4\pi g_{s}}{h^{3}} \int _{\lambda} ^{\Lambda} \frac{p}{e^{\beta(E - \mu)} \pm 1} \, dp \end{equation}

As I mentioned in the comments, in addition to the upper cutoff $\Lambda$ (which you may freely now assume to be very large/approaching infinity), you impose a lower cutoff $\lambda$ which restricts the momentum of your particles to energies that would qualify as in the relativistic regime. The result of the integral had we taken no limits would be:

\begin{equation} \int \frac{p}{e^{\beta p} \pm 1} \, dp = \begin{cases} &\frac{p \, \log{(1 - e^{-\beta p})}}{\beta } - \frac{\mathrm{Li}_{2}(e^{-\beta p})}{\beta ^{2} } \, , \qquad \text{(EB)} \\ &\frac{\mathrm{Li} _{2}(e^{-\beta p})}{\beta ^{2} } - \frac{p \, \log{(1 + e^{-\beta p})}}{\beta } \, , \qquad \text{(FD)} \\ \end{cases} \end{equation}

where $\mathrm{Li}_{n}$ is the polylogarithmic function. For the upper limit $\Lambda \rightarrow +\infty$:

\begin{equation} \mathrm{Li}_{2}(0) = 0 \end{equation}

\begin{equation} \lim _{p \rightarrow +\infty} {p \, \log{(1 \pm e^{-\beta p})}} = 0 \end{equation}

So all is left are the terms with the lower limit supplanted:

\begin{equation} \int _{\lambda}^{\Lambda} \frac{p}{e^{\beta p} \pm 1} \, dp = \begin{cases} &-\frac{\lambda \, \log{(1 - e^{-\beta \lambda})}}{\beta } + \frac{\mathrm{Li}_{2}(e^{-\beta \lambda})}{\beta ^{2} } \, , \qquad \text{(EB)} \\ &-\frac{\mathrm{Li} _{2}(e^{-\beta \lambda})}{\beta ^{2} } + \frac{\lambda \, \log{(1 + e^{-\beta \lambda})}}{\beta } \, , \qquad \text{(FD)} \\ \end{cases} \end{equation}

For relatively large $\lambda$, Fermi-Dirac statistics leads to $R \approx 0$, but negative due to having a higher power of $\beta$ in the denominator in the first term (remember that $T$ in the radiation-dominated era is much higher than in the matter-dominated one). For Bose-Einstein statistics, had the lower limit been 0, the logarithm would tend towards $-\infty$, which is where I presume the divergence you found occurs. But for relatively large values $\lambda$, the logarithm is now regulated and is expected to be something close to 0. The first term will also be very small, so to determine whether $R$ is positive or negative you'd have to input more precise values of your parameters.

  • $\begingroup$ Yes, I also got a negative answer which was precisely what was tripping me since the text says: "This is the sum, over all particles in a unit 3-volume, of the inverse energy. " $\endgroup$ Commented Jun 21, 2022 at 14:00
  • $\begingroup$ The phrasing is just misguiding, I believe. The phrase "over all particles" implies some particle statistics. This is precisely the part that turns the integral negative. The exposition earlier on that it is a scalar field should tip you off that negative values aren't necessarily prohibited. But I do agree that the way it was presented is confusing. $\endgroup$
    – rhomaios
    Commented Jun 21, 2022 at 14:15
  • $\begingroup$ Also why is it that I get an integral that does not converge when I do not resort to non-relativistic approximation? $\endgroup$ Commented Jun 21, 2022 at 14:52
  • $\begingroup$ Did you impose an appropriate lower cutoff to your integral? For relativistic particles the lower integral limit of the radial coordinate in the phase space cannot be 0 or something very small, otherwise it would belie the idea that your particles are relativistic. $\endgroup$
    – rhomaios
    Commented Jun 21, 2022 at 16:23

After pondering for a while, I thought it might be worth sharing a different angle to this problem. So let's start with the stress energy tensor $T^{\mu \nu}$ of a perfect fluid:

$$T^{\mu \nu} = \left(\rho + \frac{p}{c^2} \right) U^{\mu} U^\nu + p g^{\mu \nu}, $$

where $p$ is pressure, $\rho$ is density, $U^\mu$ is four velocity and $g^{\mu \nu}$ is the metric. Consider the trace of this quantity:

$$T \equiv T^{i}_i = (\rho + \frac{p}{c^2}) c^2+ 4p =3p + \rho c^2 $$

Now, we know for a particular regime this parameter is approximately constant. Hence, we have a conserved quantity or equation of state:

$$ 3p + \rho c^2 = \text{constant}$$

However, we also know from the linked text's equation $(3.31 c)$ in the question. The stress energy tensor is given by:

$$ T^{\mu \nu} = \int \mathcal{N}(x,p,t) p^\mu \otimes p^\nu \frac{d V_p}{E}$$

Taking the trace as mentioned before we get and using $m$ as mass:

$$ \text{constant}=3p + \rho c^2= m^2 c^4 \int \mathcal{N}(x,p,t) \frac{d V_p}{E} $$


$$ \int \mathcal{N}(x,p,t) \frac{d V_p}{E} =\frac{1}{m^2c^4} (3p + \rho{c^2} ) = \frac{\text{constant}}{m^2c^4} $$

  • $\begingroup$ This is a good idea, but you made a mistake in your third equation. The covariant derivative of the stress-energy tensor is 0, but that doesn't mean the trace is a constant. In fact, one of the Friedmann equations reads $\frac{\ddot{a}}{a} = -\frac{4\pi G}{3} (\rho +3P)$. You can argue the deceleration is approximately constant for some regime, but in general it is not since $a \propto t^{2/3}$ for the matter-dominated era and $a \propto t^{1/2}$ for the radiation-dominated one. But since $\ddot{a} >0$ and there's a negative sign included, it does give negative values as we have found so far. $\endgroup$
    – rhomaios
    Commented Jun 23, 2022 at 10:17
  • $\begingroup$ Another minor correction, but in your fifth and sixth equations, there should be $m^{2} c^{4}$ as derived from the trace $p^{\mu} p_{\mu}$. You can see it from the dimensional analysis as well ($R$ should have dimensions of $[\text{mass}]^{2}$ in natural units). The major problem however is when you want to deal with photons or other massless bosons, since then the aforementioned trace is 0, so you can't evaluate $R$ with this method. $\endgroup$
    – rhomaios
    Commented Jun 23, 2022 at 10:22
  • $\begingroup$ Thanks I'll add the corrections :) $\endgroup$ Commented Jun 23, 2022 at 10:25
  • $\begingroup$ If anything I can see why the lecturer side stepped this problematic invariant > _ < $\endgroup$ Commented Jun 23, 2022 at 10:40
  • $\begingroup$ Not so problematic as it is physically useless; it even says so in the text that it has very little application in physics. To me it's kind of logical that the results seem so strange and unintuitive, since it doesn't give you any quantity you'd reproduce or otherwise come across in realistic physical calculations. $\endgroup$
    – rhomaios
    Commented Jun 23, 2022 at 11:02

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.