# Neutrino implementation¶

SWIFT can also accurately model the effects of massive neutrinos in cosmological simulations. At the background level, massive neutrinos and other relativistic species can be included by specifying their number and masses in the cosmology section of the parameter file (see Cosmology).

At the perturbation level, neutrinos can be included as a separate particle
species (`PartType6`

). To facilitate this, SWIFT implements the
\(\delta f\) method for shot noise suppression (Elbers et al. 2020). The method
works by statistically weighting the particles during the simulation,
with weights computed from the Liouville equation using current and
initial momenta. The method can be activated by specifying
`Neutrino:use_delta_f`

in the parameter file.

The implementation of the \(\delta f\) method in SWIFT assumes a
specific method for generating the initial neutrino momenta (see below).
This makes it possible to reproduce the initial momentum when it is
needed without increasing the memory footprint of the neutrino particles.
If perturbed initial conditions are not needed, the initial momenta can
be generated internally by specifying `Neutrino:generate_ics`

in the
parameter file. This will assign `PartType6`

particles to each
neutrino mass specified in the cosmology and generate new velocities
based on the homogeneous (unperturbed) Fermi-Dirac distribution. In
this case, placeholder neutrino particles should be provided in the
initial conditions with arbitrary masses and velocities, distributed
uniformly in the box.

## Relativistic Drift¶

At high redshift, neutrino particles move faster than the speed of light if the usual Newtonian expressions are used. To rectify this, SWIFT implements a relativistic drift correction. In this convention, the internal velocity variable (see theory/Cosmology) is \(v^i=a^2u^i=a^2\dot{x}^i\gamma^{-1}\), where \(u^i\) is the spatial part of the 4-velocity, \(a\) the scale factor, and \(x^i\) a comoving position vector. The conversion factor to the coordinate 3-velocity is \(\gamma=ac/\sqrt{a^2c^2+v^2}\). This factor is applied to the neutrino particles throughout the simulation.

## Generating Fermi-Dirac momenta¶

The implementation of the \(\delta f\) method in SWIFT assumes that
neutrinos were initially assigned a Fermi-Dirac momentum using the following
method. Each particle has a fixed 64-bit unsigned integer \(\ell\) given
by the particle ID [1] (plus an optional seed: `Neutrino:neutrino_seed`

).
This number is transformed into a floating point number \(u\in(0,1)\),
using the following pseudo-code based on splitmix64:

```
m = l + 0x9E3779B97f4A7C15
m = (m ^ (m >> 30)) * 0xBF58476D1CE4E5B9;
m = (m ^ (m >> 27)) * 0x94D049BB133111EB;
m = m ^ (m >> 31);
u = (m + 0.5) / (UINT64_MAX + 1)
```

This is subsequently transformed into a Fermi-Dirac momentum \(q = F^{-1}(u)\) by evaluating the quantile function. To generate neutrino particle initial conditions with perturbations, one first generates momenta from the unperturbed Fermi-Dirac distribution using the above method and then applies perturbations in any suitable manner.

When using the \(\delta f\) method, SWIFT also assumes that `PartType6`

particles are assigned to all \(N_\nu\) massive species present in the
cosmology, such that the particle with fixed integer \(\ell\) corresponds
to species \(i = \ell\; \% \;N_\nu\in[0,N_\nu-1]\).

The sampled Fermi-Dirac speeds and neutrino masses are written into the
snapshot files as `SampledSpeeds`

and `MicroscopicMasses`

.

[1] | Currently, it is not guaranteed that a particle ID is unique. |