Yingda Cheng

Wei Guo, Yingda Cheng

In this paper, we develop a sparse grid discontinuous Galerkin (DG) scheme for transport equations and applied it to kinetic simulations. The method uses the weak formulations of traditional Runge-Kutta DG (RKDG) schemes for hyperbolic problems and is proven to be $L^2$ stable and convergent. A major advantage of the scheme lies in its low computational and storage cost due to the employed sparse finite element approximation space. This attractive feature is explored in simulating Vlasov and Boltzmann transport equations. Good performance in accuracy and conservation is verified by numerical tests in up to four dimensions.

Yingda Cheng, Irene M. Gamba, Philip J. Morrison

In this paper we consider Runge-Kutta discontinuous Galerkin (RKDG) schemes for Vlasov-Poisson systems that model collisionless plasmas. One-dimensional systems are emphasized. The RKDG method, originally devised to solve conservation laws, is seen to have excellent conservation properties, be readily designed for arbitrary order of accuracy, and capable of being used with a positivity-preserving limiter that guarantees positivity of the distribution functions. The RKDG solver for the Vlasov equation is the main focus, while the electric field is obtained through the classical representation by Green's function for the Poisson equation. A rigorous study of recurrence of the DG methods is presented by Fourier analysis, and the impact of different polynomial spaces and the positivity-preserving limiters on the quality of the solutions is ascertained. Several benchmark test problems, such as Landau damping, two-stream instability and the KEEN (Kinetic Electrostatic Electron Nonlinear) wave, are given.

Jose Morales-Escalante, Irene M. Gamba, Yingda Cheng et al.

The purpose of this work is to incorporate numerically, in a discontinuous Galerkin (DG) solver of a Boltzmann-Poisson model for hot electron transport, an electronic conduction band whose values are obtained by the spherical averaging of the full band structure given by a local empirical pseudopotential method (EPM) around a local minimum of the conduction band for silicon, as a midpoint between a radial band model and an anisotropic full band, in order to provide a more accurate physical description of the electron group velocity and conduction energy band structure in a semiconductor. This gives a better quantitative description of the transport and collision phenomena that fundamentally define the behaviour of the Boltzmann - Poisson model for electron transport used in this work. The numerical values of the derivatives of this conduction energy band, needed for the description of the electron group velocity, are obtained by means of a cubic spline interpolation. The EPM-Boltzmann-Poisson transport with this spherically averaged EPM calculated energy surface is numerically simulated and compared to the output of traditional analytic band models such as the parabolic and Kane bands, numerically implemented too, for the case of 1D $n^+-n-n^+$ silicon diodes with 400nm and 50nm channels. Quantitative differences are observed in the kinetic moments related to the conduction energy band used, such as mean velocity, average energy, and electric current (momentum).

Olivier Bokanowski, Yingda Cheng, Chi-Wang Shu

We study semi-Lagrangian discontinuous Galerkin (SLDG) and Runge-Kutta discontinuous Galerkin (RKDG) schemes for some front propagation problems in the presence of an obstacle term, modeled by a nonlinear Hamilton-Jacobi equation of the form $\min(u_t + c u_x, u - g(x))=0$, in one space dimension. New convergence results and error bounds are obtained for Lipschitz regular data. These "low regularity" assumptions are the natural ones for the solutions of the studied equations.

Yingda Cheng, Andrew J. Christlieb, Xinghui Zhong

In this paper, we generalize the idea in our previous work for the Vlasov-Ampère (VA) system \cite{cheng_va} and develop energy-conserving discontinuous Galerkin (DG) methods for the Vlasov-Maxwell (VM) system. The VM system is a fundamental model in the simulation of collisionless magnetized plasmas. Compared to \cite{cheng_va}, additional care needs to be taken for both the temporal and spatial discretizations to achieve similar type of conservation when the magnetic field is no longer negligible. Our proposed schemes conserve the total particle number and the total energy at the same time, and therefore can obtain accurate, yet physically relevant solutions. The main components of our methods include second order and above, explicit or implicit energy-conserving temporal discretizations, and DG methods for Vlasov and Maxwell's equations with carefully chosen numerical fluxes. Benchmark numerical tests such as the streaming Weibel instability are provided to validate the accuracy and conservation of the schemes.

Zhanjing Tao, Wei Guo, Yingda Cheng

In this paper, we develop sparse grid discontinuous Galerkin (DG) schemes for the Vlasov-Maxwell (VM) equations. The VM system is a fundamental kinetic model in plasma physics, and its numerical computations are quite demanding, due to its intrinsic high-dimensionality and the need to retain many properties of the physical solutions. To break the curse of dimensionality, we consider the sparse grid DG methods that were recently developed in \cite{guo2016sparse,guo2017adaptive} for transport equations. Such methods are based on multiwavelets on tensorized nested grids and can significantly reduce the numbers of degrees of freedom. We formulate two versions of the schemes: sparse grid DG and adaptive sparse grid DG methods for the VM system. Their key properties and implementation details are discussed. Accuracy and robustness are demonstrated by numerical tests, with emphasis on comparison of the performance of the two methods, as well as with their full grid counterparts.

Yingda Cheng, Irene M. Gamba, Fengyan Li et al.

Discontinuous Galerkin methods are developed for solving the Vlasov-Maxwell system, methods that are designed to be systematically as accurate as one wants with provable conservation of mass and possibly total energy. Such properties in general are hard to achieve within other numerical method frameworks for simulating the Vlasov-Maxwell system. The proposed scheme employs discontinuous Galerkin discretizations for both the Vlasov and the Maxwell equations, resulting in a consistent description of the distribution function and electromagnetic fields. It is proven, up to some boundary effects, that charge is conserved and the total energy can be preserved with suitable choices of the numerical flux for the Maxwell equations and the underlying approximation spaces. Error estimates are established for several flux choices. The scheme is tested on the streaming Weibel instability: the order of accuracy and conservation properties of the proposed method are verified.

Vrushali A. Bokil, Yingda Cheng, Yan Jiang et al.

The propagation of electromagnetic waves in general media is modeled by the time-dependent Maxwell's partial differential equations (PDEs), coupled with constitutive laws that describe the response of the media. In this work, we focus on nonlinear optical media whose response is modeled by a system of first order nonlinear ordinary differential equations (ODEs), which include a single resonance linear Lorentz dispersion, and the nonlinearity comes from the instantaneous electronic Kerr response and the residual Raman molecular vibrational response. To design efficient, accurate, and stable computational methods, we apply high order discontinuous Galerkin discretizations in space to the hybrid PDE-ODE Maxwell system with several choices of numerical fluxes, and the resulting semi-discrete methods are shown to be energy stable. Under some restrictions on the strength of the nonlinearity, error estimates are also established. When we turn to fully discrete methods, the challenge to achieve provable stability lies in the temporal discretizations of the nonlinear terms. To overcome this, novel strategies are proposed to treat the nonlinearity in our model within the framework of the second-order leap-frog and implicit trapezoidal time integrators. The performance of the overall algorithms are demonstrated through numerical simulations of kink and antikink waves, and third-harmonic generation in soliton propagation.

Anqi Chen, Fengyan Li, Yingda Cheng

In this paper, we develop an ultra-weak discontinuous Galerkin (DG) method to solve the one-dimensional nonlinear Schrödinger equation. Stability conditions and error estimates are derived for the scheme with a general class of numerical fluxes. The error estimates are based on detailed analysis of the projection operator associated with each individual flux choice. Depending on the parameters, we find out that in some cases, the projection can be defined element-wise, facilitating analysis. In most cases, the projection is global, and its analysis depends on the resulting $2\times2$ block-circulant matrix structures. For a large class of parameter choices, optimal $\textit{a priori}$ $L^2$ error estimates can be obtained. Numerical examples are provided verifying theoretical results.

Yuan Liu, Yingda Cheng, Shanqin Chen et al.

Computational costs of numerically solving multidimensional partial differential equations (PDEs) increase significantly when the spatial dimensions of the PDEs are high, due to large number of spatial grid points. For multidimensional reaction-diffusion equations, stiffness of the system provides additional challenges for achieving efficient numerical simulations. In this paper, we propose a class of Krylov implicit integration factor (IIF) discontinuous Galerkin (DG) methods on sparse grids to solve reaction-diffusion equations on high spatial dimensions. The key ingredient of spatial DG discretization is the multiwavelet bases on nested sparse grids, which can significantly reduce the numbers of degrees of freedom. To deal with the stiffness of the DG spatial operator in discretizing reaction-diffusion equations, we apply the efficient IIF time discretization methods, which are a class of exponential integrators. Krylov subspace approximations are used to evaluate the large size matrix exponentials resulting from IIF schemes for solving PDEs on high spatial dimensions. Stability and error analysis for the semi-discrete scheme are performed. Numerical examples of both scalar equations and systems in two and three spatial dimensions are provided to demonstrate the accuracy and efficiency of the methods. The stiffness of the reaction-diffusion equations is resolved well and large time step size computations are obtained.

Yingda Cheng, Andrew J. Christlieb, Wei Guo et al.

In plasma simulations, where the speed of light divided by a characteristic length is at a much higher frequency than other relevant parameters in the underlying system, such as the plasma frequency, implicit methods begin to play an important role in generating efficient solutions in these multi-scale problems. Under conditions of scale separation, one can rescale Maxwell's equations in such a way as to give a magneto static limit known as the Darwin approximation of electromagnetics. In this work, we present a new approach to solve Maxwell's equations based on a Method of Lines Transpose (MOL$^T$) formulation, combined with a fast summation method with computational complexity $O(N\log{N})$, where $N$ is the number of grid points (particles). Under appropriate scaling, we show that the proposed schemes result in asymptotic preserving methods that can recover the Darwin limit of electrodynamics.

Anqi Chen, Yingda Cheng, Yong Liu et al.

We analyze the superconvergence properties of ultra-weak discontinuous Galerkin (UWDG) methods with various choices of flux parameters for one-dimensional linear Schrödinger equation. In our previous work [10], stability and optimal convergence rate are established for a large class of flux parameters. Depending on the flux choices and if the polynomial degree $k$ is even or odd, in this paper, we prove $2k$ or $(2k-1)$-th order superconvergence rate for cell averages and numerical flux of the function, as well as $(2k-1)$ or $(2k-2)$-th order for numerical flux of the derivative. In addition, we prove superconvergence of $(k+2)$ or $(k+3)$-th order of the DG solution towards a special projection. At a class of special points, the function values and the first and second order derivatives of the DG solution are superconvergent with order $k+2, k+1, k$, respectively. The proof relies on the correction function techniques initiated in [8], and applied to [6] for direct DG (DDG) methods for diffusion problems. Compared with [6], Schrödinger equation poses unique challenges for superconvergence proof because of the lack of the dissipation mechanism from the equation. One major highlight of our proof is that we introduce specially chosen test functions in the error equation and show the superconvergence of the second derivative and jump across the cell interfaces of the difference between numerical solution and projected exact solution. This technique was originally proposed in [12] and is essential to elevate the convergence order for our analysis. Finally, by negative norm estimates, we apply the post-processing technique and show that the accuracy of our scheme can be enhanced to order $2k.$ Theoretical results are verified by numerical experiments.

Wei Guo, Yingda Cheng

In this paper, we develop an adaptive multiresolution discontinuous Galerkin (DG) scheme for time-dependent transport equations in multi-dimensions. The method is constructed using multiwavlelets on tensorized nested grids. Adaptivity is realized by error thresholding based on the hierarchical surplus, and the Runge-Kutta DG (RKDG) scheme is employed as the reference time evolution algorithm. We show that the scheme performs similarly to a sparse grid DG method when the solution is smooth, reducing computational cost in multi-dimensions. When the solution is no longer smooth, the adaptive algorithm can automatically capture fine local structures. The method is therefore very suitable for deterministic kinetic simulations. Numerical results including several benchmark tests, the Vlasov-Poisson (VP) and oscillatory VP systems are provided.

Zhanjing Tao, Anqi Chen, Mengping Zhang et al.

In this paper, we develop sparse grid central discontinuous Galerkin (CDG) scheme for linear hyperbolic systems with variable coefficients in high dimensions. The scheme combines the CDG framework with the sparse grid approach, with the aim of breaking the curse of dimensionality. A new hierarchical representation of piecewise polynomials on the dual mesh is introduced and analyzed, resulting in a sparse finite element space that can be used for non-periodic problems. Theoretical results, such as $L^2$ stability and error estimates are obtained for scalar problems. CFL conditions are studied numerically comparing discontinuous Galerkin (DG), CDG, sparse grid DG and sparse grid CDG methods. Numerical results including scalar linear equations, acoustic and elastic waves are provided.

Peter DelMastro, Daniel Appelö, Yingda Cheng

We propose a family of low-rank, completely positive and trace preserving schemes for the Lindblad equation, a common model for open quantum systems. Low-rank representation is employed at two levels: the density matrix is factorized into the product of tall-skinny matrices, and the columns of these matrices are further represented using the tensor train (TT) format, also know as matrix product states (MPS). This two-level low-rank format fits naturally into our existing Kraus is King scheme (arXiv:2409.08898v2 [math.NA]) for the Lindblad equation, whose underlying operations are arithmetic on the columns of the tall-skinny matrices. We show how these operations can be performed efficiently in the TT/MPS format, with particular emphasis on density matrix rank-truncation. We conclude with extensive numerical experiments demonstrating the convergence of this scheme and its efficiency in simulating systems with up to $10^{19}$ degrees of freedom using only modest compute resources.

Juntao Huang, Yingda Cheng, Andrew J. Christlieb et al.

This is the third paper in a series in which we develop machine learning (ML) moment closure models for the radiative transfer equation (RTE). In our previous work \cite{huang2021gradient}, we proposed an approach to learn the gradient of the unclosed high order moment, which performs much better than learning the moment itself and the conventional $P_N$ closure. However, while the ML moment closure has better accuracy, it is not able to guarantee hyperbolicity and has issues with long time stability. In our second paper \cite{huang2021hyperbolic}, we identified a symmetrizer which leads to conditions that enforce that the gradient based ML closure is symmetrizable hyperbolic and stable over long time. The limitation of this approach is that in practice the highest moment can only be related to four, or fewer, lower moments. In this paper, we propose a new method to enforce the hyperbolicity of the ML closure model. Motivated by the observation that the coefficient matrix of the closure system is a lower Hessenberg matrix, we relate its eigenvalues to the roots of an associated polynomial. We design two new neural network architectures based on this relation. The ML closure model resulting from the first neural network is weakly hyperbolic and guarantees the physical characteristic speeds, i.e., the eigenvalues are bounded by the speed of light. The second model is strictly hyperbolic and does not guarantee the boundedness of the eigenvalues. Several benchmark tests including the Gaussian source problem and the two-material problem show the good accuracy, stability and generalizability of our hyperbolic ML closure model.

Juntao Huang, Yingda Cheng, Andrew J. Christlieb et al.

This is the second paper in a series in which we develop machine learning (ML) moment closure models for the radiative transfer equation (RTE). In our previous work \cite{huang2021gradient}, we proposed an approach to directly learn the gradient of the unclosed high order moment, which performs much better than learning the moment itself and the conventional $P_N$ closure. However, the ML moment closure model in \cite{huang2021gradient} is not able to guarantee hyperbolicity and long time stability. We propose in this paper a method to enforce the global hyperbolicity of the ML closure model. The main idea is to seek a symmetrizer (a symmetric positive definite matrix) for the closure system, and derive constraints such that the system is globally symmetrizable hyperbolic. It is shown that the new ML closure system inherits the dissipativeness of the RTE and preserves the correct diffusion limit as the Knunsden number goes to zero. Several benchmark tests including the Gaussian source problem and the two-material problem show the good accuracy, long time stability and generalizability of our globally hyperbolic ML closure model.

Juntao Huang, Yingda Cheng, Andrew J. Christlieb et al.

In this paper, we take a data-driven approach and apply machine learning to the moment closure problem for radiative transfer equation in slab geometry. Instead of learning the unclosed high order moment, we propose to directly learn the gradient of the high order moment using neural networks. This new approach is consistent with the exact closure we derive for the free streaming limit and also provides a natural output normalization. A variety of benchmark tests, including the variable scattering problem, the Gaussian source problem with both periodic and reflecting boundaries, and the two-material problem, show both good accuracy and generalizability of our machine learning closure model.

Yan Jiang, Puttha Sakkaplangkul, Vrushali A. Bokil et al.

In this paper, we consider Maxwell's equations in linear dispersive media described by a single-pole Lorentz model for electronic polarization. We study two classes of commonly used spatial discretizations: finite difference methods (FD) with arbitrary even order accuracy in space and high spatial order discontinuous Galerkin (DG) finite element methods. Both types of spatial discretizations are coupled with second order semi-implicit leap-frog and implicit trapezoidal temporal schemes studied in our previous research [5,6]. By performing detailed dispersion analysis for the semi-discrete and fully discrete schemes, we obtain rigorous quantification of the dispersion error for Lorentz dispersive dielectrics. In particular, comparisons of dispersion error can be made taking into account the model parameters, and mesh sizes in the design of the two types of schemes. The results for the numerical dispersion analysis can guide us in the optimal choice of discretization parameters for the more complicated and nonlinear models. The numerical dispersion analysis of the fully discrete FD and DG schemes, for the dispersive Maxwell model considered in this paper, clearly indicate the dependence of the numerical dispersion errors on spatial and temporal discretizations, their order of accuracy, mesh discretization parameters and model parameters. The results obtained here cannot be arrived at by considering discretizations of Maxwell's equations in free space. In particular, our results contrast the advantages and disadvantages of using high order FD or DG schemes and leap-frog or trapezoidal time integrators over different frequency ranges using a variety of measures of numerical dispersion errors. Finally, we highlight the limitations of the second order accurate temporal discretizations considered.

Zixuan Wang, Qi Tang, Wei Guo et al.

This paper constitutes our initial effort in developing sparse grid discontinuous Galerkin (DG) methods for high-dimensional partial differential equations (PDEs). Over the past few decades, DG methods have gained popularity in many applications due to their distinctive features. However, they are often deemed too costly because of the large number of degrees of freedom of the approximation space, which are the main bottleneck for simulations in high dimensions. In this paper, we develop sparse grid DG methods for elliptic equations with the aim of breaking the \emph{curse of dimensionality}. Using a hierarchical basis representation, we construct a sparse finite element approximation space, reducing the degrees of freedom from the standard {$O(h^{-d})$ to $O(h^{-1}|\log_2 h|^{d-1})$} for $d$-dimensional problems, where $h$ is the uniform mesh size in each dimension. Our method, based on the interior penalty (IP) DG framework, can achieve accuracy of $O(h^{k}|\log_2 h|^{d-1})$ in the energy norm, where $k$ is the degree of polynomials used. Error estimates are provided and confirmed by numerical tests in multi-dimensions.

Yingda Cheng, Irene M. Gamba, Armando Majorana et al.

In this paper, we present results of a discontinuous Galerkin (DG) scheme applied to deterministic computations of the transients for the Boltzmann-Poisson system describing electron transport in semiconductor devices. The collisional term models optical-phonon interactions which become dominant under strong energetic conditions corresponding to nano-scale active regions under applied bias. The proposed numerical technique is a finite element method using discontinuous piecewise polynomials as basis functions on unstructured meshes. It is applied to simulate hot electron transport in bulk silicon, in a silicon $n^+$-$n$-$n^+$ diode and in a double gated 12nm MOSFET. Additionally, the obtained results are compared to those of a high order WENO scheme simulation and DSMC (Discrete Simulation Monte Carlo) solvers.