# Characteristic Length and Time Scales of the Highly Forward Scattering of Photons in Random Media
## 2 Materials and Methods
### Time-Domain Photon Transport Models
#### 2.1.3 \(M\)-th Order Delta-Eddington Equation (dEM)  
The dEA decomposes the highly forward-peaked phase function into a purely forward-peaked component expressed by the delta function and the other components [16; 17]:  
\[p_{\delta}^{M}(\mathbf{\Omega}\cdot\mathbf{\Omega}^{\prime})=\frac{1}{2\pi}h\delta(1- \mathbf{\Omega}\cdot\mathbf{\Omega}^{\prime})+(1-h)p_{\delta\delta}^{M}(\mathbf{\Omega} \cdot\mathbf{\Omega}^{\prime}), \tag{4}\]  
where \(h\) is a coefficient of the decomposition. \(p_{\delta 2}^{M}\) is a phase function excluding the delta-function component and expanded in the finite series of Legendre polynomials up to the order of \(M\):  
\[p_{\delta 2}^{M}(\mathbf{\Omega}\cdot\mathbf{\Omega}^{\prime})=\sum_{l=0}^{M}\frac{2l+1} {4\pi}\sigma_{dEl}P_{l}(\mathbf{\Omega}\cdot\mathbf{\Omega}^{\prime}). \tag{5}\]  
In the case of the HG phase function, the expansion coefficient \(\sigma_{dEl}\) and the decomposition coefficient \(h\) are determined so as to satisfy the moment conditions up to the order of \(M+1\):  
\[\sigma_{dEl}:=\int_{\mathbb{S}^{2}}d\mathbf{\Omega}^{\prime}p_{\delta 2}^{M}(\mathbf{\Omega}\cdot\mathbf{\Omega}^{ \prime})P_{l}(\mathbf{\Omega}\cdot\mathbf{\Omega}^{\prime})=\frac{g^{l}-h}{1-h},\ \ \ \ \ h=g^{M+1}. \tag{6}\]  
The case of \(l=0\) in Equation (6) corresponds to the normalization condition of \(p_{\delta 2}^{M}\cdot\int_{\mathbb{S}^{2}}d\mathbf{\Omega}^{\prime}p_{\delta 2}^{M}(\mathbf{ \Omega}\cdot\mathbf{\Omega}^{\prime})=1\). By using \(p_{\delta}^{M}\) (Equation (4)), the RTE (Equation (1)) is approximated to  
\[\left[\frac{\partial}{v\partial t}+\mathbf{\Omega}\cdot\
abla+\mu_{a}(\mathbf{r})+\mu _{s}^{M}(\mathbf{r})\right]I(\mathbf{r},\mathbf{\Omega},t)=\mu_{s}^{M}(\mathbf{r})\int_{ \mathbb{S}^{2}}d\mathbf{\Omega}^{\prime}p_{\delta 2}^{M}(\mathbf{\Omega}\cdot\mathbf{\Omega}^{ \prime})I(\mathbf{r},\mathbf{\Omega}^{\prime},t)+q(\mathbf{r},\mathbf{\Omega},t), \tag{7}\]  
where \(\mu_{s}^{M}(\mathbf{r})=(1-h)\mu_{s}(\mathbf{r})=(1-g^{M+1})\mu_{s}(\mathbf{r})\). In this paper, Equation (7) is called as the \(M\)-th order delta-Eddington equation (dEM). Comparing the RTE (Equation (1)) with the dEM, \(\mu_{s}\) is modified to \(\mu_{s}^{M}\) and \(p\) is approximated to \(p_{\delta 2}^{M}\). The modification of \(\mu_{s}\) to \(\mu_{s}^{M}\) means a change of the scattering length scale. Figure 2a,b plot a scattering length ratio of \(\mu_{s}^{M}/\mu_{s}\) and a relative error \(E_{phase}^{M}\) of the phase functions between the HG function (Equation (2)) and the \(M\)-th order dEA (Equation (5)) at different \(g\)-values of 0.3, 0.6, 0.9, and 0.95. Here, the error \(E_{phase}^{M}\) is defined by  
\[E_{phase}^{M}=\int_{-1}^{1}d\mu\left|\frac{p_{\delta 2}^{M}(\mu)-p_{HG}(\mu)}{p_{HG }(\mu)}\right|,\ \ \ \ \ \mu=\mathbf{\Omega}\cdot\mathbf{\Omega}^{\prime}. \tag{8}\]  
Figure 1: Henyey-Greenstein phase function as a function of \(\mathbf{\Omega}\cdot\mathbf{\Omega}^{\prime}\) for isotropic scattering (\(g=0.0\)), moderately forward scattering (\(g=0.6\)), and highly forward scattering (\(g=0.9\) and 0.95), in a logarithmic scale.  
As the expansion order \(M\) increases, the ratio \(\mu_{s}^{M}/\mu_{s}\) approaches to unity and the error \(E_{phase}^{M}\) decreases, meaning the scattering length scale and the phase function for the dEM converge to those for the RTE. At \(g=0.9\) and 0.95 (highly forward scattering), the convergences are slower than the cases of \(g=0.3\) and 0.6 (moderate forward scattering).