Many Body in General Relativity: A Thermal Equivalence Principle

1. Introduction

An interesting problem in physics is to study the properties of a (quantum) many body system at low (non-zero) temperature on a curved surface. For example colloidal particles may be adsorbed or confined on a substrate with nonzero curvature, be it the wall of a porous material, or a membrane, a vesicle, a micelle for example made of ampiphilic surfactant molecules such as lipids, or a biological membrane, or the surface of a large solid particle, or an interface in an oil-water emulsion [1]. For a fluid of ⁴He atoms it would be interesting to study the superfluid fraction. For a fluid of electrons it would be interesting to study the superfluidity. Moreover it would be interesting to study the properties of the electron plasma on a sphere in presence of a magnetic field.

One important point to discuss is whether the space in which the particles live is exactly two dimensional, as it happens in the satirical novella of Edwin Abbott Abbott [2], or if it can be treated as quasi two dimensional. There is a profound difference between the two scenarios to the point that the form of the interaction between the particles also changes. For example for colloidal particles one may choose the polarizable hard sphere pair interaction or for the fluid of helium atoms one may use the Lennard Jones pair potential, but the distance between the two interacting particles may be chosen either as the geodesic distance between them or the Euclidean distance in the three dimensional space where the surface is embedded. For the electron gas the Coulomb pair potential as a solution to the Poisson equation has different forms in two or three dimensions and in general depends on the metric of the curved surface.

These properties can be studied exactly with the path integral (Monte Carlo) method and these studies certainly enrich the knowledge on many bodies in (quantum) general relativity [3,4,5,6,7]. Not even the two body problem can be treated analytically in general relativity [8]. The problem of gravitating many bodies should be separated by the problem of many bodies with non-gravitational interactions in general relativity. In fact mass curves spacetime through the Einstein field equations and gravitating bodies will behave as free particles on that curved spacetime, whereas non-gravitational interactions produce particles accelerations on the spacetime. So being able to treat many (quantum) bodies on a curved surface would be an important step forward for the much more complicated problem of gravitating many (quantum) bodies in general relativity.

We find it of fundamental importance issuing a bridge between the two scientific communities of the exact simulations of a many body (quantum) system and of general relativity. We foresee an important progress in the physics of (quantum) gravitating many body systems beyond the simple ideal gases or hydrodynamic systems that are usually treated [9,10]. We here review the physics of many bodies in the context of general relativity. Starting from the stress energy tensor for one body, for a swarm of bodies, for a perfect fluid, we review relativistic hydrodynamics, kinetic theory, and statistical physics of N identical bodies. We conclude our excursion with a thermal equivalence principle in physics.

In this work we consider spacetime as a smooth manifold $\mathcal{M}$ of dimension d and metric tensor $\mathit{g}\mathit{g}$ with covariant components $g_{\alpha \beta }$. We will denote with an arrow over a bold face letter the corresponding 4-vector and with just the bold face symbol the corresponding 3-dimensional vector. Greek indexes run over the d spacetime dimensions. Roman indexes run only over the $d-1$ space dimensions. We use Einstein summation convention of tacitly assuming a sum over repeated indexes. We will assume the speed of light $c=1$ throughout.

One Particle

For one body of mass m we have a self gravitating system with a stress energy tensor given by

$$\begin{array}{c}\hfill T^{\alpha \beta }\left(\overrightarrow{\mathit{x}}\right)=m\int u^{\alpha }{u}^{\beta }{\delta }^{\left(4\right)}(\overrightarrow{\mathit{x}}-\overrightarrow{\mathit{z}}\left(\tau \right))d\tau ,\end{array}$$

(1)

where $\tau$ is the body proper time, $d\overrightarrow{\mathit{z}}/d\tau =\overrightarrow{\mathit{u}}=(\gamma ,\gamma \mathit{v})$ with $u^0=dt/d\tau =\gamma ={(1-v^2)}^{-1/2}$ and

$$\begin{array}{c}\hfill T^{\alpha \beta }\left(\overrightarrow{\mathit{x}}\right)=m{\displaystyle \frac{u^{\alpha }{u}^{\beta }}{u^0}}{\delta }^{\left(3\right)}(\mathit{x}-\mathit{z}\left(t\right)),\end{array}$$

(2)

where the body is at $\mathit{z}\left(t\right)$ with velocity $\mathit{v}\left(t\right)$ at time t.

Swarm of Particles

For a swarm of N bodies all of the same mass m and $\mathit{v}$

$$\begin{array}{ccc}\hfill T^{\alpha \beta }\left(\overrightarrow{\mathit{x}}\right)& =& m{u}^{\alpha }{u}^{\beta }\sum _{i=1}^N\int {\delta }^{\left(4\right)}(\overrightarrow{\mathit{x}}-{\overrightarrow{\mathit{z}}}_i\left({\tau }_i\right))d{\tau }_i\hfill \\ & =& m{\displaystyle \frac{u^{\alpha }{u}^{\beta }}{u^0}}\sum _{i=1}^N{\delta }^{\left(3\right)}(\mathit{x}-{\mathit{z}}_i\left(t\right))\hfill \\ & =& m{u}^{\alpha }{u}^{\beta }n,\hfill \end{array}$$

(3)

where

$$\begin{array}{c}\hfill n={\displaystyle \frac{1}{u^0}}\sum _{i=1}^N{\delta }^{\left(3\right)}(\mathit{x}-{\mathit{z}}_i\left(t\right)),\end{array}$$

(4)

is the proper number density of bodies measured in a comoving frame where $\overrightarrow{\mathit{u}}=(1,\mathbf{0})$.

Perfect Fluid

For a perfect fluid of proper number density n of non interacting bodies all of the same mass m and $v=\left|\mathit{v}\right|$ but isotropic velocity profile $\mathit{v}=v\mathit{n}$

$$\begin{array}{c}\hfill T^{\alpha \beta }=\chi {\langle u^{\alpha }{u}^{\beta }{\rangle }_{\mathit{n}}}_{,}\end{array}$$

(5)

so that $T^{\alpha \beta }=0$ for $\alpha \ne \beta$ and $T^{00}=\chi {\gamma }^2$. Since $T^{00}=\rho =n\left(\gamma m\right)$ is the energy density of the fluid we require $\chi =mn/\gamma$. Then

$$\begin{array}{ccc}\hfill T^{ij}& =& \chi {\gamma }^2{v}^2\langle n^i{n}^j{\rangle }_{\mathit{n}}\hfill \\ & =& \chi {\gamma }^2{v}^2{\displaystyle \frac{1}{3}}{\delta }^{ij}\hfill \\ & =& n\left(\gamma m\right)v^2{\displaystyle \frac{1}{3}}{\delta }^{ij}\hfill \\ & =& p{\delta }^{ij},\hfill \end{array}$$

(6)

where $\delta$ is a Kronecker delta and in the second equality we used isotropy of n and

$$\begin{array}{c}\hfill \left\{\begin{array}{c}\rho =n\left(m\gamma \right)\hfill \\ p={\displaystyle \frac{1}{3}}\rho v^2\hfill \end{array}\right.\end{array}$$

(7)

are respectively the mass density and pressure in the isotropic frame of the fluid. Summarizing

$$\begin{array}{c}\hfill T^{\alpha \beta }=(\rho +p)u^{\alpha }{u}^{\beta }+p{\eta }^{\alpha \beta },\end{array}$$

(8)

where $\left|\right|{\eta }^{\alpha \beta }\left|\right|=\mathrm{diag}\{-1,1,1,1\}$ is the metric in Minkowski spacetime. For photons $v=1$ and $p=\rho /3$. For $v\ll 1$, $\rho =nm(1+v^2/2+\dots )$, and $p\approx nm{v}^2/3=(2/3)(\rho -nm)=(2/3)ϵ$, where $ϵ=(3/2)k_{B}T$ is the internal energy of a monatomic ideal gas in thermal equilibrium at a temperature T, $k_B$ is Boltzmann constant, and $p=n{k}_{B}T$ is the ideal gas equation of state.

2. Hydrodynamics

Hydrodynamics concerns itself with the study of the motion of fluids (liquids and gases). Since the phenomena considered in fluid dynamics are macroscopic, a fluid is regarded as a continuous medium. Therefore when we speak of the “point” of a fluid (or of an infinitesimal volume of it) we mean not a single molecule of the fluid but a volume element still containing very many molecules but yet small compared with the volume of the whole fluid.

2.1. Newtonian

A mathematical description of the state of a moving fluid consists in specifying the fluid velocity $\mathit{v}=\mathit{v}(t,\mathit{x})$ and any two thermodynamic functions pertaining to the fluid, for instance the pressure $p=p(t,\mathit{x})$ and the density $\rho =\rho (t,\mathit{x})$, from which one can determine all other thermodynamic quantities. These 5 quantities are functions of the coordinates $\mathit{x}=(x,y,z)$ and of time t. Once again we stress that a point $\mathit{r}$ in space at a given time t refers to a fixed point and not to specific particles of the fluid. From Chapter I of Ref. [11] we find

$$\begin{array}{ccc}& & {\displaystyle \frac{\partial \rho }{\partial t}}+\nabla \left(\rho \mathit{v}\right)=0,\hfill \end{array}$$

(9)

$$\begin{array}{ccc}& & {\displaystyle \frac{\partial \mathit{v}}{\partial t}}+(\mathit{v}·\nabla )\mathit{v}=-{\displaystyle \frac{1}{\rho }}\nabla p,\hfill \end{array}$$

(10)

$$\begin{array}{ccc}& & {\displaystyle \frac{\partial s}{\partial t}}+(\mathit{v}·\nabla )s=0,\hfill \end{array}$$

(11)

where the first equation is the continuity equation, the second is the Euler equation, and the third one is the equation for the adiabatic flow in which $s=s(t,\mathit{x})$ is the entropy per particle.

From the first law of thermodynamics follows

$$\begin{array}{ccc}& & dϵ=Tds-pd(m/\rho ),\hfill \end{array}$$

(12)

$$\begin{array}{ccc}& & ϵ=ϵ(\rho ,s),\hfill \end{array}$$

(13)

$$\begin{array}{ccc}& & p={\displaystyle \frac{{\rho }^2}{m}}{\left.{\displaystyle \frac{\partial ϵ}{\partial \rho }}\right|}_s,\hfill \end{array}$$

(14)

$$\begin{array}{ccc}& & T={\left.{\displaystyle \frac{\partial ϵ}{\partial s}}\right|}_{\rho },\hfill \end{array}$$

(15)

where $ϵ$ is the internal energy per particle. Equations (14) and (15) can be considered as algebraic relations for the right hand side of Equations (10) and (11) respectively.

For an ideal gas $ϵ=ϵ\left(T\right)$ and for a monatomic gas

$$\begin{array}{c}\hfill s=k_{B}ln(T^{3/2}m/\rho )+constant.\end{array}$$

(16)

2.2. Relativistic

We will work in a Local Lorentz Frame (LLF). Recalling that the stress energy tensor is divergenceless, from the stress energy tensor of a perfect fluid (8) we find

$$\begin{array}{ccc}\hfill 0=T^{\alpha \beta }{,}_{\beta }& =& {(\rho +p)}_{,\beta }{u}^{\alpha }{u}^{\beta }+(\rho +p)u^{\alpha }{}_{,\beta }u^{\beta }+(\rho +p)u^{\alpha }{u}^{\beta }{}_{,\beta }+p_{,\beta }{\eta }^{\alpha \beta }\hfill \\ & =& {\displaystyle \frac{d(\rho +p)}{d\tau }}{u}^{\alpha }+(\rho +p)a^{\alpha }+(\rho +p)u^{\alpha }{u}^{\beta }{}_{,\beta }+p_{,}{}^{\alpha },\hfill \end{array}$$

(17)

where the comma stands for a partial derivative. Multiplying by $\overrightarrow{\mathit{u}}$ and recalling that $\overrightarrow{\mathit{u}}·\overrightarrow{\mathit{a}}=0$ we find

$$\begin{array}{c}\hfill {\displaystyle \frac{d\rho }{d\tau }}=-(\rho +p)\overrightarrow{\nabla }\overrightarrow{\mathit{u}},\end{array}$$

(18)

which is the relativistic continuity expression which extends Equation (9).

To find the extension of the Euler equation we introduce the projector tensor

$$\begin{array}{ccc}\hfill P^{\alpha \beta }& =& {\eta }^{\alpha \beta }+u^{\alpha }{u}^{\beta }for\overrightarrow{\mathit{u}}\mathrm{timelike}\overrightarrow{\mathit{u}}·\overrightarrow{\mathit{u}}=-1\hfill \\ \hfill P^{\alpha \beta }& =& {\eta }^{\alpha \beta }-n^{\alpha }{n}^{\beta }for\overrightarrow{\mathit{u}}\mathrm{spacelike}\overrightarrow{\mathit{u}}·\overrightarrow{\mathit{u}}=+1\hfill \end{array}$$

Then

$$\begin{array}{c}\hfill 0=P_{\alpha \gamma }{T}^{\alpha \beta }{}_{,\beta }=(\rho +p)a_{\gamma }+P_{\alpha \gamma }{p}_{,}{}^{\alpha },\end{array}$$

(19)

or

$$\begin{array}{c}\hfill (\rho +p)\overrightarrow{\mathit{a}}=-\overrightarrow{\nabla }p-\overrightarrow{\mathit{u}}{\displaystyle \frac{dp}{d\tau }},\end{array}$$

(20)

which is the relativistic Euler equation which extends Equation (9).

It is easy to see that in the non relativistic limit $\overrightarrow{\mathit{u}}=(\gamma ,\gamma \mathit{v})\approx (1,\mathit{v})$ with $v\ll 1$ and $p\ll \rho$, Equation (18) reduces to Equation (9) and Equation (20) reduces to Equation (10) [12] 1

Let us now discuss the continuity Equation (18). First of all we observe that the mass density is not conserved $d\rho /d\tau \ne 0$. But the baryon, lepton, charge, … numbers are conserved. For example if we call $n=N/V$ the baryon number density in the rest frame of the fluid with N baryons in a volume V, N is certainly constant but V will change, so that

$$\begin{array}{c}\hfill 0={\displaystyle \frac{dN}{d\tau }}={\displaystyle \frac{d\left(nV\right)}{d\tau }},\end{array}$$

(21)

but $(dV/d\tau )/V=\overrightarrow{\nabla }\overrightarrow{\mathit{u}}$ (see Ex. 22.1 in Ref. [10]). So

$$\begin{array}{ccc}\hfill 0& =& {\displaystyle \frac{1}{V}}{\displaystyle \frac{d\left(nV\right)}{d\tau }}\hfill \\ & =& {\displaystyle \frac{dn}{d\tau }}+n\overrightarrow{\nabla }\overrightarrow{\mathit{u}}\hfill \\ & =& \overrightarrow{\mathit{u}}·\overrightarrow{\nabla }n+n\overrightarrow{\nabla }\overrightarrow{\mathit{u}}\hfill \\ & =& \overrightarrow{\nabla }\left(n\overrightarrow{\mathit{u}}\right),\hfill \end{array}$$

(22)

where we may define the divergenceless current density

$$\begin{array}{c}\hfill \overrightarrow{\mathit{J}}=n\overrightarrow{\mathit{u}}.\end{array}$$

(23)

Let us now discuss the thermodynamics. The second law tells that $ds/d\tau \ge 0$ where s is the entropy per baryon. The first law becomes

$$\begin{array}{c}\hfill d(\rho /n)=-pd(1/n)+Tds,\end{array}$$

(24)

or

$$\begin{array}{c}\hfill d\rho ={\displaystyle \frac{\rho +p}{n}}dn+nTds,\end{array}$$

(25)

which is the relativistic extension of Equation (12). In this equation the differential d can be substituted either with an exterior derivative $\tilde{\mathit{d}}$, with a gradient $\overrightarrow{\nabla }$, or with a directional derivative ${\overrightarrow{\nabla }}_{\overrightarrow{\mathit{u}}}=u^{\alpha }\partial /\partial x^{\alpha }=d/d\tau$. Given an equation of state $\rho =\rho (n,s)$ we will have

$$\begin{array}{ccc}& & p=n{\left.{\displaystyle \frac{\partial \rho }{\partial n}}\right|}_s-\rho ,\hfill \end{array}$$

(26)

$$\begin{array}{ccc}& & T={\displaystyle \frac{1}{n}}{\left.{\displaystyle \frac{\partial \rho }{\partial s}}\right|}_n,\hfill \end{array}$$

(27)

which are the relativistic extensions of Equations (14) and (15).

It is easy to show that a perfect fluid flow is adiabatic. From the relativistic continuity Equations (18) and (22) follows

$$\begin{array}{c}\hfill {\displaystyle \frac{d\rho }{d\tau }}={\displaystyle \frac{\rho +p}{n}}{\displaystyle \frac{dn}{d\tau }}.\end{array}$$

(28)

Then from the relativistic first thermodynamic Equation (25) follows

$$\begin{array}{c}\hfill {\displaystyle \frac{ds}{d\tau }}=0.\end{array}$$

(29)

2.2.1. Shock Wave

Consider a homogeneous, static, perfect fluid. A sound wave in the fluid is an adiabatic perturbation. The speed of sound is

$$\begin{array}{c}\hfill v_s^2={\left.{\displaystyle \frac{\partial p}{\partial \rho }}\right|}_s\end{array}$$

(30)

Expand

$$\begin{array}{ccc}\hfill \rho & =& {\rho }_0+{\rho }_1,\hfill \\ \hfill p& =& p_0+p_1,\hfill \\ \hfill n& =& n_0+n_1,\hfill \end{array}$$

where ${\rho }_0,p_0,n_0$ are constant in space (uniform fluid) and in time (static fluid) and ${\rho }_1,p_1,n_1$ are small perturbations. Taking $\overrightarrow{\mathit{u}}=(1,{\mathit{v}}_1)$ with $v_1\ll 1$ we find from the continuity Equation (18)

$$\begin{array}{c}\hfill {\displaystyle \frac{\partial {\rho }_1}{\partial t}}=-({\rho }_0+p_0)\nabla {\mathit{v}}_1,\end{array}$$

(31)

and from the spatial part of Euler Equation (20)

$$\begin{array}{c}\hfill ({\rho }_0+p_0){\displaystyle \frac{\partial {\mathit{v}}_1}{\partial t}}=-\nabla p_1,\end{array}$$

(32)

where we neglect the last term $\overrightarrow{\mathit{u}}dp/d\tau =\overrightarrow{\mathit{u}}{u}^{\alpha }\partial p/\partial x^{\alpha }$ because an infinitesimal of second order and $\partial p_0/\partial t=0$. Therefore putting together Equations (31) and (32) we find

$$\begin{array}{c}\hfill {\displaystyle \frac{{\partial }^2{\rho }_1}{d{t}^2}}=-({\rho }_0+p_0)\nabla {\displaystyle \frac{\partial {\mathit{v}}_1}{\partial t}}={\nabla }^2{p}_1.\end{array}$$

(33)

In a perfect fluid $p=p(\rho ,T)$ so that $p({\rho }_0+{\rho }_1,T)=p({\rho }_0,T)+\partial p({\rho }_0,T){/\partial \rho |}_s{\rho }_1=p_0+p_1$ with $p_1=v_c^2{\rho }_1$ and we finally find

$$\begin{array}{c}\hfill {\displaystyle \frac{{\partial }^2{\rho }_1}{d{t}^2}}=v_s^2{\nabla }^2{\rho }_1,\end{array}$$

(34)

which is the shock wave equation.

2.2.2. Bernoulli Equation

Consider a steady, adiabatic flow of a perfect fluid. Since in a steady state $\partial p/\partial t=0$ from the relativistic Euler Equation (20) follows

$$\begin{array}{c}\hfill (\rho +p){\displaystyle \frac{d{u}^0}{d\tau }}=-u^0{\displaystyle \frac{dp}{d\tau }},\end{array}$$

(35)

So

$$\begin{array}{ccc}\hfill {\displaystyle \frac{d}{d\tau }}\left(u^0{\displaystyle \frac{\rho +p}{n}}\right)& =& {\displaystyle \frac{d{u}^0}{d\tau }}{\displaystyle \frac{\rho +p}{n}}+u^0{\displaystyle \frac{d}{d\tau }}\left({\displaystyle \frac{\rho +p}{n}}\right)\hfill \\ & =& -{\displaystyle \frac{u^0}{n}}{\displaystyle \frac{dp}{d\tau }}+u^0{\displaystyle \frac{d}{d\tau }}\left({\displaystyle \frac{\rho +p}{n}}\right)\hfill \\ & =& {\displaystyle \frac{u^0}{n}}\left[{\displaystyle \frac{d\rho }{d\tau }}-{\displaystyle \frac{\rho +p}{n}}{\displaystyle \frac{dn}{d\tau }}\right]=0,\hfill \end{array}$$

(36)

where in the second equality we used Equation (35) and in the last equality we used the relativistic first law of thermodynamics Equation (25). We then conclude that

$$\begin{array}{c}\hfill u^0{\displaystyle \frac{\rho +p}{n}}isconstantalongthefluidflowlines.\end{array}$$

(37)

In the Newtonian limit $\overrightarrow{\mathit{u}}=(\gamma ,\gamma \mathit{v})\approx (1,\mathit{v})$ with $v\ll 1$, $u^0=\gamma ={(1-v^2)}^{-1/2}\approx 1+v^2/2$, and $P\ll \rho$. We can take $\rho ={\rho }_0(1+\pi )$ and $(\rho +p)/n\approx {\rho }_0(1+\pi +p/{\rho }_0)/n$, so that

$$\begin{array}{c}\hfill {\displaystyle \frac{1}{2}}{v}^2+\pi +p/{\rho }_{0}isconstantalongthefluidflowlines,\end{array}$$

(38)

where the sum of the last two terms is the enthalpy.

3. Kinetic Theory Approach

The kinetic theory approach is based on a one body distribution function [10,13].

Distribution Function

We will construct a Lorentz invariant phase space distribution function as a number density of particles in phase space

$$\begin{array}{c}\hfill f={\displaystyle \frac{d\mathcal{N}}{d\mathit{x}d\mathit{p}}},\int fd\mathit{x}d\mathit{p}=N,\end{array}$$

(39)

for a fluid of N bodies, where $d\mathit{x}=d{x}^{1}d{x}^{2}d{x}^3$ and $d\mathit{p}=d{p}^{1}d{p}^{2}d{p}^3$. So that $f/N$ can be considered as a probability distribution function. We will now prove that f as defined above is a Lorentz invariant distribution. We start defining a proper 3-volume. The 4-volume $d^4\mathsf{\Omega }=d{x}^{0}d{x}^{1}d{x}^{2}d{x}^3$ is invariant under a Lorentz transformation. Dividing by $d\tau$ we find another Lorentz invariant

$$\begin{array}{c}\hfill dV=u^{0}d{x}^{1}d{x}^{2}d{x}^3.\end{array}$$

(40)

Then we want to define a 3-volume element in momentum space. The 4-volume $d^{4}p=d{p}^{0}d{p}^{1}d{p}^{2}d{p}^3$ is invariant. Since $p^0=\sqrt{{\mathit{p}}^2+m^2}$ we will define

$$\begin{array}{ccc}\hfill d\mathsf{\Pi }& =& \int d^{4}p\delta \left(\sqrt{-\overrightarrow{\mathit{p}}·\overrightarrow{\mathit{p}}}-m\right)\hfill \\ & =& {\displaystyle \frac{\sqrt{{\left(p^0\right)}^2-{\mathit{p}}^2}}{p^0}}d{p}^{1}d{p}^{2}d{p}^3\hfill \\ & =& {\displaystyle \frac{m}{p^0}}d{p}^{1}d{p}^{2}d{p}^3.\hfill \end{array}$$

(41)

And

$$\begin{array}{c}\hfill dVd\mathsf{\Pi }=d{x}^{1}d{x}^{2}d{x}^{3}d{p}^{1}d{p}^{2}d{p}^3,\end{array}$$

(42)

is Lorentz invariant.

We will now prove conservation of volume in phase space in curved spacetime (see Liouville theorem in BOX 22.6 of Ref. [10]). Consider a very small bundle of identical particles that move through curved spacetime on a neghboring geodesics. We want to prove that $d\left(dVd\mathsf{\Pi }\right)/d\lambda =0$ where $\lambda$ is an affine parameter along the central geodesic of the bundle. Given any function of phase space $g(\overrightarrow{\mathit{x}},\overrightarrow{\mathit{p}})$, if $m\ne 0$ 2 take $\tau =a\lambda +b$ for arbitrary a and b. Then

$$\begin{array}{c}\hfill {\displaystyle \frac{dg(\overrightarrow{\mathit{x}},\overrightarrow{\mathit{p}})}{d\lambda }}={\displaystyle \frac{\partial g}{\partial x^{\alpha }}}{\displaystyle \frac{d{x}^{\alpha }}{d\lambda }}+{\displaystyle \frac{\partial g}{\partial p^{\alpha }}}{\displaystyle \frac{d{p}^{\alpha }}{d\lambda }},\end{array}$$

(43)

on a geodesic $d{p}^{\alpha }/d\tau =0$ so

$$\begin{array}{c}\hfill {\displaystyle \frac{dg(\overrightarrow{\mathit{x}},\overrightarrow{\mathit{p}})}{d\lambda }}={\displaystyle \frac{\partial g}{\partial x^{\alpha }}}{p}^{\alpha }{\displaystyle \frac{a}{m}},\end{array}$$

(44)

and for $g=dVd\mathsf{\Pi }$ and $p^{\alpha }=md{x}^{\alpha }/d\tau$

$$\begin{array}{c}\hfill {\displaystyle \frac{dg}{d\lambda }}=a{\displaystyle \frac{dg}{d\lambda }},\end{array}$$

(45)

for any a. so $dg/d\lambda =0$. Since $d\mathcal{N}$ and $dVd\mathsf{\Pi }$ are unchanged then also $f=d\mathcal{N}/dVd\mathsf{\Pi }$ is unchanged

$$\begin{array}{c}\hfill {\displaystyle \frac{df}{d\lambda }}=0.\end{array}$$

(46)

This equation is at the heart of the collisionless Boltzmann equation and the Vlasov equation. All these approximate theories are valid at sufficiently low density. Whereas for the full Boltzmann equation [13,14] one has

$$\begin{array}{c}\hfill {\displaystyle \frac{df}{d\lambda }}={\left({\displaystyle \frac{\partial f}{\partial \lambda }}\right)}_{\mathrm{collisions}}.\end{array}$$

(47)

This is the most famous of all Kinetic equations and was obtained by Boltzmann more than a century ago.

The phase space probability density of a system in thermodynamic equilibrium at an inverse temperature $\beta =1/k_{B}T$ with $k_B$ Boltzmann constant, is not an explicit function of proper time. We shall use the symbol $f_0$ to denote the equilibrium probability density.

Ideal Gas

For an ideal gas, i.e. a fluid of non interacting many identical bodies, from §37 & 52 & 53 of Ref. [15] we know that, on a comoving frame with $\overrightarrow{\mathit{u}}=(1,\mathbf{0})$, we can write

$$\begin{array}{c}\hfill f_0(\overrightarrow{\mathit{x}},\overrightarrow{\mathit{p}})={\displaystyle \frac{d\mathcal{N}}{dVd\mathsf{\Pi }}}={\displaystyle \frac{g}{h^3}}{\displaystyle \frac{1}{e^{-\beta (\overrightarrow{\mathit{p}}·\overrightarrow{\mathit{u}}+\mu )}-\epsilon }},\end{array}$$

(48)

where h is Planck constant, $\mu$ is the chemical potential, $g=2J+1$ is the spin J degeneracy (2 polarizations for photons), and

$$\begin{array}{c}\hfill \epsilon =\left\{\begin{array}{cc}+1\hfill & Bose-Einsteinstatistics\hfill \\ 0\hfill & Maxwell-Boltzmannstatistics\hfill \\ -1\hfill & Fermi-Diracstatistics\hfill \end{array}\right.\end{array}$$

(49)

Where $\epsilon =0$ takes care of the statistics for classical distinguishable bodies at sufficiently high temperature. At sufficiently low temperature a different statistics must be devised, in which the mean occupation number of the various quantum states of bodies are not assumed small. The statistics, however, differs according to the type of many body wave function by which the gas is described. These functions must be either symmetrical or antisymmetrical with respect to interchange of any pair of particles (see §61 in Ref. [16]). The former case occurring for bodies with integral spin, bosons $\epsilon =+1$, and the latter case for those of half-integral spin, fermions $\epsilon =-1$.

Moments of the distribution function of the ideal gas

Next we can take moments of $f_0$ with respect to $\overrightarrow{\mathit{p}}$

$$\begin{array}{ccc}\hfill \int f_0{p}^{\mu }d\mathsf{\Pi }& =& J^{\mu },\hfill \end{array}$$

(50)

$$\begin{array}{ccc}\hfill \int f_0{p}^{\mu }{p}^{\nu }d\mathsf{\Pi }& =& T^{\mu \nu },\hfill \end{array}$$

(51)

where here $d\mathsf{\Pi }=d{}_{\stackrel{++}{=}}/p^0$ with $p^0=\sqrt{{\mathit{p}}^2+m^2}$. For $\overrightarrow{\mathit{u}}=(1,\mathbf{0})$, from Equations (23) and (8), we must have

$$\begin{array}{ccc}\hfill J^{\mu }& =& n{u}^{\mu },\hfill \end{array}$$

(52)

$$\begin{array}{ccc}\hfill T^{\mu \nu }& =& (\rho +p)u^{\mu }{u}^{\nu }+p{\eta }^{\mu \nu }.\hfill \end{array}$$

(53)

So

$$\begin{array}{ccc}\hfill n& =& -J^{\mu }{u}_{\mu }=-\int f_0{p}^{\mu }{u}_{\mu }d\mathsf{\Pi }=\int f_{0}d\mathit{p}\hfill \\ & =& {\displaystyle \frac{g}{h^3}}{\int }_0^{\infty }{\displaystyle \frac{4\pi p^{2}dp}{e^{\beta [\sqrt{{\mathit{p}}^2+m^2}-\mu ]}-\epsilon }}.\hfill \end{array}$$

(54)

Introduce the following change of variables

$$\begin{array}{c}\hfill \left\{\begin{array}{c}p=msinh\chi \hfill \\ \overline{\beta }=m\beta \hfill \end{array}\right.\end{array}$$

(55)

so that from Equation (54) we find

$$\begin{array}{c}\hfill n={\displaystyle \frac{4\pi g{m}^3}{h^3}}{\int }_0^{\infty }{\displaystyle \frac{{sinh}^2\chi cosh\chi d\chi }{e^{[\overline{\beta }cosh\chi -\beta \mu ]}-\epsilon }}.\end{array}$$

(56)

For the pressure

$$\begin{array}{ccc}\hfill p& =& {\displaystyle \frac{1}{3}}(u_{\mu }{u}_{\nu }+{\eta }_{\mu \nu })T^{\mu \nu }\hfill \\ & =& {\displaystyle \frac{1}{3}}\int f_0{\mathit{p}}^{2}d\mathsf{\Pi }\hfill \\ & =& {\displaystyle \frac{1}{3}}\int f_0{\mathit{p}}^2{\displaystyle \frac{d\mathit{p}}{p^0}}\hfill \\ & =& {\displaystyle \frac{4\pi g{m}^4}{3h^3}}{\int }_0^{\infty }{\displaystyle \frac{{sinh}^4\chi d\chi }{e^{[\overline{\beta }cosh\chi -\beta \mu ]}-\epsilon }}.\hfill \end{array}$$

(57)

Also

$$\begin{array}{ccc}\hfill \rho -3p& =& -T^{\alpha }{}_{\alpha }=m^2\int f_0{\displaystyle \frac{d\mathit{p}}{p^0}}\hfill \\ & =& {\displaystyle \frac{4\pi g{m}^4}{h^3}}{\int }_0^{\infty }{\displaystyle \frac{{sinh}^2\chi d\chi }{e^{[\overline{\beta }cosh\chi -\beta \mu ]}-\epsilon }}.\hfill \end{array}$$

(58)

Maxwell-Boltzmann statistics ($\epsilon =0$) [17]

From Ref. [18] we learn that

$$\begin{array}{ccc}\hfill K_n\left(\overline{\beta }\right)& =& {\displaystyle \frac{{\overline{\beta }}^n}{(2n-1)!!}}{\int }_0^{\infty }d\chi {sinh}^{2n}\chi e^{-\overline{\beta }cosh\chi }\hfill \end{array}$$

(59)

$$\begin{array}{ccc}& =& {\displaystyle \frac{{\overline{\beta }}^{n-1}}{(2n-3)!!}}{\int }_0^{\infty }d\chi {sinh}^{2n-2}\chi cosh\chi e^{-\overline{\beta }cosh\chi },\hfill \end{array}$$

(60)

where $K_n$ is a modified Bessel function of the second kind and in the second equality we performed an integration by parts. The asymptotic behaviors of the modified Bessel function are as follows

$$\begin{array}{ccc}\hfill K_n\left(\overline{\beta }\right)& =& \sqrt{{\displaystyle \frac{\pi }{2\overline{\beta }}}}{e}^{-\overline{\beta }}\left[1+{\displaystyle \frac{4n^2-1}{8\overline{\beta }}}+O\left({\overline{\beta }}^{-2}\right)\right]\overline{\beta }\gg 1,\hfill \end{array}$$

(61)

$$\begin{array}{ccc}\hfill K_n\left(\overline{\beta }\right)& =& {\displaystyle \frac{(n-1)!}{{\overline{\beta }}^n}}\left[2^{n-1}-{\displaystyle \frac{2^{n-3}{\overline{\beta }}^2}{n-1}}+O\left({\overline{\beta }}^3\right)\right]\overline{\beta }\ll 1.\hfill \end{array}$$

(62)

We then find

$$\begin{array}{ccc}\hfill n& =& a{K}_2\left(\overline{\beta }\right)/\overline{\beta },\hfill \end{array}$$

(63)

$$\begin{array}{ccc}\hfill p& =& am{K}_2\left(\overline{\beta }\right)/{\overline{\beta }}^2=n{k}_{B}T,\hfill \end{array}$$

(64)

$$\begin{array}{ccc}\hfill \rho -3p& =& am{K}_1\left(\overline{\beta }\right)/\overline{\beta },\hfill \end{array}$$

(65)

where $a=4\pi g{m}^3{e}^{\beta \mu }/h^3$. Note that the ideal gas equation of state (64) is a relativistic invariant.

For the internal energy per particle we then find

$$\begin{array}{c}\hfill u\left(T\right)={\displaystyle \frac{\rho }{n}}=m{\displaystyle \frac{K_1\left(\overline{\beta }\right)}{K_2\left(\overline{\beta }\right)}}+3k_{B}T=\left\{\begin{array}{cc}{\displaystyle m\left[1+\frac{3}{2}\frac{k_{B}T}{m}+\dots \right]}\hfill & \overline{\beta }\gg 1,\hfill \\ 3k_{B}T\hfill & \overline{\beta }\ll 1,\hfill \end{array}\right.\end{array}$$

(66)

where we used the asymptotic expansions (61) and (62).

For the ratio of the specific heats $\gamma \left(T\right)=c_p/c_v$ we then find

$$\begin{array}{ccc}\hfill \gamma \left(T\right)& =& {\displaystyle \frac{{\left.\frac{du}{dT}\right|}_p}{{\left.\frac{du}{dT}\right|}_v}}=1+{\displaystyle \frac{k_B}{{\left.\frac{du}{dT}\right|}_v}}=\left\{\begin{array}{cc}{\displaystyle 5/3}\hfill & \overline{\beta }\gg 1,\hfill \\ 4/3\hfill & \overline{\beta }\ll 1.\hfill \end{array}\right.\hfill \end{array}$$

(67)

Quantum Statistics ($\epsilon =\pm 1$) [19]

For the Bose-Einstein and the Fermi-Dirac statistics and the zero temperature, $\beta \to \infty$, limit see Ref. [19].

4. Statistical Mechanics Approach

The statistical mechanics approach is based on a many body distribution function [10,13].

Thermal Equilibrium of the Many Bodies

The aim of equilibrium statistical mechanics is to calculate observable properties of a system of interest either as averages over a phase trajectory (the method of Boltzmann), or as averages over an ensemble of systems, each of which is a replica of the system of interest (the method of Gibbs). In Gibbs formulation of statistical mechanics the equilibrium probability distribution for the systems of N identical bodies of the ensemble is described by ${\rho }_0^{\left(N\right)}$, a phase space probability density, in a $6N$ dimensional phase space ${\prod }_{i=1}^{N}d{\mathit{x}}_{i}d{\mathit{p}}_i=d^{N}V{d}^N\mathsf{\Pi }$, in the classical case or a density matrix in the quantum case. Here we will only consider the more general quantum case that reduces to the classical case at high temperature. We will then have [13,20,21]

$$\begin{array}{ccc}\hfill {\rho }_0^{\left(N\right)}& =& exp\left(-{\int }_{\partial \mathsf{\Omega }}{T}^{\mu \nu }{\beta }_{\nu }d{S}_{\mu }\right),\hfill \end{array}$$

(68)

$$\begin{array}{ccc}\hfill Z_N& =& \mathrm{tr}\left({\rho }_0^{\left(N\right)}\right),\hfill \end{array}$$

(69)

where $\partial \mathsf{\Omega }$ is a general, arbitrary, spacelike hypersurface bounding the 4-volume $\mathsf{\Omega }$ and $\overrightarrow{\beta }\left(\mathit{x}\right)$ is a 4-vector such that $\beta =\sqrt{{\beta }_{\mu }{\beta }^{\mu }}$ and as usual $1/k_B\beta \left(\mathit{x}\right)=T\left(\mathit{x}\right)$ the invariant absolute temperature, i.e. the temperature measured by a comoving thermometer. $Z_N$ is the canonical partition function where $\mathrm{tr}(\dots )$ denotes a trace that requires a path integral in position representation [20] foot:identical. And for example $\int T^{00}dV=\mathcal{H}$ with $\mathcal{H}=\mathcal{K}+\mathcal{V}$ the Hamiltonian operator of the fluid where $\mathcal{K}$ is the kinetic energy operator of the N bodies. The covariant form of Equation (68) of the equilibrium statistical operator was first used by Weldon [22] for the Belinfante symmetrized stress energy tensor.

It is then possible to define the n-body reduced equilibrium distribution functions as

$$\begin{array}{ccc}\hfill f_0^{\left(n\right)}({\overrightarrow{\mathit{x}}}_1,\dots ,{\overrightarrow{\mathit{x}}}_n)& =& 〈{\left[\prod _{i=1}^n\sum _{j=1}^N\delta ({\overrightarrow{\mathit{x}}}_i-{\overrightarrow{\mathit{r}}}_j)\right]}_{\mathrm{DP}}〉\hfill \end{array}$$

(70)

$$\begin{array}{ccc}& =& {\displaystyle \frac{1}{Z_N}}{\left[\prod _{i=1}^n\sum _{j=1}^N\delta ({\overrightarrow{\mathit{x}}}_i-{\overrightarrow{\mathit{r}}}_j)\right]}_{\mathrm{DP}}{\rho }_0^{\left(N\right)}(\left\{{\overrightarrow{\mathit{r}}}_k\right\},\left\{{\overrightarrow{\mathit{r}}}_k\right\};\beta )d^{4N}\overrightarrow{\mathit{r}},\hfill \end{array}$$

(71)

where the thermal average of an operator $\mathcal{O}$ is $\langle \mathcal{O}\rangle =\mathrm{tr}\left({\rho }_0^{\left(N\right)}\mathcal{O}\right)/Z_N$, $d^{4N}\overrightarrow{\mathit{r}}={\prod }_{i=1}^{N}d{\overrightarrow{\mathit{r}}}_i$, ${\rho }_0^{\left(N\right)}(\left\{{\overrightarrow{\mathit{r}}}_k\right\},\left\{{\overrightarrow{\mathit{r}}}_k^{\prime }\right\};\beta )$ is the position representation of the density matrix (68) at an inverse temperature $\beta$ that results from a path integral [20], 3 and the subscript DP means that only the products of Dirac delta functions relative to Different Particles should be considered. Here $f_0^{\left(1\right)}=\int f_{0}d\mathsf{\Pi }$ where $f_0$ is the one body distribution function of the previous Section 3.

In a recent project [23] we studied an electron gas at low temperatures, the Jellium, on the surface of a sphere through the path integral Monte Carlo method. A unit sphere is the surface 4 of constant positive scalar curvature 2. In particular we noticed as the simulation “speed” of the path in a neighborhood of the poles diminishes. This is a consequence of the hairy ball theorem, according to which her Euler class is the obstruction to her tangent planes, the tangent bundle, 5 having always a non vanishing fiber, or hair, for any section 6. The theorem was first proven by Henri Poincaré for the sphere in 1885 [25], and extended to higher even dimensions in 1912 by Luitzen Egbertus Jan Brouwer [26]. The theorem has been expressed colloquially as “you can’t comb a hairy ball flat without creating a cowlick” or “you can’t comb the hair on a coconut”. If z is a continuous function that assigns a vector in the three dimensional space to every point $\mathcal{P}$ on a sphere such that $z\left(\mathcal{P}\right)$ is always tangent to the sphere at $\mathcal{P}$, then there is at least one pole, a point where the field vanishes, i.e. a $\mathcal{P}$ such that $z\left(\mathcal{P}\right)=0$. Every zero of a vector field has a (non-zero) index 7, and it can be shown that the sum of all of the indexes at all of the zeros must be two, because the Euler characteristic of the sphere is two. Therefore, there must be at least one zero. This is a consequence of the Poincaré-Hopf theorem. The theorem was proven for two dimensions by Henri Poincaré and later generalized to higher dimensions by Heinz Hopf [28]. In particular we see how, even a single free particle have a path which will be subject to some anisotropy due to the effective potential induced by the curvature of the sphere. This effect was studied in Refs. [23,29].

Thermal Equilibrium of the Metric Tensor

A different story is to move the temperature from the stress energy tensor to the metric tensor as is done in Refs. [3,7] and in the trilogy [4,5,6] also applied to study the vacuum in cosmic space in Ref. [30]. That is, to move the statistical physics description from the right hand side of Einstein field equations

$$\begin{array}{c}\hfill G_{\mu \nu }={\displaystyle \frac{8\pi G}{c^4}}{T}_{\mu \nu },\end{array}$$

(72)

to the left hand side. Of course the two descriptions has to give the same picture. In Ref. [4] we took the statistical average of the trace of Einstein field equations

$$\begin{array}{c}\hfill {\langle -R\rangle }_g=\kappa {\langle T_{\mu }{}^{\mu }\rangle }_t.\end{array}$$

(73)

where $\kappa =8\pi G/c^4$ and $R=-{G_{\mu }}^{\mu }$ is the scalar curvature, ${\langle \dots \rangle }_g$ is a statistical average on the metric tensor, and ${\langle \dots \rangle }_t$ is a time average. On the right hand side, replacing the time average with an ensemble average, we find [21]

$$\begin{array}{c}\hfill {\langle z_{\mu }{T}^{\mu \nu }\rangle }_t={\displaystyle \frac{\mathrm{tr}\left({\rho }_0^{\left(N\right)}{z}_{\mu }{T}^{\mu \nu }\right)}{Z_N}}=-{\displaystyle \frac{\delta }{\delta {\beta }_{\nu }\left(\mathit{x}\right)}}ln{Z}_N.\end{array}$$

(74)

In the above formula, while the left hand side depends on a arbitrary vector $\overrightarrow{\mathit{z}}$, the right hand side is not manifestly dependent on it. In fact, the functional derivative of $Z_N$ of Equation (69) includes a hidden dependence on the normal vector as the functional derivation implies the choice of a measure, hence of a hypersurface and a corresponding normal vector. We will also have

$$\begin{array}{c}\hfill {\langle T_{\mu }{}^{\mu }\rangle }_t=-{\displaystyle \frac{\delta }{\delta \beta }}ln{Z}_N.\end{array}$$

(75)

Then the virial theorem of Equation (73) can be rewritten as

$$\begin{array}{c}\hfill {\langle R\rangle }_g=\kappa {\displaystyle \frac{\delta }{\delta \beta }}ln{Z}_N,\end{array}$$

(76)

where $\mathrm{tr}(\dots )$ denotes a trace, and $\mathsf{\Lambda }=\sqrt{4\pi \lambda \beta }$ is the de Broglie thermal wavelength, with $\lambda ={\hslash }^2/2m$.

For an ideal gas, where the bodies are non interacting, $\mathcal{V}=0$, we immediately find $Z_N=V^N/{\mathsf{\Lambda }}^{-3N}N!$ where $\mathsf{\Lambda }=\sqrt{4\pi \lambda \beta }$ is the de Broglie thermal wavelength, with $\lambda ={\hslash }^2/2m$. Then

$$\begin{array}{c}\hfill {\langle R\rangle }_g=\kappa {\displaystyle \frac{1}{V}}{\displaystyle \frac{\partial }{\partial \beta }}ln\left[{\displaystyle \frac{1}{N!}}{\left({\displaystyle \frac{V}{{\mathsf{\Lambda }}^3}}\right)}^N\right]=-3n\kappa {\displaystyle \frac{\partial }{\partial \beta }}ln\mathsf{\Lambda }=-{\displaystyle \frac{3n}{2\beta }}\kappa ,\end{array}$$

(77)

where the functional derivative has been replaced by a partial derivative and $V/{\mathsf{\Lambda }}^3$ is the single particle translational partition function, familiar from elementary statistical mechanics.

In Ref. [4] we defined a virial inverse temperature $\tilde{\beta }$ stemming from the thermal fluctuations of the metric tensor, as 8

$$\begin{array}{c}\hfill {\tilde{\beta }}^{-1}\left(\mathit{x}\right)={\displaystyle \frac{\tilde{v}}{4}}{\langle T_{\mu }{}^{\mu }\rangle }_t,\end{array}$$

(78)

where $\tilde{v}$ is a positive constant. Therefore we find the following equivalence

$$\begin{array}{c}\hfill {\tilde{\beta }}^{-1}\left(\mathit{x}\right)={\displaystyle \frac{3n\tilde{v}}{8}}{\beta }^{-1}\left(\mathit{x}\right),\end{array}$$

(79)

or

$$\begin{array}{c}\hfill \tilde{T}\left(\mathit{x}\right)=T\left(\mathit{x}\right),\end{array}$$

(80)

with

$$\begin{array}{c}\hfill {\tilde{k}}_B={\displaystyle \frac{3n\tilde{v}}{8}}{k}_B,\end{array}$$

(81)

where $n\tilde{v}$ is an intensive quantity, if fact $\tilde{v}$ is a local volume [4]. In Ref. [21] it was also shown that at thermodynamic equilibrium ${\beta }^{\mu }\left(\mathit{x}\right)$ must be a killing vector of the manifold so it must be a constant four vector. Then T should be independent of $\mathit{x}$ at thermodynamic equilibrium.

Thermal equivalence principle

This equivalence proves that Einstein field equations offer a symmetric way to study statistical physics where one can either work at the level of the many body system encoded in the stress energy tensor on the right hand side of (72) or at the level of the thermal fluctuations of the metric tensor on the left hand side of (72). The two descriptions are equivalent. This is a thermal equivalence principle in physics: “given a many body system in general relativity its thermal equilibrium properties derived from its statistical physics description are equivalent to the properties of a statistical physics description of the metric of the spacetime that it influences and viceversa”. In other words: “A statistical ensemble of bodies goes into thermal equilibrium with the spacetime it occupies”. This is established by Equation (80).

5. Conclusions

In this work we determined a thermal equivalence principle for a statistical theory of gravitation. First of all it is important to realize that at low temperatures a statistical theory of gravity will necessarily put together the quantum world with our Universe ruled by general relativity. Outer space, or simply space, is the expanse that exists beyond Earth’s atmosphere and between celestial bodies. It contains very low particle densities, constituting a near perfect vacuum of predominantly hydrogen and helium plasma, permeated by electromagnetic radiation, cosmic rays, neutrinos, magnetic fields and dust. The baseline temperature of outer space, as set by the background radiation from the Big Bang, is $\approx 2.7$ K. Intergalactic space takes up most of the volume of the universe, but even galaxies and star systems consist almost entirely of empty space. Most of the remaining mass-energy in the observable universe is made up of an unknown form, dubbed dark matter (60% of the Universe) and dark energy (27% of the Universe) [30].

The program of constructing a well defined statistical theory of our Universe is one of the greatest challenges of contemporary physics which had been foreseen by Einstein in his renown iconic phrase “God doesn’t play dice”. From the point of view of the challenge that it offers to mathematics one needs a way to create a bridge between the variational theory of functional integrals or more specifically path integrals and differential geometry or more specifically Riemannian geometry. From this point of view it seems natural to predict that differential topology will play a crucial role. Recently we carried out some path integral (Monte Carlo) simulations for Jellium (an electron plasma at low temperature) on the surface of a sphere, probably the simplest of all curved smooth manifolds. And already in that study we found important topological effects on the electrons paths. It is important to realize that this kind of calculations can be considered as toy simulations for a many body system on a more complex smooth manifold as needed by spacetime in general relativity.

In this work we show that the temperature as defined in this kind of statistical physics studies of many bodies plays the same role as the one that can be defined by a path integral on the spacetime metric that was introduced in Ref. [4] 9. This is perfectly natural from the point of view of the linear constraint given by the Einstein field equations. This symmetry between the statistical physics of many body matter in the Universe and a statistical physics theory of the metric tensor where matter lives, offers naturally a thermal equivalence principle stating that the material and the spacetime are in thermal equilibrium one another.