Radiation and Neutrino Transport

Radiation and neutrino transport couple a distribution in position, direction, energy, and species to relativistic matter. Full Boltzmann evolution is usually too expensive for merger-scale calculations, so leakage, moment closures, Monte Carlo packets, or hybrids make different approximations. No method is uniformly best across opaque diffusion, transparent free streaming, intersecting beams, and rapidly moving matter. The finite claim must name regime, interactions, spectral treatment, closure, resolution, and observable.

Helpful background. Formulations of General Relativity supplies spacetime coupling, Ordinary Differential Equations supplies stiff source integration, Relativistic Matter and Compact Objects supplies microphysics, Computational Relativity supplies evidence standards, and Fluid and Kinetic Theory supplies kinetic moments and collision terms.

Transport representations

The distribution function obeys a curved-spacetime Boltzmann equation. Moment schemes evolve radiation energy and momentum and require a closure for higher moments. M1 analytic closure is deterministic and efficient, with correct diffusion and free-stream limits in simple geometries, but cannot faithfully represent crossing beams because it combines them into one local flux direction.

Monte Carlo methods sample phase-space packets and naturally represent angular and spectral structure, but introduce sampling noise and can be inefficient in optically thick matter. Leakage estimates local emission and diffusion times without transporting the full radiation field. Hybrids can use Monte Carlo information to close moments. Foucart 2018, §2 and §§3.1–3.4, Figures 4–11 develops one moment–Monte Carlo closure and explicitly notes coupling limitations.

Inputs include species, energy groups or gray approximation, opacities and reaction-rate versions, fluid-frame conventions, closure, implicit–explicit coupling, packet count and random seed, floors, boundaries, and gravity treatment. Energy and lepton-number exchange must be equal and opposite between matter and radiation up to declared discretization and sampling error.

Why two moments cannot determine the angular distribution

In a local inertial frame with $c=1$, let the energy-integrated intensity be $I(\mathbf n)$ for unit propagation direction $\mathbf n$. Choose its normalization so

$$E=\int I\,d\Omega,\qquad F^i=\int n^iI\,d\Omega,\qquad P^{ij}=\int n^in^jI\,d\Omega.$$

In flat, source-free transport, $\partial_tE+\partial_iF^i=0$ and $\partial_tF^j+\partial_iP^{ij}=0$. Evolving only $E$ and $F^i$ therefore requires an additional prescription for $P^{ij}$. Isotropic radiation has $F^i=0$ and $P^{ij}=E\delta^{ij}/3$. Two equal, opposite pencil beams along $x$ have the same $E$ and $F^i=0$, but $P^{xx}=E$ and $P^{yy}=P^{zz}=0$. No local algebraic closure depending only on $(E,F^i)$ can represent both states. Refining a spatial grid cannot recover angular information absent from the evolved variables. This elementary example is why a crossing-beam failure is useful evidence about closure scope rather than a reason to hide the test.

Comparison of methods and results

The evidence cutoff is 2026-08-18.

Evidence	Method comparison	Supported scope	Important limitation
2018	Analytic moment closure versus Monte Carlo-informed closure	Idealized relativistic transport tests	Stability and cost of full coupling remained open
2023 review	Leakage, moment, and Monte Carlo methods	Merger and remnant use cases	Important interactions and full fidelity remain difficult
2024 study	Improved two-moment versus Monte Carlo and reaction-rate variants	One binary with a short-lived remnant	Some global properties changed by 10–30%; not a population-wide bound

Foucart 2023, §§4.1–4.3 and §5.2 compares algorithms and missing interactions. Foucart et al. 2024, §II.F, §III, Tables II–IV and §IV studies method and rate sensitivity for one short-lived-remnant system; its 10–30 percent differences must not be generalized to every merger.

Competing positions and benchmark ladder

One credible position favors moments for production because they are deterministic, relatively affordable, and accurate in many optically thick or single-beam regimes. Another favors Monte Carlo or hybrid information where angular crossings, spectral effects, and complex transparency make analytic closure the dominant systematic. A third accepts leakage for exploratory bulk dynamics while excluding precision composition or luminosity claims. Which position is appropriate depends on the observable and regime; the question is not settled by one benchmark.

Make the crossing-beam negative control quantitative in flat spacetime with $c=1$. Start from vacuum on $(x,y)\in[-1,1]^2$. At $x=-1$ inject a beam along $+x$ with amplitude $I_x(y)=\exp[-(y/0.1)^2]$; at $y=-1$ inject a beam along $+y$ with $I_y(x)=\exp[-(x/0.1)^2]$. Set every other incoming angular component to zero, use outflow at $x=1$ and $y=1$, and set absorption, emission, and scattering opacities to zero. At $t=3$, after both fronts have crossed the domain, the independent kinetic reference is

$$E_{\rm ex}=I_x(y)+I_y(x), \qquad F^x_{\rm ex}=I_x(y), \qquad F^y_{\rm ex}=I_y(x),$$

with $P^{xx}_{\rm ex}=I_x$, $P^{yy}_{\rm ex}=I_y$, and $P^{xy}_{\rm ex}=0$. Run $N^2$ spatial grids with $N=64,128,256$ and CFL factors $0.4,0.2$. For discrete ordinates use $16,32,64$ uniformly spaced directions including the coordinate axes. For Monte Carlo use $10^5$, $4\times10^5$, and $1.6\times10^6$ packets with eight fixed, recorded seeds at each level. Run M1 on the same spatial and temporal sequence as the deliberate closure-negative control.

The pencil beams are angular delta distributions. In a discrete-ordinate implementation, normalize their populated ordinate intensities by the quadrature weights so the prescribed integrated beam energy remains $I_x$ or $I_y$ when the angular grid changes. Because both beam directions are included exactly, this test primarily probes streaming, crossing, closure, and spatial error; it does not establish accuracy for a general angular distribution. Add an oblique or finite-width beam for that purpose. In a packet comparison, report both per-cell noise and the uncertainty of integrated diagnostics; a larger number of spatial cells does not provide more packets per cell.

Define the moment error by

$$E_{\rm mom}^2= \frac{\|E-E_{\rm ex}\|_2^2 +\|F^x-F^x_{\rm ex}\|_2^2 +\|F^y-F^y_{\rm ex}\|_2^2} {\|E_{\rm ex}\|_2^2 +\|F^x_{\rm ex}\|_2^2 +\|F^y_{\rm ex}\|_2^2},$$

using cell-area weights. Also compute the time-integrated energy-budget residual from stored incoming and outgoing boundary fluxes, normalized by total injected energy. Accept the finest angular-transport result or the eight-seed Monte Carlo ensemble only if $E_{\rm mom}\leq2\times10^{-2}$ and the energy-budget residual is at most $10^{-3}$. The Monte Carlo standard error must decrease within $25\%$ of the expected packet-count scaling, and halving the CFL factor may change $E_{\rm mom}$ by no more than $10\%$. M1 is expected to fail the moment-error threshold near the crossing; that failure must be reported, not tuned away.

Observed discrepancy: not available—specification only.

Complete the ladder with homogeneous diffusion, free streaming along geodesics, radiation in moving matter, and absorption–emission equilibrium. Compare moment and Monte Carlo solutions with independently refined spatial, temporal, energy-group, angular, or packet resolution. For Monte Carlo, report seed ensembles and uncertainty. Track energy, momentum, and lepton-number budgets.

Then test a production-like snapshot or dynamical remnant while varying closure and reaction rates. Report luminosities, spectra, composition, heating, and fluid feedback with numerical and model uncertainty kept separate. The validity domain is the tested matter, optical depth, energy range, interactions, and duration.

Reproducing the calculation and testing its robustness

Pin transport equations, species, groups, opacity and reaction-rate datasets with checksums, closure or packet algorithm, random seeds, source integrator and tolerances, grid, boundaries, fluid coupling, precision, and code revision. A clean interface could be transport-benchmark --config crossing-beams.yaml --output run/, emitting fields, packets or moments, budgets, convergence or ensemble sweeps, failures, and environment metadata. This specifies a rerun record, not an executed program.

Reject a result that violates equilibrium, free-stream, diffusion, or exchange budgets; hides packet uncertainty; or applies an M1 beam result outside its closure regime. Review is triggered by new opacity or reaction-rate data, added oscillation or inelastic physics, changed closure, materially different remnant lifetime or optical-depth regime, or an independent comparison shifting a claimed observable beyond uncertainty. The frontier assessment must then be dated again.

Diagnostic

An M1 code converges cleanly on diffusion and single-beam tests but merges two crossing beams into one. A merger application’s heating region contains opposing beams. May spatial refinement alone support the heating rate?

Solution

No. The crossing-beam error is a closure limitation, not spatial truncation, so refinement converges to the wrong moment solution. Compare with angular transport or Monte Carlo, vary the closure, and include the resulting spread as model uncertainty. Narrow the claim if no independent angular calculation covers the relevant regime.

References

• Foucart, F. (2018). “Monte-Carlo closure for moment-based transport schemes in general relativistic radiation hydrodynamics simulations.” Monthly Notices of the Royal Astronomical Society 475, 4186–4207. DOI.
• Foucart, F. (2023). “Neutrino transport in general relativistic neutron star merger simulations.” Living Reviews in Computational Astrophysics 9, 1. DOI.
• Foucart, F., Cheong, P. C.-K., Duez, M. D., et al. (2024). “Robustness of neutron star merger simulations to changes in neutrino transport and neutrino-matter interactions.” Physical Review D 110, 083028. arXiv.