Midwest Numerical Analysis Day

Invited Talks

Michele Benzi, Emory University
Numerical Linear Algebra Problems Arising in Network Analysis

Network science is a rapidly growing interdisciplinary area at the intersection of mathematics, physics, computer science, and a multitude of disciplines ranging from the life sciences to the social sciences and even the humanities. Network analysis methods are now widely used in biology, the study of social networks (both human and animal), finance, ecology, bibliometric studies, archeology, the evolution of cities, and a host of other fields.

After giving a brief overview of network science, I will discuss some basic mathematical and computational problems arising in network analysis, with a focus on the fundamental notions of centrality, communicability, and robustness. Quantitative, walk-based measures of these notions will be introduced and motivated. I will then show how these measures can be efficiently computed, even for large networks, using state-of-the-art numerical linear algebra techniques. Heuristics for edge manipulation aimed at obtaining robust networks will also be discussed. I will conclude with some open challenges in computational network analysis. Most of the talk is intended to be accessible to a broad audience.

Susanne Brenner, Louisiana State University 
C0 Interior Penalty Methods

C0 interior penalty methods are discontinuous Galerkin methods for fourth order problems that are based on standard Lagrange finite element spaces for second order problems. In this talk we will discuss the apriori and a posteriori error analyses of these methods for fourth order elliptic boundary value problems and elliptic variational inequalities.

Yoichiro Mori, University of Minnesota
Analysis of the Dynamics of Immersed Elastic Filaments in Stokes Flow

Problems in which immersed elastic structures interact with the surrounding fluid abound in science and engineering. Despite their scientific importance, analysis and numerical analysis of such problems are scarce or non-existent. In this talk, we consider the problem of an elastic filament immersed in a 2D or 3D Stokes fluid. We first discuss our recent results on the analysis of the immersed filament problem in a 2D Stokes fluid (the Peskin problem). We prove well-posedness and immediate regularization of the elastic filament configuration, and discuss the implication of these results for numerical analysis. We will then discuss the immersed filament problem in a 3D Stokes fluid (the Slender Body problem). Here, it has not even been clear what the appropriate mathematical formulation of the problem should be. We propose a mathematical formulation for the Slender Body problem and discuss well- posedness for the stationary version of this problem. Furthermore, we prove that the Slender Body approximation, introduced by Keller and Rubinow in the 1980's and used widely in the fluid-structure interaction community, provides an approximation to the Slender Body problem with some error bound. This is joint work with Analise Rodenberg, Laurel Ohm and Dan Spirn.

Zhimin Zhang, Beijing Computational Science Research Center and Wayne State University 
Why Spectral Methods are Preferred in PDE Eigenvalue Computations—in Some Cases?

When approximating PDE eigenvalue problems by numerical methods such as finite difference and finite element, it is common knowledge that only a small portion of numerical eigenvalues are reliable. As a comparison, spectral methods may perform extremely well in some situation, especially for 1-D problems. In addition, we demonstrate that spectral methods can outperform traditional methods and the state-of-the-art method in 2-D problems even with singularities.

Contributed Talks

Mahboub Baccouch, University of Nebraska, Omaha
Asymptotically Exact a Posteriori Error Estimates for the Local Discontinuous Galerkin Method for Nonlinear Convection-Diffusion Problems

In this talk, we present and analyze implicit a posteriori error estimates for the local discontinuous Galerkin (LDG) method applied to nonlinear convection-diffusion problems in one space dimension. Optimal a priori error estimates for the solution and for the auxiliary variable that approximates the first-order derivative are derived in the L2-norm for the semi- discrete formulation. More precisely, we identify special numerical fluxes and a suitable projection of the initial condition for the LDG scheme to achieve p+1 order of convergence for the solution and its spatial derivative in the L2-norm, when piecewise polynomials of degree at most p are used. We further prove that the LDG solution is superconvergent with order p+3/2 towards a special Gauss-Radau projection of the exact solution. Our superconvergence results allow us to show that the leading error term on each element is proportional to the (p+1)-degree right Radau polynomial. We use these results to construct asymptotically exact a posteriori error estimator. The proposed a posteriori LDG error estimate converges at a fixed time to the true spatial error in the L2-norm at O(h{p+3/2}) rate. Finally, several numerical examples are given to validate the theoretical results.

Majid Bani-Yaghoub, University of Missouri-Kansas City
Numerical Simulations of a Delayed Hyperbolic–Parabolic Population Model

We numerically explore a class of nonlocal delayed hyperbolic–parabolic model. The existence and uniqueness of a traveling wavefront of the model has already been established. In the present work, we numerically show that the solution of the model may converge to a non-monotonic traveling wavefront. This is a joint work with Dr. Chunhua Ou ( Memorial University of Newfoundland) and Dr. Guangming Yao (Clarkson University).

Paul Cazeaux, University of Kansas
Relaxation and Moiré Domain Formation in Incommensurate 2D Heterostructures

The recent discovery of a whole family of two-dimensional crystalline materials such as graphene, hexagonal boron nitride (h-BN) and many others leads to study the properties of their combinations, particularly by stacking a number of layers vertically. Such structures are generally non-periodic, with interesting geometric properties such as Moiré effects which have been observed experimentally as highly regular domain-wall patterns. In this talk, we will discuss the origin of this relaxation and study a novel model based on an alternative description in configuration space, bypassing the need for supercell approximations. The numerical discretization of this model leads to efficient simulations that give new insights as to the elastic mechanisms at play in the Moiré patterns.

Gang Chen, Missouri University of Science and Technology
Analysis of Hybridizable Discontinuous Galerkin Finite Element Method for Time-Harmonic Maxwell Equations

We consider the time-harmonic Maxwell equations with Dirichlet boundary condition and zero frequency. The formulation of the Maxwell equations we consider here is a mixed curl-curl formulation, where the divergence condition is imposed by introducing a Lagrange multiplier. In the first part, we generalize regularity of Maxwell equations. In the second part, we introduce and analyze a hybridizable discontinuous Galerkin (HDG) method for the time- harmonic Maxwell equations. All the unknowns, as long as the unknowns at inter-element faces, of the underlying system are approximated by discontinuous polynomials. For the low regularity case, special care has been done to approach the data on boundary. This method is shown be be stable and optimal convergent, under high regularity and minimal regularity, respectively. Numerical experiments with both smooth and singular analytical solutions are presented.

Vani Cheruvu, University of Toledo
Spectra of Boundary Integral Operators Defined on the Unit Sphere for the Modified Laplace Equation

We consider a modified Laplace equation on a unit sphere. Spherical harmonics are used for the expansion of the unknown function. We show that on the unit sphere, both modified Laplace single and double layer operators diagonalize in spherical harmonic basis. The analytic expressions for evaluating the operators away from the boundary are also derived. Currently, we are working on the numerical aspects. In this talk, we present both the analytical and numerical results of our work.

Hailiang Liu, Iowa State University
Invariant-Region-Preserving DG Methods for Multi-Dimensional Hyperbolic Conservation Law Systems, with an Application to Compressible Euler Equations

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. This is a joint work with Yi Jiang (Iowa State University).

Suzanne Shontz, University of Kansas
A Novel Algorithm for High-Order Triangular and Tetrahedral Mesh Generation

In this talk, we present our method for generating second and third-order Lagrange triangular and tetrahedral meshes based on affine combinations of nodal positions. Our mesh generation approach consists of the following steps. First, for each interior mesh node, an optimization problem is solved to calculate an affine combination of nodal positions that relate each interior node to its neighbors. Second, a deformation is applied to the high-order boundary nodes to move them onto the true boundary. Third, new positions for the interior nodes are calculated by solving a linear system of equations using the weights and the new boundary node positions from the first two steps, respectively. To reconstruct each of the true boundaries in the second step, a Radial Basis Function (RBF) is fit to give an implicit representation of each boundary. The boundary deformation then consists of moving each high-order boundary node to the position that minimizes the squared value of the implicit function; this minimum occurs when the node reaches the true boundary. We present several numerical examples in two and three dimensions that demonstrate the viability of our method, as well as discuss the distortion of the elements in the resulting meshes. We will also compare this method to our previous high-order mesh generation method that uses convex combinations of nodal positions to relate each interior node to its neighbors.

Jimmy Vogel, Purdue University
Divide-and-Conquer Eigensolvers for Computational Physics

In this talk I illustrate how hierarchical semiseparable (HSS) matrix structure can be used to efficiently and accurately solve several eigenvalue problems from active research areas in computational physics. The talk will include some basic background to HSS matrices and relevant physics, algorithm descriptions, discussion of various numerical stability and efficiency issues, and some representative test results. This is joint work with my advisor Prof. Jianlin Xia.

Lei Wang, University of Wisconsin-Milwaukee 
A Fast Treecode Algorithm for 3D Stokes Flow

A large number of problems in fluid dynamics are modeled as many-particle interactions in Stokes flows, for example, simulations of falling jets of particles in viscous fluids, microfluidic crystals, and vesicle flows. The formulation is often based on fundamental solutions. The Stokeslet and the Stresslet are the kernels in the single and double layer potentials, respectively. Many situations (e.g., through superposition or discretization of boundary integrals) involve sums of Stokeslets and Stresslets, which is an example of an N-body problem and the direct sum requires O(N2) operations. This can make the numerical calculation prohibitively expensive. A Barnes-Hut tree treecode algorithm is developed for speeding up the computation. The particles are restructured recursively into a tree, and the particle-particle interactions are replaced with particle-cluster interactions computed by either a far-field expansion or a direct summation. Numerical results exhibit the promising performance of the algorithm both in serial and parallel machine.

Jue Yan, Iowa State University
Recent Developments on Direct Discontinuous Galerkin Methods

We first introduce the direct discontinuous Galerkin (DDG) method and its variations, namely the DDGIC and symmetric DDG methods. Compared to the leading diffusion DG method solvers like the interior penalty method (SIPG), we find out our diffusion solver the DDG methods have many advantages. Under the topic of maximum principle, DDG methods numerical solution can be proved to satisfy strict maximum principle even on unstructured triangular meshes with at least third order of accuracy.

Recently we develop DDG methods to solve Keller-Segel Chemotaxis equations. Different to available DG methods or other numerical methods in literature, we introduce no extra variable to approximate the chemical density gradients and we solve the system directly. With Pk polynomial approximations, we observe no order loss for the density variable. The reason behind is that the DDGIC or symmetric DDG methods have the hidden super convergence property on its approximation to the solution gradients. With Fourier (Von Neumann) analysis technique, we prove the DDG solution’s spatial derivative is super convergent with at least (k+1)th order under moment format or in the weak sense. Notice that we do not have super convergence for the SIPG method. We show the cell density approximations are strictly positive with at least third order of accuracy.

He Yang, Augusta University
Local Discontinuous Galerkin Methods for the Khokhlov–Zabolotskaya–Kuznetzov Equation

Khokhlov–Zabolotskaya–Kuznetzov (KZK) equation is a model that describes the propagation of the ultrasound beams in the thermoviscous fluid. It contains a nonlocal diffraction term, an absorption term and a nonlinear term. Accurate numerical methods to simulate the KZK equation are important to its broad applications in medical ultrasound simulations. In this paper, we propose a local discontinuous Galerkin method to solve the KZK equation. We prove the L2 stability of our scheme and conduct a series of numerical experiments including the focused circular short tone burst excitation and the propagation of unfocused sound beams, which show that our scheme leads to accurate solutions and performs better than the benchmark solutions in the literature.

Yufei Yu, University of Kansas
Effects of Permanent Charge and Boundary Condition on Ionic Flow via a Quasi-1D Poisson- Nernst-Planck Model

Ionic channels—large proteins on cell membrane—are a major way for ions to transport through cell membrane that carries electric signals for cells to communicate with each other. The permanent charge of an ion channel is the crucial structure for ionic flow properties of the channel. The effects of permanent charges interacting with boundary conditions have been studied analytically via the Quasi-1D Poisson-Nernst-Planck (PNP) model for small permanent charge and for large permanent charge. In this talk, we will present results of numerical investigation to bridge between the two extrema. As expected, our numerical results verify the analytical predictions for small and large permanent charges. On the other hand, non-trivial behavior emerges as one varies the permanent charge from small to large, in particular, bifurcations are revealed, showing the rich phenomena of permanent charge effects by the power of combining the analytical and numerical studies. An adaptive moving mesh finite element method has been applied which is critical due to the presence of Debye layers at the interface between the permanent charge regions and uncharged regions of ion channels. This talk is based on a joint work with W. Huang and W. Liu, both from KU.

Yangwen Zhang, Missouri University of Science and Technology
An HDG Method for Tangential Boundary Control of Stokes Equations

We propose a hybridizable discontinuous Galerkin(HDG) method to approximate the solution of a tangential Dirichlet boundary control problem governed by Stokes equation. Dirichlet boundary control problems are naturally studied in the L2-setting, and therefore well-known to be challenging due to its non variational form and low regularity, even for simpler elliptic PDEs. Although there are many works in the literature on Dirichlet boundary control problems for fluids flow, the authors are not aware of any existing theoretical or numerical analysis for Dirichlet boundary control of Stokes or Navier-Stokes equation in such setting. To avoid the extreme low regularity of solutions for Stokes Dirichlet boundary control, this work gives a first attempt to the analysis of tangential boundary control of fluid flows with controls in L2. The contribution of this paper is twofold. First, we obtain well-posedness and regularity results for the tangential Dirichlet control problem. Second, under certain assumptions on the domain and the target state, we obtain optimal a priori error estimates for the control in 2D for the HDG method. We present numerical experiments to demonstrate the performance of the HDG method.