Hailiang Liu

Hailiang Liu, Zhongming Wang

We design an arbitrary-order free energy satisfying discontinuous Galerkin (DG) method for solving time-dependent Poisson-Nernst-Planck systems. Both the semi-discrete and fully discrete DG methods are shown to satisfy the corresponding discrete free energy dissipation law for positive numerical solutions. Positivities of numerical solutions are enforced by an accuracy-preserving limiter in reference to positive cell averages. Numerical examples are presented to demonstrate the high resolution of the numerical algorithm and to illustrate the proven properties of mass conservation, free energy dissipation, as well as the preservation of steady states.

Hailiang Liu, Olof Runborg, Nicolay M. Tanushev

Gaussian beams are asymptotically valid high frequency solutions to hyperbolic partial differential equations, concentrated on a single curve through the physical domain. They can also be extended to some dispersive wave equations, such as the Schrödinger equation. Superpositions of Gaussian beams provide a powerful tool to generate more general high frequency solutions that are not necessarily concentrated on a single curve. This work is concerned with the accuracy of Gaussian beam superpositions in terms of the wavelength $ε$. We present a systematic construction of Gaussian beam superpositions for all strictly hyperbolic and Schrödinger equations subject to highly oscillatory initial data of the form $Ae^{iΦ/ε}$. Through a careful estimate of an oscillatory integral operator, we prove that the $k$-th order Gaussian beam superposition converges to the original wave field at a rate proportional to $ε^{k/2}$ in the appropriate norm dictated by the well-posedness estimate. In particular, we prove that the Gaussian beam superposition converges at this rate for the acoustic wave equation in the standard, $ε$-scaled, energy norm and for the Schrödinger equation in the $L^2$ norm. The obtained results are valid for any number of spatial dimensions and are unaffected by the presence of caustics. We present a numerical study of convergence for the constant coefficient acoustic wave equation in $\Real^2$ to analyze the sharpness of the theoretical results.

Hailiang Liu, Zhongming Wang

We propose a high order discontinuous Galerkin (DG) method for solving nonlinear Fokker-Planck equations with a gradient flow structure. For some of these models it is known that the transient solutions converge to steady-states when time tends to infinity. The scheme is shown to satisfy a discrete version of the entropy dissipation law and preserve steady-states, therefore providing numerical solutions with satisfying long-time behavior. The positivity of numerical solutions is enforced through a reconstruction algorithm, based on positive cell averages. For the model with trivial potential, a parameter range sufficient for positivity preservation is rigorously established. For other cases, cell averages can be made positive at each time step by tuning the numerical flux parameters. A selected set of numerical examples is presented to confirm both the high-order accuracy and the efficiency to capture the large-time asymptotic.

Yi Jiang, Hailiang Liu

An invariant-region-preserving (IRP) limiter for multi-dimensional hyperbolic conservation law systems is introduced, as long as the system admits a global invariant region which is a convex set in the phase space. It is shown that the order of approximation accuracy is not destroyed by the IRP limiter, provided the cell average is away from the boundary of the convex set. Moreover, this limiter is explicit, and easy for computer implementation. A generic algorithm incorporating the IRP limiter is presented for high order finite volume type schemes. For arbitrarily high order discontinuous Galerkin (DG) schemes to hyperbolic conservation law systems, sufficient conditions are obtained for cell averages to remain in the invariant region provided the projected one-dimensional system shares the same invariant region as the full multi-dimensional hyperbolic system {does}. The general results are then applied to both one and two dimensional compressible Euler equations so to obtain high order IRP DG schemes. Numerical experiments are provided to validate the proven properties of the IRP limiter and the performance of IRP DG schemes for compressible Euler equations.

Hailiang Liu, James Ralston, Olof Runborg et al.

In this work we construct Gaussian beam approximations to solutions of the high frequency Helmholtz equation with a localized source. Under the assumption of non-trapping rays we show error estimates between the exact outgoing solution and Gaussian beams in terms of the wave number $k$, both for single beams and superposition of beams. The main result is that the relative local $L^2$ error in the beam approximations decay as {$k^{-N/2}$ independent of dimension and presence of caustics, for $N$-th order beams.

Hailiang Liu, Xuping Tian

We introduce a novel algorithm for gradient-based optimization of stochastic objective functions. The method may be seen as a variant of SGD with momentum equipped with an adaptive learning rate automatically adjusted by an 'energy' variable. The method is simple to implement, computationally efficient, and well suited for large-scale machine learning problems. The method exhibits unconditional energy stability for any size of the base learning rate. We provide a regret bound on the convergence rate under the online convex optimization framework. We also establish the energy-dependent convergence rate of the algorithm to a stationary point in the stochastic non-convex setting. In addition, a sufficient condition is provided to guarantee a positive lower threshold for the energy variable. Our experiments demonstrate that the algorithm converges fast while generalizing better than or as well as SGD with momentum in training deep neural networks, and compares also favorably to Adam.

Hailiang Liu, Robin Gröpler, Gerald Warnecke

We design, analyze and numerically validate a novel discontinuous Galerkin method for solving the coagulation-fragmentation equations. The DG discretization is applied to the conservative form of the model, with flux terms evaluated by Gaussian quadrature with $Q=k+1$ quadrature points for polynomials of degree $k$. The positivity of the numerical solution is enforced through a simple scaling limiter based on positive cell averages. The positivity of cell averages is propagated by the time discretization provided a proper time step restriction is imposed.

Yi Jiang, Hailiang Liu

In this paper, we introduce an invariant-region-preserving (IRP) limiter for the p-system and the corresponding viscous p-system, both of which share the same invariant region. Rigorous analysis is presented to show that for smooth solutions the order of approximation accuracy is not destroyed by the IRP limiter, provided the cell average stays away from the boundary of the invariant region. Moreover, this limiter is explicit, and easy for computer implementation. A generic algorithm incorporating the IRP limiter is presented for high order finite volume type schemes as long as the evolved cell average of the underlying scheme stays strictly within the invariant region. For any high order discontinuous Galerkin (DG) scheme to the p-system, sufficient conditions are obtained for cell averages to stay in the invariant region. For the viscous p-system, we design both second and third order IRP DG schemes. Numerical experiments are provided to test the proven properties of the IRP limiter and the performance of IRP DG schemes.

Yi Jiang, Hailiang Liu

We introduce an explicit invariant-region-preserving limiter applied to DG methods for compressible Euler equations. The invariant region considered consists of positivity of density and pressure and a maximum principle of a specific entropy. The modified polynomial by the limiter preserves the cell average, lies entirely within the invariant region and does not destroy the high order of accuracy for smooth solutions. Numerical tests are presented to illustrate the properties of the limiter. In particular, the tests on Riemann problems show that the limiter helps to damp the oscillations near discontinuities.

Hailiang Liu, Xuping Tian

In this paper, we propose SGEM, Stochastic Gradient with Energy and Momentum, to solve a large class of general non-convex stochastic optimization problems, based on the AEGD method that originated in the work [AEGD: Adaptive Gradient Descent with Energy. arXiv: 2010.05109]. SGEM incorporates both energy and momentum at the same time so as to inherit their dual advantages. We show that SGEM features an unconditional energy stability property, and derive energy-dependent convergence rates in the general nonconvex stochastic setting, as well as a regret bound in the online convex setting. A lower threshold for the energy variable is also provided. Our experimental results show that SGEM converges faster than AEGD and generalizes better or at least as well as SGDM in training some deep neural networks.

Hailiang Liu, James Ralston, Peimeng Yin

Gaussian beams are asymptotically valid high frequency solutions concentrated on a single curve through the physical domain, and superposition of Gaussian beams provides a powerful tool to generate more general high frequency solutions to PDEs. We present a superposition of Gaussian beams over an arbitrary bounded set of dimension $m$ in phase space, and show that the tools recently developed in [ H. Liu, O. Runborg, and N. M. Tanushev, Math. Comp., 82: 919--952, 2013] can be applied to obtain the propagation error of order $k^{1- \frac{N}{2}- \frac{d-m}{4}}$, where $N$ is the order of beams and $d$ is the spatial dimension. Moreover, we study the sharpness of this estimate in examples.

Tianxiang Gao, Xiaokai Huo, Hailiang Liu et al.

Neural networks with wide layers have attracted significant attention due to their equivalence to Gaussian processes, enabling perfect fitting of training data while maintaining generalization performance, known as benign overfitting. However, existing results mainly focus on shallow or finite-depth networks, necessitating a comprehensive analysis of wide neural networks with infinite-depth layers, such as neural ordinary differential equations (ODEs) and deep equilibrium models (DEQs). In this paper, we specifically investigate the deep equilibrium model (DEQ), an infinite-depth neural network with shared weight matrices across layers. Our analysis reveals that as the width of DEQ layers approaches infinity, it converges to a Gaussian process, establishing what is known as the Neural Network and Gaussian Process (NNGP) correspondence. Remarkably, this convergence holds even when the limits of depth and width are interchanged, which is not observed in typical infinite-depth Multilayer Perceptron (MLP) networks. Furthermore, we demonstrate that the associated Gaussian vector remains non-degenerate for any pairwise distinct input data, ensuring a strictly positive smallest eigenvalue of the corresponding kernel matrix using the NNGP kernel. These findings serve as fundamental elements for studying the training and generalization of DEQs, laying the groundwork for future research in this area.

Hailiang Liu, Yan-Han Chen

In recent years, the machine learning community has increasingly embraced the optimal transport (OT) framework for modeling distributional relationships. In this work, we introduce a sample-based neural solver for computing the Wasserstein geodesic between a source and target distribution, along with the associated velocity field. Building on the dynamical formulation of the optimal transport (OT) problem, we recast the constrained optimization as a minimax problem, using deep neural networks to approximate the relevant functions. This approach not only provides the Wasserstein geodesic but also recovers the OT map, enabling direct sampling from the target distribution. By estimating the OT map, we obtain velocity estimates along particle trajectories, which in turn allow us to learn the full velocity field. The framework is flexible and readily extends to general cost functions, including the commonly used quadratic cost. We demonstrate the effectiveness of our method through experiments on both synthetic and real datasets.

Tianxiang Gao, Siyuan Sun, Hailiang Liu et al.

Neural Ordinary Differential Equations (ODEs) have been successful in various applications due to their continuous nature and parameter-sharing efficiency. However, these unique characteristics also introduce challenges in training, particularly with respect to gradient computation accuracy and convergence analysis. In this paper, we address these challenges by investigating the impact of activation functions. We demonstrate that the properties of activation functions, specifically smoothness and nonlinearity, are critical to the training dynamics. Smooth activation functions guarantee globally unique solutions for both forward and backward ODEs, while sufficient nonlinearity is essential for maintaining the spectral properties of the Neural Tangent Kernel (NTK) during training. Together, these properties enable us to establish the global convergence of Neural ODEs under gradient descent in overparameterized regimes. Our theoretical findings are validated by numerical experiments, which not only support our analysis but also provide practical guidelines for scaling Neural ODEs, potentially leading to faster training and improved performance in real-world applications.

Tianxiang Gao, Hailiang Liu, Jia Liu et al.

Implicit deep learning has received increasing attention recently due to the fact that it generalizes the recursive prediction rules of many commonly used neural network architectures. Its prediction rule is provided implicitly based on the solution of an equilibrium equation. Although a line of recent empirical studies has demonstrated its superior performances, the theoretical understanding of implicit neural networks is limited. In general, the equilibrium equation may not be well-posed during the training. As a result, there is no guarantee that a vanilla (stochastic) gradient descent (SGD) training nonlinear implicit neural networks can converge. This paper fills the gap by analyzing the gradient flow of Rectified Linear Unit (ReLU) activated implicit neural networks. For an $m$-width implicit neural network with ReLU activation and $n$ training samples, we show that a randomly initialized gradient descent converges to a global minimum at a linear rate for the square loss function if the implicit neural network is \textit{over-parameterized}. It is worth noting that, unlike existing works on the convergence of (S)GD on finite-layer over-parameterized neural networks, our convergence results hold for implicit neural networks, where the number of layers is \textit{infinite}.

Hailiang Liu, Xuping Tian

We propose AEGD, a new algorithm for first-order gradient-based optimization of non-convex objective functions, based on a dynamically updated energy variable. The method is shown to be unconditionally energy stable, irrespective of the step size. We prove energy-dependent convergence rates of AEGD for both non-convex and convex objectives, which for a suitably small step size recovers desired convergence rates for the batch gradient descent. We also provide an energy-dependent bound on the stationary convergence of AEGD in the stochastic non-convex setting. The method is straightforward to implement and requires little tuning of hyper-parameters. Experimental results demonstrate that AEGD works well for a large variety of optimization problems: it is robust with respect to initial data, capable of making rapid initial progress. The stochastic AEGD shows comparable and often better generalization performance than SGD with momentum for deep neural networks.

Xiaobing Feng, Hailiang Liu, Shu Ma

In this paper, we propose a family of time-stepping schemes for approximating general nonlinear Schrödinger equations. The proposed schemes all satisfy both mass conservation and energy conservation. Truncation and dispersion error analyses are provided for each proposed scheme. Efficient fixed-point iterative solvers are also constructed to solve the resulting nonlinear discrete problems. As a byproduct, an efficient one-step implementation of the BDF schemes is obtained as well. Extensive numerical experiments are presented to demonstrate the convergence and the capability of capturing the blow-up phenomenon of the proposed schemes.