Precision calculation of $N_{\text{eff}}$ with Neutrino Direct Simulation Monte Carlo

TL;DR

The paper introduces $ u$DSMC, a Direct Simulation Monte Carlo method tailored to the expanding, thermal plasma of the early Universe, to solve neutrino Boltzmann dynamics with finite electron mass, three-flavor oscillations, and finite-temperature QED corrections. By representing EM energy with a temperature field and neutrinos as simulation particles, the approach avoids momentum-space discretization and preserves energy conservation in collisions, enabling precise neutrino decoupling studies within $ Lambda$CDM. The authors validate $ u$DSMC against established integrated and unintegrated Boltzmann solvers and demonstrate a high-precision calculation of $N_{ ext{eff}} = 3.0439 \,\pm\, 0.0006$, consistent with state-of-the-art results and providing a robust platform for exploring non-standard scenarios with non-thermal injections. The work highlights $ u$DSMC as an independent cross-check tool and a flexible framework for future explorations of beyond-$ Lambda$CDM physics in the neutrino sector, including higher-order QED corrections and exotic injection channels.

Abstract

Neutrino Direct Simulation Monte Carlo ($ν$DSMC) is a Monte Carlo method for solving the neutrino Boltzmann equation in the early Universe, designed to track the evolution of cosmic neutrinos across a wide range of cosmological scenarios. We develop a complete $ν$DSMC solver that consistently incorporates the effects of the electron mass, three-flavour neutrino oscillations, and finite-temperature QED corrections to the thermodynamics of the electromagnetic plasma. As a first application, we perform a high-precision calculation of neutrino decoupling in the standard cosmological model and obtain $N_{\text{eff}} = 3.0439 \pm 0.0006$, in excellent agreement with state-of-the-art results.

Figures (2)

• Figure 1: The modification of the No-Time-Counter scheme used in $\nu$DSMC to simulate interactions of neutrinos and EM particles inside cells. First, we perform neutrino oscillations based on the stationary probabilities \ref{['eq:oscillation-probabilities']}. Then, we determine $N_{\text{pairs}}$ pairs of particles to interact -- the sum of Eqs. \ref{['eq:Npairs-category-1-optimized']}-\ref{['eq:Npairs-category-3-optimized']}. Then, based on the same equations, iteratively, we decide what type of pair to sample -- $\nu\nu$, $\nu e$, or $ee$. For each pair, we compute its interaction weight $(\sigma v)_{\text{pair}}$ and make an intermediate decision on whether it will interact using the criterion \ref{['eq:interaction-acceptance-optimized']}. Then, we sample the kinematics of the interacting particles, generate the final states resulting from the collision, and make the final decision of whether the interaction takes place from the Pauli principle \ref{['eq:pauli-blocking']}. Finally, we update the local properties of the plasma: the EM plasma temperature and the number of EM particles, as well as neutrino flavor distributions by the oscillation probabilities, Eq. \ref{['eq:oscillation-probabilities']}.
• Figure 2: The evolution of the Universe in the temperature domain $T \in (30\text{ keV},10\text{ MeV})$ in $\Lambda$CDM as obtained within the $\nu$DSMC framework. Top panel: the temperature evolution of the ratio $8/7 \cdot (11/4)^{\frac{4}{3}}\rho_{\nu}/\rho_{\text{EM}}$, which tends to $N_{\text{eff}}$ in the limit of small temperatures $T\ll m_{e}$. Bottom panel: the ratio of the electron neutrino distribution function to the Fermi-Dirac distribution in terms of the comoving momentum $y = p \cdot (a(T)/a(T = 10\text{ MeV}))$. The Fermi-Dirac function is evaluated at temperatures $T_{\text{eff}}$ being the solution $E_{\nu_{e},\text{total}}/V_{\text{system}} = 3/2\zeta(3)/\pi^{2} T_{\text{eff}}^{3}$. The wiggles in the domain $y \lesssim 10$ are caused by limited statistics in this region (corresponding to $p/T \lesssim 0.02$).