(local) backward-Euler methods. The neutrino energy bins are coupled through (1) frame-dependent energy shifts since the material neutrinos interact with is moving, (2) gravitational redshift, and (3) energy transfer in scatterings off of electrons and nucleons. Neutrino-matter interaction rates are usually precomputed and stored in dense multi-D tables within which simulations interpolate.

    Full 6+1 D general-relativistic Boltzmann neutrino-radiation hydrodynamics is exceedingly challenging and has so far not been possible to included in core-collapse supernova simulations. 3+1 D (1D in space, 2D in momentum space) (e.g., [@liebendoerfer:05]), 5+1 D (2D in space, 3D in momentum space) simulations [@ott:08] and static 6D simulations [@sumiyoshi:15] have been carried out.

    Most (spatially) multi-D simulations treat neutrino transport in some dimensionally-reduced approximation. The most common is an expansion of the radiation field into angular moments. The $n$-th moment of this expansion requires information about the $(n+1)$-th moment (and in some cases also about the $(n+2)$-th moment). This necessitates a closure relation for the moment at which the expansion is truncated. Multi-group flux-limited diffusion evolves the $0$-th moment (the radiation energy density). The flux limiter is the closure that interpolates between diffusion and free streaming. The disadvantages of this method are its very diffusive nature that 