## 3 Simulation Results  
In this section, to demonstrate the validity for the numerical approach, we first compare the results for the DOC obtained with our numerical approach with those derived from an analytical expression. Then, we present numerical examples for the evolution of the DOC with the propagation distance, using our numerical method. In this numerical case, it is impossible to obtain the analytical expression of the GPCB during propagation. In the following simulation, the number of elements \(N_{1}\times N_{1}\) and \(N\) are chosen to be \(512\times 512\) and \(10^{4}\), respectively. Under this circumstance, it takes about 3.7 min for each plot. The computer is from Dell, Texas, USA, whose configurations used in the simulation are i7-8700 CPU @ 3.2 GHz and 16 GB RAM. Here, we want to stress that the advantage of our numerical approach is that it can deal with partially coherent beams that cannot be described by an analytical expression. In order to demonstrate the validity of our simulation method, a typical Schell-model beam, named a partially coherent elegant Hermite-Gaussian (PCEHG) beam [42], is used to compare the evolution of DOC from the analytical results with that from the simulation results. The CSD function of the PCEHG beam in the source plane is expressed as  
\[ \begin{array}{l}W_{0}(\mathbf{r}_{1},\mathbf{r}_{2})=\tau_{0}(\mathbf{r}_{1 })\tau_{0}^{*}(\mathbf{r}_{2})\mu_{0}(\mathbf{r}_{1}-\mathbf{r}_{2})\\ =H_{n}\Big{(}\frac{\gamma_{1}}{\omega_{0}}\Big{)}H_{m}\Big{(}\frac{\gamma_{1} }{\omega_{0}}\Big{)}H_{n}\Big{(}\frac{\gamma_{2}}{\omega_{0}}\Big{)}H_{m}\Big{(} \frac{\gamma_{2}}{\omega_{0}}\Big{)}\\ \times\exp\Big{(}-\frac{\mathbf{r}_{1}^{2}+\mathbf{r}_{2}^{2}}{\omega_{0}^{2}} \Big{)}\exp\Big{(}-\frac{\mathbf{r}_{1}-\mathbf{r}_{2}^{2}}{2\sigma_{g}^{2}} \Big{)}.\end{array}  \tag{13}\]  
\(H_{n}\) is the Hermite polynomial of order \(n\), \(\omega_{0}\) and \(\sigma_{g}\) are the beam waist size and the transverse coherence width, respectively. The analytical expression for the CSD function of the PCEHG beam at distance z in free-space propagation can be found in [42] (see Equation (17) in [42], where the elements of the transfer matrix are \(A=1\), \(B=z\), and \(D=1\)). Figure 1a1-e1 illustrates the evolution of the modulus of the square of the DOC with the propagation distance in free space calculated from the analytical formula. The beam parameters in the calculation are chosen to be \(n=2\), \(m=0\), \(\lambda=532\)nm, \(\omega_{0}=1\) mm, and \(\sigma_{g}=1\) mm. For comparison, the corresponding results attained by our numerical approach are shown in Figure 1a2-e2. It can be seen that the simulation results agree reasonably well with those obtained from the analytical formula, implying that our method is valid and effective in simulating the evolution of the DOC during propagation. As we state in Section 2, the advantage of the numerical approach is that it can deal with the evolution of the DOC of a GPCB whose CSD function during propagation is impossible to acquire by an analytical expression. In the following example, we assume that the function \(p(\mathbf{v})\) (Fourier  
Figure 1: Evolution of the square of the modulus of the degree of coherence of the partially coherent elegant Hermite–Gaussian beam with the propagation distance in free space. (**a1–e1**) Results obtained from the analytical formula; (**a2–e2**) the corresponding results obtained using our numerical approach. transform of the DOC) is a soccer pattern, as shown in Figure 2a; the corresponding modulus of the square of the DOC is plotted in Figure 2c. The initial amplitude function \(\tau_{0}\) we choose is the elegant Hermite-Gaussian distribution, as follows  
\[\tau_{0}(x,y)=H_{n}(x/\omega_{0})H_{m}(y/\omega_{0})\exp\Bigl{[}-\Bigl{(}x^{2}+ y^{2}\Bigr{)}/\omega_{0}^{2}\Bigr{]}. \tag{14}\]  
Figure 2b,d presents the source intensity distribution of such GPCB and its modulus of the square of the Fourier transform, respectively.