## 2 Mathematical Formulation
### Discretization of Transient Radiative Transfer Equation (TRTE)  
Equation (1) can be expressed (with linear anisotropy) as  
\[\frac{1}{\beta\cdot c_{med}}\frac{\partial l_{d}^{l}}{\partial t}+\frac{1}{\beta} \frac{\partial l_{d}^{l}}{\partial\mathsf{s}}\ =\ -l_{d}^{l}+(1-\omega)I_{b}+\frac{\omega}{4\pi}\sum_{l^{ \prime}}^{L}l_{d}^{l}\Phi^{l^{\prime}l}\Delta\Omega^{l}+S_{c}^{l} \tag{5}\]  
where \(L\) defines the total number of discrete solid angles. In the present work, the scattering phase function is expressed as a series in terms of [PERSON]'s polynomial [23]. Thus, the linearized TRTE is given by  
\[\frac{1}{\beta\cdot c_{md}}\frac{\partial l_{d}^{l}}{\partial t}+\frac{1}{ \beta}\frac{\partial l_{d}^{l}}{\partial\mathsf{s}}\ =\ -\beta_{m}^{l}l_{d}^{l}+S_{m}^{l} \tag{6}\]  
where the modified extinction coefficient \(\beta_{m}\) is given as  
\[\beta_{m}^{l}\ =\ \left(1-\frac{\omega}{4\pi}\Delta\Omega^{l}\right)-\frac{ \omega}{4\pi}g\cos^{2}\theta^{l}\Delta\Omega^{l}-g\sin^{2}\theta^{l}\Delta \Omega^{l} \tag{7}\]  
The modified source term is given by  
\[S_{m}^{l}\ =\ (1-\omega)I_{b}+S_{c}^{l}+\left[\begin{array}{ccc}\frac{ \omega}{4\pi}&\sum\limits_{l}^{L}&l_{d}^{l^{\prime}}\Delta\Omega^{l^{\prime}}+ g\cos\theta&\sum\limits_{l}^{L}&l_{d}^{l^{\prime}}\cos\theta^{l^{\prime}}\Delta \Omega^{l^{\prime}}+\\ l^{\prime}&=\ 1&&l^{\prime}\ =\ 1\\ l^{\prime}\
eq l&&l^{\prime}\
eq l\end{array}\right] \tag{8}\]  
The FTn FVM used \(n(n+2)\) control angles where the polar angle was uniformly subdivided into an even number (n). However, a non-uniform subdivision of the azimuthal angle was considered. In fact, for each level of the polar angle, the azimuthal angle was subdivided into the numbers of the following arrangement: 4,8,12... 2n,2n... 12,8,4 [14; 15; 16; 19].  
Different spatial discretization schemes may be applied to relate radiation intensities on the interfaces to the radiation intensity at the grid nodes. The curved-line advection method (CLAM) scheme differs from the STEP scheme by the interpolation of intensities in a specified interface [9; 10; 11; 12; 13; 14; 15; 16]. These differences have been analyzed in [9; 10; 11; 12; 13; 14; 15; 16].  
The forward Euler scheme (or the so-called the explicit scheme) is a first-order numerical procedure in which the temporal dependent variable in the time instant (\(t+\Delta t\)) at a given point (e.g., the radiation intensity at time instant \(t+\Delta t\): \(l^{t+\Delta t}\) for the TRTE) can be calculated explicitly from the current time instant (\(t\)) solution values (e.g., \(l^{t}\) for the TRTE) according to (for the TRTE example):  
\[\frac{1}{c_{med}}\frac{\partial I}{\partial t}\approx\frac{l^{t+\Delta t}-l^{ t}}{c_{mal}\Delta t}\ =\ F\Big{(}t,l^{t}\Big{)} \tag{9}\]  
where \(\Delta t\) is the time step and \(F\) is the temporal derivative of \(I\) at instant \(t\), divided by the propagation speed of light within the medium, \(C_{med}\). This temporal scheme is easier to apply and computationally more efficient than the backward Euler (implicit) scheme, where the solution \(l^{t+\Delta t}\) is determined by resolving a set of equations containing both the current and latter state. It is also noted that the explicit scheme needs a smaller time step \(\Delta t\), by comparison to the implicit scheme, in order to reach numerical stability.