$Y_e$, entropy), the source term of Eq.(1) ($\mathcal{S}[{j}, {\chi}, {\Sigma}]$) can be determined. For closure, we use a prescribed relation between the radiation energy and flux as ($\mathcal{F}^{\rm s}/\mathcal{E}^{\rm s}= \frac{1}{2}(1 +\sqrt{1-
[{R_{\nu}}/\max\left(r,R_\nu\right)]^2})$ with $R_\nu$ being the radius of an energy-dependent scattering sphere (see Eq.(11) in @idsa). Since the cell-centered value of the flux, $\mathcal{F}^{\rm s}$, is obtained by the prescribed relation, the cell-interface value is estimated by the first-order upwind scheme assuming that the flux is out-going along the radial direction. With the numerical flux, the transport equation of $\mathcal{E}^{\rm s}$ (Eq.(1)) now expressed in a hyperbolic form is numerically solved. This modification does not produce any significant changes in the numerical results (see Takiwaki et al. in preparation for more details), however, the computational cost becomes more than 10 times smaller than that in the previous treatment. The velocity dependent terms ($O(v/c)$) are only included (up to the leading order) in the trapped part of the distribution function (Eq. (15) in @idsa). Concerning heavy-lepton neutrinos ($\nu_x = 
\nu_{\mu},\nu_{\tau},\bar{\nu}_{\mu}, \bar{\nu}_{\tau}$), we employ a leakage scheme to include the $\nu_x$ cooling via pair-, photo, plasma processes (e.g., @ross03 [@itoh89]). We apply the so-called ray-by-ray approach in which the neutrino transport 