What is the correct relativistic distribution function?

Question

General Statement and Questions

I am trying to figure out the proper way to model a velocity/momentum distribution function that is correct in the relativistic limit. I would like to determine/know two things:

1. Is there an analytical form for an anisotropic relativistic momentum distribution (i.e., the relativistic analog of the bi-Maxwellian distribution)?
2. What does temperature (i.e., kinetic temperature) mean in the relativistic limit?
   • Temperature cannot be a Lorentz invariant, can it?
   • It certainly cannot be invariant if the average particle thermal energies correspond to relativistic thermal velocities, correct?
     • So how can just a simple scalar temperature be a proper normalization factor as a Lagrange multiplier in the Maxwell-Jüttner distribution, for instance?
3. Extra: Is there an appropriate relativistic version of the $\kappa$ distribution (see this arXiv PDF for reference, or e-print number 1003.3532 if you do not trust links)?

Background

I am aware of the Maxwell-Jüttner distribution of particle species $\alpha$, given by: $$ f_{\alpha} \left( p \right) = \Lambda \ exp \left[ - \Theta_{0} \ m c^{2} \ \gamma \left( p \right) - \sum_{i=1}^{3} \Theta_{i} \ c \ p^{i} \right] $$ where $\Lambda$, $\Theta_{0}$, and $\Theta_{i}$ are Lagrange multipliers, $p$ is the relativistic momentum, and $\gamma\left( p \right)$ is the Lorentz factor. The $\Theta_{\nu}$ terms are 4-vector components with units of inverse energy.

In the isotropic limit, one can set all $\Theta_{i}$ $\rightarrow$ 0. This leads to the canonical form of the isotropic relativistic momentum distribution function is given by: $$ f_{\alpha} \left[ \gamma \left( p \right) \right] = \Lambda e^{ - \Theta_{0} \ m c^{2} \ \gamma \left( p \right) } $$ where $\Theta_{0}$ was shown¹ to be the inverse of a temperature.

The Issue

The definition of $\Lambda$, however, has led to multiple results, as stated by Treumann et al. [2011]:

the correct (non-angular-dependent part of the) relativistic thermal-equilibrium distribution should become the modified-Jüttner distribution. (The ordinary Maxwell-Jüttner distribution function was derived by F. Jüttner, 1911, who obtained it imposing translational invariance in momentum space only.)

An attempt⁴ was made to derive $\Lambda$ by imposing Lorentz invariance on only momentum space, ignoring the spatial coordinates of the volume integral. However, Treumann et al. [2011] note that:

This is either not justified at all or it is argued that the particles are all confined to a fixed box which is unaffected by the Lorentz transformation and invariance. However, the momentum and configuration space volume elements the product of which forms the phase-space volume element, are not independent, as we have demonstrated above. Even in this case of a fixed outer box, the particle's proper spaces experience linear Lorentz contractions when seen from the stationary frame of the observer, i.e., from the box-frame perspective. The consequence is that the extra proper Lorentz factor $\gamma\left( p \right)$ in the phase-space volume element cancels thereby guaranteeing and restoring Lorentz invariance...

They go on to show that the correct Lagrange multipliers are: $$ \Theta_{0} = \frac{ 1 }{ T } \\ \Lambda = \frac{ N^{0} }{ 4 \pi \ m^{2} T^{2} } \left[ 3 K_{2}\left( \frac{ m c^{2} }{ T } \right) + \frac{ m c^{2} }{ T } K_{1}\left( \frac{ m c^{2} }{ T } \right) \right]^{-1} $$ where $N^{0}$ is the scalar part of the particle current density 4-vector (i.e., number density), $K_{i}(x)$ is the second order modified Bessel function, and $T$ is a scalar temperature. Notice there is an additional term (i.e., $K_{1}(x)$) in the normalization factor $\Lambda$, which is why they called this the modified Maxwell-Jüttner distribution. This accounts for Lorentz invariance in the phase-space element, not just momentum-space.

What I am looking for...

Regardless of its accuracy, the distribution function in Treumann et al. [2011] still only assumes an isotropic distribution and I am still a bit confused how the temperature is just a scalar. In plasma physics, it is more appropriate to think of it as kind of a pseudotensor derived from the pressure tensor or 2nd moment of the distribution function. So am I supposed to interpret relativistic temperatures through the energy-momentum tensor or something else? See more details about velocity moments here: https://physics.stackexchange.com/a/218643/59023.

In many situations, plasmas can be described as either a bi-Maxwellian or bi-kappa [e.g., Livadiotis, 2015] velocity distribution functions. The bi-Maxwellian is given by: $$ f\left( v_{\parallel}, v_{\perp} \right) = \frac{ 1 }{ \pi^{3/2} \ V_{T \parallel} \ V_{T \perp}^{2} } \ exp\left[ - \left( \frac{ v_{\parallel} - v_{o, \parallel} }{ V_{T \parallel} } \right)^{2} - \left( \frac{ v_{\perp} - v_{o, \perp} }{ V_{T \perp} } \right)^{2} \right] $$ where $\parallel$($\perp$) refer to directions parallel(perpendicular) with respect to a quasi-static magnetic field, $\mathbf{B}_{o}$, $V_{T_{j}}$ is the $j^{th}$ thermal speed (actually the most probable speed), and $v_{o, j}$ is the $j^{th}$ component of the bulk drift velocity of the distribution (i.e., from the 1st velocity moment).

The bi-kappa distribution function is given by: $$ f\left( v_{\parallel}, v_{\perp} \right) = A \left[ 1 + \left( \frac{ v_{\parallel} - v_{o, \parallel} }{ \sqrt{ \kappa - 3/2 } \ \theta_{\parallel} } \right)^{2} + \left( \frac{ v_{\perp} - v_{o, \perp} }{ \sqrt{ \kappa - 3/2 } \ \theta_{\perp} } \right)^{2} \right]^{- \left( \kappa + 1 \right) } $$ where the amplitude is given by: $$ A = \left( \frac{ \Gamma\left( \kappa + 1 \right) }{ \left( \pi \left( \kappa - 3/2 \right) \right)^{3/2} \ \theta_{\parallel} \ \theta_{\perp}^{2} \ \Gamma\left( \kappa - 1/2 \right) } \right) $$ and where $\theta_{j}$ is the $j^{th}$ thermal speed (also the most probable speed), $\Gamma(x)$ is the complete gamma function and we can show that the average temperature is just given by: $$ T = \frac{ 1 }{ 3 } \left( T_{\parallel} + 2 \ T_{\perp} \right) $$ if we assume a gyrotropic distribution (i.e., shows symmetry about $\mathbf{B}_{o}$ so that the two perpendicular components of a diagonalized pressure tensor are equal).

In summary, I would prefer a relativistically consistent bi-kappa distribution but would be very happy with the bi-Maxwellian version as well.

Update

After several conversations with the R. Treumann, he and his colleague decided to look into an anisotropic Maxwell-Jüttner distribution. I also referred him to this page and he decided to try and remain consistent with the original Maxwell-Jüttner distribution normalization to avoid further confusion.

His new results can be found in the arXiv paper with e-print number 1512.04015.

Summary of Results
One of the interesting things noted by Treumann and Baumjohann is that one cannot simply take the expression for energy and split the momentum into parallel and perpendicular terms as has occasionally been done in the past. Part of the issue is that the normalizations factors, i.e., temperature-like quantities, are not relativistically invariant. The temperature in this case is more akin to a pseudotensor than a scalar (Note: I use pseudotensor very lightly/carelessly here).

They use the Dirac tensor, from the Klein-Gordon approach, to define the energies. They treat the pressure as a proper tensor, with an assumed inverse, to define what they call the temperature tensor.

Unfortunately, the equation cannot be reduced analytically, but it is useful none-the-less given the alternative is to assume the unrealistic case of an isotropic velocity distribution in a relativistic plasma.

References

1. Israel, W. "Relativistic kinetic theory of a simple gas," J. Math. Phys. 4, 1163-1181, doi:10.1063/1.1704047, 1963.
2. Treumann, R.A., R. Nakamura, and W. Baumjohann "Relativistic transformation of phase-space distributions," Ann. Geophys. 29, 1259-1265, doi:10.5194/angeo-29-1259-2011, 2011.
3. Jüttner, F. "Das Maxwellsche Gesetz der Geschwindigkeitsverteilung in der Relativtheorie," Ann. Phys. 339, 856-882, doi:10.1002/andp.19113390503, 1911.
4. Dunkel, J., P. Talkner, and P. Hänggi "Relative entropy, Haar measures and relativistic canonical velocity distributions," New J. Phys. 9, 144-157, doi:10.1088/1367-2630/9/5/144, 2007.
5. Livadiotis, G. "Introduction to special section on Origins and Properties of Kappa Distributions: Statistical Background and Properties of Kappa Distributions in Space Plasmas," J. Geophys. Res. Space Physics 120, 1607-1619, doi:10.1002/2014JA020825, 2015.

Answer 1

I think you make it sound much more mysterious than it is. The relativistic distribution function is $$ f_p = \frac{1}{(2\pi)^3}\exp(-(\mu-u\cdot p)/T)\, $$ where $u_\alpha$ is the 4-velocity of the fluid, $p_\alpha$ is the 4-momentum of the particle, $T$ is temperature, and $\mu$ is the chemical potential. This is sometimes called the Juttner distribution, and it has obvious generalizations to Bose-Einstein and Fermi-Dirac statistics.

This expression is obviously Lorentz-invariant ($u\cdot p$ is a scalar, and so are $\mu$ and $T$), it reduces to Boltzmann in the rest frame of the fluid, and it is Gallilean invariant for slowly moving fluids.

The funny Bessel functions appear if you try to determine the fugacity $e^{\mu/T}$ in terms of the density, $n=N_0$ with $$ N_\mu = \int \frac{d^3p}{p^0} p_\mu f_p\, , $$ because now your normalization factor contains integrals over $\exp(-\sqrt{p^2+m^2}/T)$

Additional Remarks: After some prodding by the OP I looked at the Treumann et al. paper (it is available on the arxive, http://arxiv.org/abs/1105.2120). I initially thought that this was just a needlessly complicated rederivation of well-known results, but this is not the case. The paper is just wrong. (Frankly, it is not a good sign if a paper that points out a major flaw in relativistic kinetic theory is published in the physics section of the arxive, and a geophysics journal.)

The paper starts out o.k., observing that since $d^3xd^3p$ is Lorentz invariant, $f_p$ must be a scalar. However, they write down a distribution function which is clearly not a scalar unless I assume that $T$ is the zeroth component of a vector. Even if that could be arranged, the result is nor right because it is not Galilean invariant for small velocities. He then writes down a non-covariant expression for $N_\mu$. If $f_p$ is a scalar, $d^3p p_\mu f_p$ is not a vector, and as a result his particle density $N_0$ does not transform as a density (I wrote the correct expression above). Since $N_0$ is wrong, the normalization $\Lambda$ is also wrong.

More Remarks: After more prodding, an attempt to answer the original questions:

1) The original Boltzmann distribution is anisotropic if the fluid velocity $\vec{u}$ is not zero (the distribution is isotropic in the fluid rest frame), but presumably you are looking for something more general. The distribution function of a viscous fluid is anisotropic already in the rest frame of the fluid, $f_p=f_p^0 +\delta f_p$ with $$\delta f_p \sim \eta\sigma_{\alpha\beta}p^\alpha p^\beta,$$ where $\sigma_{\alpha\beta}$ is the relativistic strain tensor. Even more generally, we can parameterize the distribution function in terms of the currents and the stresses. This is the relativistic version of Grad's 13 moment method, see, for example, http://arxiv.org/abs/1301.2912. People have also written down models for highly anisotropic distribution functions, see for example http://arxiv.org/abs/1007.0130 .

2) The temperature is a scalar. It is related to the kinetic energy per particle in the rest frame. The kinetic energy density is the 00 component of a tensor, and transforms accordingly.

Answer 2

Disclaimer

The following answer is largely taken from the arXiv paper with e-print number 1512.04015 by Treumann and Baumjohann (from here on I abbreviate references to this paper as TB15).

Background

It is well known that the Maxwell-Jüttner distribution works well for a momentum/energy distribution with an isotropic, scalar temperature, $T$. This distribution is generally applied to hot, relativistic plasmas. Herein lies the problem and the source of the question. Hot plasmas, even in non-relativistic cases, are rarely isotropic.

As an aside and to clarify definitions, TB15 state:

The meaning of classical here is not that quantum effects are excluded which, in a hot plasma is anyway the case, because of the tinniness of the thermal quantum length $\lambda_{q} = \sqrt{2 \pi \hbar^{2} / m_{e} T}$. It implies that particle numbers are conserved. This inhibits pair creation and annihilation and thus restricts to temperatures below, roughly, $T_{e} < 2 m_{e} c^{2} \approx 1 \ MeV$.

For non-relativistic plasmas, dealing with a temperature anisotropy is easy because the kinetic energy is additive, as shown in the example bi-Maxwellian and bi-kappa equations in the above question. In a relativistic plasma, the energy has an additional momentum-dependent factor, $\gamma\left( \mathbf{p} \right)$, which is the Lorentz factor, defined as: $$ \gamma\left( \mathbf{p} \right) = \sqrt{ 1 + \left( \frac{ \mathbf{p} }{ m \ c } \right)^{2} } \tag{1} $$ where $\mathbf{p}$ is the momentum, $m$ is the particle rest mass, and $c$ is the speed of light in vacuum. Then the particle's total energy can be defined as: $$ \epsilon\left( \mathbf{p} \right) = \gamma\left( \mathbf{p} \right) \ m \ c^{2} \tag{2} $$ The problem is that $\epsilon\left( \mathbf{p} \right)$ is no longer additive. It is often the case that anisotropic temperatures are introduced in the following form: $$ \epsilon\left( \mathbf{p} \right) = m \ c^{2} \sqrt{ 1 + \left( \alpha_{\perp} \bar{p}_{\perp} \right)^{2} + \left( \alpha_{\parallel} \bar{p}_{\parallel} \right)^{2} } \tag{3} $$ where $\bar{p}_{\perp (\parallel)} = c \ p_{\perp (\parallel)}/T_{\perp (\parallel)}$ and $\alpha_{\perp (\parallel)} = T_{\perp (\parallel)}/mc^{2}$. There are several issues produced by this format, as TB15 state:

This inhibits a simple definition of a distribution expressed solely in $\gamma$ enforcing the use of the momentum distribution. Moreover, in the exponent of the relativistic distribution the normalisation to temperature becomes arbitrary... A further aspect concerns the question, which $\gamma$ is meant. If the observer is dealing from his inertial system with a moving extended gaseous or plasma body with relativistic bulk momentum $\mathbf{P}$ the problem is different from the one where the observer is embedded at relative rest into the volume with all particles moving relative to him at relativistic momenta $\mathbf{p}$. In the latter case the relativistic effects are intrinsic to the system which is in contrast to the above simpler extrinsic bulk effect...

Answer

They use the Klein-Gordon approach because particle spins are ignored. Thus, squaring Equation 2 results in: $$ \begin{align} \epsilon^{2} & = \left[ c \ \mathbf{p} + i \ m \ c^{2} \ \mathbf{e} \right] \cdot \left[ c \ \mathbf{p} - i \ m \ c^{2} \ \mathbf{e} \right] \tag{4a} \\ & = \left( c \ p_{\nu} + i \ m \ c^{2} \ e_{\nu} \right) \delta_{\mu}^{\nu} \left( c \ p^{\mu} - i \ m \ c^{2} \ e^{\mu} \right) \tag{4b} \end{align} $$ where $\mathbf{p} = \left( p_{\perp} \cos{\phi}, p_{\perp} \sin{\phi}, p_{\parallel} \right)$, $\delta_{\mu}^{\nu}$ is the Dirac tensor, and $\mathbf{e}$ is a unit vector that can be arbitrarily chosen. TB15 make several points about the form of Equation 4b:

Of course, the energy maintains its scalar property which is reflected in the scalar dot-product of the two bracketed expressions. The appearance of the imaginary element has no effect on the energy. In particular, it does not imply any damping of particle motion! However, by the above factorisation, the Hamiltonian can be split into two Hamiltonians which are linear in the vector $\mathbf{p}$. This is important to realise. It permits the consistent introduction of a pressure anisotropy...

They then construct a pressure tensor assuming a fixed scalar number density, $N$, and write it in terms of the temperature as: $$ \begin{align} \mathbb{P} & = P_{\perp} \delta_{\mu}^{\nu} + \left( P_{\parallel} - P_{\perp} \right) \delta_{3}^{3} \tag{5a} \\ & = N \left[ T_{\perp} \delta_{\mu}^{\nu} + \left( T_{\parallel} - T_{\perp} \right) \delta_{3}^{3} \right] \tag{5b} \end{align} $$ where direction $3$ is assumed to be along the quasi-static magnetic field direction, $\mathbf{b} = \mathbf{B}/B$ and Equation 5b applies for an ideal gas. Here is where things become difficult because the temperature is not merely a scalar quantity anymore, thus we cannot simply divide by it as was done in Equation 3 above. If the pressure tensor is invertible, which it should be, then we have: $$ \begin{align} \mathbb{P}^{-1} & = N^{-1} \left[ T_{\perp}^{-1} \delta_{\nu}^{\mu} + \left( T_{\parallel}^{-1} - T_{\perp}^{-1} \right) \delta_{3}^{3} \right] \tag{6a} \\ & = N^{-1} \Theta \tag{6b} \end{align} $$ where $\Theta$ is the "inverse temperature tensor." Now we can rewrite Equation 4b as: $$ \begin{align} \epsilon^{2} & = \left[ c \ p_{\nu} + i \ m \ c^{2} \ e_{\nu} \right] \Theta_{\mu}^{\nu} \Theta_{\lambda}^{\mu} \left[ c \ p^{\lambda} - i \ m \ c^{2} \ e^{\lambda} \right] \tag{7a} \\ & = \frac{ m^{2} c^{4} }{ T_{\perp}^{2} } \left( \frac{ p^{2} }{ m^{2} c^{2} } + 1 \right) \left( \sin^{2}{\theta} + A^{2} \cos^{2}{\theta} \right) \tag{7b} \end{align} $$ where $A = T_{\perp}/T_{\parallel}$ and $\mathbf{p} = \left( p \sin{\theta}, p \cos{\theta} \right)$ (i.e., just reduce the momentum to parallel and perpendicular). If we define $\beta_{\perp} = m \ c^{2}/T_{\perp}$ and $\psi = A^{2} - 1$ (i.e., isotropy implies $\psi = 0$), then we can write the anisotropic, relativistic momentum distribution for an ideal gas as: $$ \begin{align} f\left( \mathbf{p} \right) & = C_{o} exp\left( - \beta_{\perp} \sqrt{ \left( 1 + \frac{ p^{2} }{ m^{2} c^{2} } \right) \left( 1 + \psi \cos^{2}{\theta} \right) } \right) \tag{8a} \\ & = C_{o} exp\left( - \beta_{\perp} \ \gamma\left( p \right) \sqrt{ \left( 1 + \psi \cos^{2}{\theta} \right) } \right) \tag{8b} \end{align} $$

Normalization

The normalization is with respect to the scalar density, $N$, and the approach is given in TB15. They start by defining $\beta = \beta_{\perp} \sqrt{ 1 + \psi \cos^{2}{\theta} }$ and then show: $$ N = 4 \ \pi \ C_{o} \ \left( m \ c \right)^{3} \int_{0}^{1} \ dx \ \frac{ K_{2}\left( \beta\left( x \right) \right) }{ \beta\left( x \right) } \tag{9} $$ where $K_{2}\left( z \right)$ is the second order modified Bessel function of argument $z$ and the integral does not have an analytical solution. One can rewrite Equation 9 in the following form for $C_{o}$ given by: $$ C_{o} = \frac{ N \ \beta_{\perp} }{ 2 \ \pi \ \left( m \ c \right)^{3} } \left[ \int_{1}^{1 + \psi} \ dz \ \frac{ K_{2}\left( \beta_{\perp} \sqrt{z} \right) }{ \sqrt{z \left( z - 1 \right)} } \right]^{-1} \tag{10} $$ which also cannot be analytically solved.

In the limit of small anisotropy (i.e., $\psi < 1$), one can approximate the integral in Equation 9 as: $$ \int_{0}^{1} \ dx \ \frac{ K_{2}\left( \beta\left( x \right) \right) }{ \beta\left( x \right) } \approx \frac{ 1 }{ 2 \ \beta_{\perp} \sqrt{ \psi } } \left\{ \left[ K_{1}\left( \frac{ \beta_{\perp} }{ 2 } \right) \right]^{2} - \left[ K_{1}\left( \frac{ 1 }{ 2 } \beta_{\perp} \sqrt{ 1 + \psi } \right) \right]^{2} \right\} \tag{11} $$

All of the above assume a perpendicular anisotropy, but analogous forms can be found for parallel anisotropies. TB15 also comments on the generalization to include an external potential fields, but I leave that to the article.

Important Note

The above expressions only apply to bosons. To include fermions, one would need to include spin in the calculations which would greatly increase the complexity of the problem. As stated by TB15:

We note that a similar calculation can also be performed for Fermionic plasmas. Then, however, one must refer to Dirac’s method of splitting the relativistic energy. This leads to the use of the Dirac matrices and substantially more involved expressions including the spin of the particles. Fermionic effects of the anisotropy will be expected only at very low temperatures and in very strong magnetic fields. Such situations might occur in the Quantum-Hall effect and when dealing with the interiors of highly magnetised objects like pulsars and magnetars... An interesting conlusion which can be drawn is that in relativistic media the temperature and its inverse, usually called $\beta = 1/T$ must be understood as vectors. This is the consequence of the pressure tensor.

As a final conclusion, they note:

One might, moreover, think of that the split of the energy would provide two different versions of the distribution. This is, however, not the case. The linearity in the momentum of the two energy-vectors has no classical meaning. Energies are scalars and do not transform like vectors. Multiplication with the inverse relativistic temperature vector does not release this property. The linearity would lead to a linear dependence of the exponent in the distribution which is not Gaussian anymore and thus loses the property of a probability. Moreover, each of the energy vectors becomes complex. This complexity resolves in quantum mechanics in the operator formalism of the Klein-Gordon and Dirac equations but has no meaning in classical physics.

Updates

Since the posting of this question and answer there have been two refereed publications on the topic, one by Treumann and Baumjohann [2016] (doi:10.5194/angeo-34-737-2016) and one by Livadiotis [2016] (doi:10.5194/angeo-34-1145-2016) (both papers are Open Access so there should be no paywall). Treumann and Baumjohann [2016] correct their arXiv result and find the same normalization factor as the original paper by Jüttner [1911].

Livadiotis [2016] is more rigorous in its treatment giving the probability distribution: $$ P\left( \mathbf{p}, \Theta \right) = C_{o} \ e^{ - \Theta^{-1} \ \sqrt{ 1 + \sum_{i=1}^{\delta} \tfrac{\Theta}{\Theta_{i}} \left( \tfrac{p_{i}}{m \ c} \right)^{2} } } $$ where $\delta$ is the number of degrees of freedom, $\Theta_{i} = \tfrac{k_{B} \ T_{i}}{m \ c^{2}}$ is the i^th component of the normalized temperature, $k_{B}$ is the Boltzmann constant, and $\Theta = \tfrac{1}{\delta} \ \sum_{i=1}^{\delta} \Theta_{i}$, and $C_{o}$ is given by: $$ C_{o} = \pi^{\left( \tfrac{ 1 - \delta }{ 2 } \right)} \ 2^{-\left( \frac{ \delta + 1 }{ 2 } \right)} \ \left( m \ c \right)^{-\delta} \ \frac{ \sqrt{ \Theta } }{ K_{\tfrac{ \delta + 1 }{ 2 }}\left( \tfrac{ 1 }{ \Theta } \right) } \ \left[ \prod_{i=1}^{\delta} \ \Theta_{i}^{-1/2} \right] $$ where $K_{n}\left( x \right)$ is the Macdonald function (it is also a modified Bessel function), $m$ is the particle mass, and $c$ is the speed of light.

In the specific case of $\delta = 3$ and $\left\{ T_{i} \right\}_{i=1}^{3} = \left( T_{\parallel}, T_{\perp} \right)$, $C_{o}$ will reduce to: $$ C_{o} \rightarrow \frac{ \sqrt{ \Theta \ A } \ \beta_{\perp} }{ 4 \pi \ \left( m \ c \right)^{3} \ K_{2}\left( \tfrac{ 1 }{ \Theta } \right) } $$ where $A = T_{\perp}/T_{\parallel}$ and $\beta_{j} = \tfrac{ m \ c^{2} }{ k_{B} \ T_{j} }$.

Answer 3

E. Lehmann (J. of Math. Physics 47 023303, 2006) attempted this in his paper Covariant equilibrium statistical mechanics (the paper is also in arXiv). He found that his fully covariant distribution became flatter than the non-covariant one. Also that while the non-covariant partition function continued to increase with temperature the fully covariant one reached a maximum at about logT~12.

Even though there is loss of particle count invariance in the fully covariant case one can still define an energy density and with that a temperature in the classical sense assuming thermodynamic equilibrium. However as the temperature increases we will reach conditions that prevailed in the universe at the time of the big bang and I doubt if anyone will claim that can be an equilibrium situation. Physical continuity would suggest that the effects would come gradually. The temperature at Lehmann's maximum occurs in core collapse supernovae.

It is possible to reproduce the principal results of Lehmann's analysis (maximum temperature and flattening distribution) using the simple argument Maxwell first used in the classical case when he was still a teenager. That is we split the 3D distribution function into a product of three identical functions and use the composition of velocities to get the form of the distribution. The classical composition is easily replaced by the relativistic using the Lorentz transformation. The resulting distribution becomes flat at logT~11 after which it becomes U-shaped. There are actually two underlying assumptions in Maxwell's approach: molecular chaos and that an equilibrium distribution must be isotropic. Isotropy might be a necessary condition but this result suggests that it is not sufficient.