12:00 PM - 12:15 PM
[MGI33-11] A study on the influence of numerical errors with discontinuous Galerkin method in atmospheric large-eddy simulations
Keywords:atmospheric boundary layer, large-eddy simulation, high-order fluid scheme, numerical errors
Introduction
One of the issues for future large-eddy simulation (LES) of global atmosphere is the discretization accuracy of the dynamical processes. Kawai and Tomita (2021, hereinafter KT2021), in the framework of the finite difference method (FDM), discussed the order of accuracy required for discretization of the advection term not to dominate the eddy viscosity term associated with the turbulence model. Based on the numerical criteria derived, it is suggested that seven to eight orders of accuracy are needed. If we attempt to achieve such a high accuracy in the conventional grid point methods, the complexity of the discretization and the deterioration of computational locality due to the extended stencil will be a problem. Thus, we now focus on the discontinuous Galerkin method (DGM), which is characterized by the simplicity of high-order strategy and computational compactness. To investigate the suitability of DGM to atmospheric calculations, we constructed a regional LES model using DGM and conducted numerical experiments of atmospheric boundary layer turbulence. It is confirmed that the vertical structure and energy spectrum well reproduced the results of SCALE-RM based on the conservative FDM using the higher-order advection scheme (reported in JpGU 2021). However, our mathematical interpretation of the numerical errors on the spectra impact in the short wavelength region was insufficient. Thus, we extend the numerical indices derived in KT2021 to the framework of DGM and investigate the effect of numerical errors in the short wavelength region of DGM. We also clarify the required order of the expansion polynomial (p) in the manner same as in KT2021. In this presentation, we will show the preliminary results.
Formulation of indices with numerical errors
An index for numerical viscosity errors is defined as the ratio of decay time with Laplacian term of the eddy viscosity term to that with the numerical viscosity (Rdiff). On the other hand, an index for numerical dispersion is defined as the ratio of the phase speed with the numerical dispersion to that with the cross term of the eddy viscosity term (Rdisp). In order for the numerical error not to contaminate the eddy viscosity term, they must be sufficiently smaller than unity. In FDM, the decay coefficient and the phase speed with the numerical error terms can be obtained explicitly because the corresponding modified equation can be easily derived, but such derivation is too difficult in high-order DGM. Thus, we adopt Fourier eigenvalue analysis, following Moura et al (2015). In DGM, multiple eigenvalues can be obtained because there are many degrees of freedom within one finite element. Considering coupled modes, we investigate the net behavior of these modes. Finally, comparing the exact solution with the numerical solution enables us to quantify amplification factor and phase error and estimate decay time constant with numerical viscosity and the phase velocity with numerical dispersion. Figs. (a,b) show that the grid spacing and wavelength dependence of Rdiff and Rdisp for DGM with p=3 and p=7 used upwind numerical fluxes . When we require that Rdiff and Rdisp < 10-1 is required at wavelengths longer than 8 grids, at least p=3 is required in O(10 m) grid spacing.
Validation of numerical criteria
To verify the numerical indices derived semi-theoretically, we conduct numerical experiments of planetary boundary layer turbulence (based on Nishizawa et al., 2015) using our developing atmospheric LES model with DGM. The effective grid spacing is fixed at 10 m, and hexahedral elements with p=3 and p=7 are used. The spatial filter length with twice the effective grid spacing results in the energy spectra obeying -5/3 power law up to about 8 grids. At higher wavenumber, both the spatial filter and the numerical viscosity contribute to the steep slope. Fig. (c) reveals that p=3 is acceptable because the spectra follow the -5/3 power law up to about 8 grids and satisfies our criterion. At shorter wavelengths, the effect of numerical diffusion becomes significant. These features are consistent to the indication by Rdiff. The numerical criteria derived here will be useful to consider the required order of polynomial and the contributions of the numerical stabilization mechanism with DGM.
One of the issues for future large-eddy simulation (LES) of global atmosphere is the discretization accuracy of the dynamical processes. Kawai and Tomita (2021, hereinafter KT2021), in the framework of the finite difference method (FDM), discussed the order of accuracy required for discretization of the advection term not to dominate the eddy viscosity term associated with the turbulence model. Based on the numerical criteria derived, it is suggested that seven to eight orders of accuracy are needed. If we attempt to achieve such a high accuracy in the conventional grid point methods, the complexity of the discretization and the deterioration of computational locality due to the extended stencil will be a problem. Thus, we now focus on the discontinuous Galerkin method (DGM), which is characterized by the simplicity of high-order strategy and computational compactness. To investigate the suitability of DGM to atmospheric calculations, we constructed a regional LES model using DGM and conducted numerical experiments of atmospheric boundary layer turbulence. It is confirmed that the vertical structure and energy spectrum well reproduced the results of SCALE-RM based on the conservative FDM using the higher-order advection scheme (reported in JpGU 2021). However, our mathematical interpretation of the numerical errors on the spectra impact in the short wavelength region was insufficient. Thus, we extend the numerical indices derived in KT2021 to the framework of DGM and investigate the effect of numerical errors in the short wavelength region of DGM. We also clarify the required order of the expansion polynomial (p) in the manner same as in KT2021. In this presentation, we will show the preliminary results.
Formulation of indices with numerical errors
An index for numerical viscosity errors is defined as the ratio of decay time with Laplacian term of the eddy viscosity term to that with the numerical viscosity (Rdiff). On the other hand, an index for numerical dispersion is defined as the ratio of the phase speed with the numerical dispersion to that with the cross term of the eddy viscosity term (Rdisp). In order for the numerical error not to contaminate the eddy viscosity term, they must be sufficiently smaller than unity. In FDM, the decay coefficient and the phase speed with the numerical error terms can be obtained explicitly because the corresponding modified equation can be easily derived, but such derivation is too difficult in high-order DGM. Thus, we adopt Fourier eigenvalue analysis, following Moura et al (2015). In DGM, multiple eigenvalues can be obtained because there are many degrees of freedom within one finite element. Considering coupled modes, we investigate the net behavior of these modes. Finally, comparing the exact solution with the numerical solution enables us to quantify amplification factor and phase error and estimate decay time constant with numerical viscosity and the phase velocity with numerical dispersion. Figs. (a,b) show that the grid spacing and wavelength dependence of Rdiff and Rdisp for DGM with p=3 and p=7 used upwind numerical fluxes . When we require that Rdiff and Rdisp < 10-1 is required at wavelengths longer than 8 grids, at least p=3 is required in O(10 m) grid spacing.
Validation of numerical criteria
To verify the numerical indices derived semi-theoretically, we conduct numerical experiments of planetary boundary layer turbulence (based on Nishizawa et al., 2015) using our developing atmospheric LES model with DGM. The effective grid spacing is fixed at 10 m, and hexahedral elements with p=3 and p=7 are used. The spatial filter length with twice the effective grid spacing results in the energy spectra obeying -5/3 power law up to about 8 grids. At higher wavenumber, both the spatial filter and the numerical viscosity contribute to the steep slope. Fig. (c) reveals that p=3 is acceptable because the spectra follow the -5/3 power law up to about 8 grids and satisfies our criterion. At shorter wavelengths, the effect of numerical diffusion becomes significant. These features are consistent to the indication by Rdiff. The numerical criteria derived here will be useful to consider the required order of polynomial and the contributions of the numerical stabilization mechanism with DGM.