A Fourth-Order Upwinding Embedded Boundary Method (UEBM) for Maxwell's Equations in Media with Material Interfaces: Part I

In this paper, we present a new fourth-order upwinding embedded boundary method (UEBM) over Cartesian grids, originally proposed in the Journal of Computational Physics [190 (2003), pp. 159-183.] as a second-order method for treating material interfaces for Maxwell's equations. In addition to the idea of the UEBM to evolve solutions at interfaces, we utilize the ghost fluid method to construct finite difference approximation of spatial derivatives at Cartesian grid points near the material interfaces. As a result, Runge-Kutta type time discretization can be used for the semidiscretized system to yield an overall fourth-order method, in contrast to the original second-order UEBM based on a Lax-Wendroff type difference. The final scheme allows time step sizes independent of the interface locations. Numerical examples are given to demonstrate the fourth-order accuracy as well as the stability of the method. We tested the scheme for several wave problems with various material interface locations, including electromagnetic scattering of a plane wave incident on a planar boundary and a two-dimensional electromagnetic application with an interface parallel to the y-axis.     

A Kernel-Free Boundary Integral Method for Variable Coefficients Elliptic PDEs

This work proposes a generalized boundary integral method for variable coefficients elliptic partial differential equations (PDEs), including both boundary value and interface problems. The method is kernel-free in the sense that there is no need to know analytical expressions for kernels of the boundary and volume integrals in the solution of boundary integral equations. Evaluation of a boundary or volume integral is replaced with interpolation of a Cartesian grid based solution, which satisfies an equivalent discrete interface problem, while the interface problem is solved by a fast solver in the Cartesian grid. The computational work involved with the generalized boundary integral method is essentially linearly proportional to the number of grid nodes in the domain. This paper gives implementation details for a second-order version of the kernel-free boundary integral method in two space dimensions and presents numerical experiments to demonstrate the efficiency and accuracy of the method for both boundary value and interface problems. The interface problems demonstrated include those with piecewise constant and large-ratio coefficients and the heterogeneous interface problem, where the elliptic PDEs on two sides of the interface are of different types.     

A Bilinear Immersed Finite Volume Element Method for the Diffusion Equation with Discontinuous Coefficient

This paper is to present a finite volume element (FVE) method based on the bilinear immersed finite element (IFE) for solving the boundary value problems of the diffusion equation with a discontinuous coefficient (interface problem). This method possesses the usual FVE method's local conservation property and can use a structured mesh or even the Cartesian mesh to solve a boundary value problem whose coefficient has discontinuity along piecewise smooth nontrivial curves. Numerical examples are provided to demonstrate features of this method. In particular, this method can produce a numerical solution to an interface problem with the usual O(h²) (in L² norm) and O(h) (in H¹ norm) convergence rates.     

Approximate Boundary Conditions for Patch Antennas Mounted on Thin Dielectric Layers

In this paper we discuss scattering problems inherent in curved microstrip structures mounted on thin dielectric structures. We provide approximate boundary conditions for such structures in the framework of integral equations.     

A Highly Scalable Boundary Integral Equation and Walk-on-Spheres (BIE-WOS) Method for the Laplace Equation with Dirichlet Data

In this paper, we study a highly scalable communication-free parallel domain boundary decomposition algorithm for the Laplace equation based on a hybrid method combining boundary integral equations and walk-on-spheres (BIE-WOS) method, which provides a numerical approximation of the Dirichlet-to-Neumann (DtN) mapping for the Laplace equation. The BIE-WOS is a local method on the boundary of the domain and does not require a structured mesh, and only needs a covering of the domain boundary by patches and a local mesh for each patch for a local BIE. A new version of the BIE-WOS method with second kind integral equations is introduced for better error controls. The effect of errors from the Feynman-Kac formula based path integral WOS method on the overall accuracy of the BIE-WOS method is analyzed for the BIEs, especially in the calculation of the right hand sides of the BIEs. For the special case of flat patches, it is shown that the second kind integral equation of BIE-WOS method can be simplified where the local BIE solutions can be given in closed forms. A key advantage of the parallel BIE-WOS method is the absence of communications during the computation of the DtN mapping on individual patches of the boundary, resulting in a complete independent computation using a large number of cluster nodes. In addition, the BIE-WOS has an intrinsic capability of fault tolerance for exascale computations. The nearly linear scalability of the parallel BIE-WOS method on a large-scale cluster with 6400 CPU cores is verified for computing the DtN mapping of exterior Laplace problems with Dirichlet data for several domains.     

An Augmented Lagrangian Deep Learning Method for Variational Problems with Essential Boundary Conditions

This paper is concerned with a novel deep learning method for variational problems with essential boundary conditions. To this end, we first reformulate the original problem into a minimax problem corresponding to a feasible augmented Lagrangian, which can be solved by the augmented Lagrangian method in an infinite dimensional setting. Based on this, by expressing the primal and dual variables with two individual deep neural network functions, we present an augmented Lagrangian deep learning method for which the parameters are trained by the stochastic optimization method together with a projection technique. Compared to the traditional penalty method, the new method admits two main advantages: i) the choice of the penalty parameter is flexible and robust, and ii) the numerical solution is more accurate in the same magnitude of computational cost. As typical applications, we apply the new approach to solve elliptic problems and (nonlinear) eigenvalue problems with essential boundary conditions, and numerical experiments are presented to show the effectiveness of the new method.     

Two Uniform Tailored Finite Point Schemes for the Two Dimensional Discrete Ordinates Transport Equations with Boundary and Interface Layers

This paper presents two uniformly convergent numerical schemes for the two dimensional steady state discrete ordinates transport equation in the diffusive regime, which is valid up to the boundary and interface layers. A five-point node-centered and a four-point cell-centered tailored finite point schemes (TFPS) are introduced. The schemes first approximate the scattering coefficients and sources by piecewise constant functions and then use special solutions to the constant coefficient equation as local basis functions to formulate a discrete linear system. Numerically, both methods can not only capture the diffusion limit, but also exhibit uniform convergence in the diffusive regime, even with boundary layers. Numerical results show that the five-point scheme has first-order accuracy and the four-point scheme has second-order accuracy, uniformly with respect to the mean free path. Therefore, a relatively coarse grid can be used to capture the two dimensional boundary and interface layers.     

Integral Equation Method for a Non-Selfadjoint Steklov Eigenvalue Problem

We propose a numerical method for a non-selfadjoint Steklov eigenvalue problem of the Helmholtz equation. The problem is formulated using boundary integrals. The Nyström method is employed to discretize the integral operators, which leads to a non-Hermitian generalized matrix eigenvalue problems. The spectral indicator method (SIM) is then applied to calculate the (complex) eigenvalues. The convergence is proved using the spectral approximation theory for (non-selfadjoint) compact operators. Numerical examples are presented for validation.     

A High-Order Accurate Gas-Kinetic Scheme for One- and Two-Dimensional Flow Simulation

This paper develops a high-order accurate gas-kinetic scheme in the framework of the finite volume method for the one- and two-dimensional flow simulations, which is an extension of the third-order accurate gas-kinetic scheme [Q.B. Li, K. Xu, and S. Fu, J. Comput. Phys., 229(2010), 6715-6731] and the second-order accurate gas-kinetic scheme [K. Xu, J. Comput. Phys., 171(2001), 289-335]. It is formed by two parts: quartic polynomial reconstruction of the macroscopic variables and fourth-order accurate flux evolution. The first part reconstructs a piecewise cell-center based quartic polynomial and a cell-vertex based quartic polynomial according to the "initial" cell average approximation of macroscopic variables to recover locally the non-equilibrium and equilibrium single particle velocity distribution functions around the cell interface. It is in view of the fact that all macroscopic variables become moments of a single particle velocity distribution function in the gas-kinetic theory. The generalized moment limiter is employed there to suppress the possible numerical oscillation. In the second part, the macroscopic flux at the cell interface is evolved in fourth-order accuracy by means of the simple particle transport mechanism in the microscopic level, i.e. free transport and the Bhatnagar-Gross-Krook (BGK) collisions. In other words, the fourth-order flux evolution is based on the solution (i.e. the particle velocity distribution function) of the BGK model for the Boltzmann equation. Several 1D and 2D test problems are numerically solved by using the proposed high-order accurate gas-kinetic scheme. By comparing with the exact solutions or the numerical solutions obtained the second-order or third-order accurate gas-kinetic scheme, the computations demonstrate that our scheme is effective and accurate for simulating invisid and viscous fluid flows, and the accuracy of the high-order GKS depends on the choice of the (numerical) collision time.     

On Integral Equation and Least Squares Methods for Scattering by Diffraction Gratings

In this paper we consider the scattering of a plane acoustic or electromagnetic wave by a one-dimensional, periodic rough surface. We restrict the discussion to the case when the boundary is sound soft in the acoustic case, perfectly reflecting with TE polarization in the EM case, so that the total field vanishes on the boundary. We propose a uniquely solvable first kind integral equation formulation of the problem, which amounts to a requirement that the normal derivative of the Green's representation formula for the total field vanish on a horizontal line below the scattering surface. We then discuss the numerical solution by Galerkin's method of this (ill-posed) integral equation. We point out that, with two particular choices of the trial and test spaces, we recover the so-called SC (spectral-coordinate) and SS (spectral-spectral) numerical schemes of DeSanto et al., Waves Random Media, 8, 315-414, 1998. We next propose a new Galerkin scheme, a modification of the SS method that we term the SS^∗ method, which is an instance of the well-known dual least squares Galerkin method. We show that the SS^∗ method is always well-defined and is optimally convergent as the size of the approximation space increases. Moreover, we make a connection with the classical least squares method, in which the coefficients in the Rayleigh expansion of the solution are determined by enforcing the boundary condition in a least squares sense, pointing out that the linear system to be solved in the SS^∗ method is identical to that in the least squares method. Using this connection we show that (reflecting the ill-posed nature of the integral equation solved) the condition number of the linear system in the SS^∗ and least squares methods approaches infinity as the approximation space increases in size. We also provide theoretical error bounds on the condition number and on the errors induced in the numerical solution computed as a result of ill-conditioning. Numerical results confirm the convergence of the SS^∗ method and illustrate the ill-conditioning that arises.     

A Novel Iterative Penalty Method to Enforce Boundary Conditions in Finite Volume POD-Galerkin Reduced Order Models for Fluid Dynamics Problems

A Finite-Volume based POD-Galerkin reduced order model is developed for fluid dynamics problems where the (time-dependent) boundary conditions are controlled using two different boundary control strategies: the lifting function method, whose aim is to obtain homogeneous basis functions for the reduced basis space and the penalty method where the boundary conditions are enforced in the reduced order model using a penalty factor. The penalty method is improved by using an iterative solver for the determination of the penalty factor rather than tuning the factor with a sensitivity analysis or numerical experimentation.The boundary control methods are compared and tested for two cases: the classical lid driven cavity benchmark problem and a Y-junction flow case with two inlet channels and one outlet channel. The results show that the boundaries of the reduced order model can be controlled with the boundary control methods and the same order of accuracy is achieved for the velocity and pressure fields. Finally, the reduced order models are 270-308 times faster than the full order models for the lid driven cavity test case and 13-24 times for the Y-junction test case.     

Numerical Study on Viscous Fingering Using Electric Fields in a Hele-Shaw Cell

We investigate the nonlinear dynamics of a moving interface in a Hele-Shaw cell subject to an in-plane applied electric field. We develop a spectrally accurate numerical method for solving a coupled integral equation system. Although the stiffness due to the high order spatial derivatives can be removed using a small scale decomposition technique, the long-time simulation is still expensive since the evolving velocity of the interface drops dramatically as the interface expands. We remove this physically imposed stiffness by employing a rescaling scheme, which accelerates the slow dynamics and reduces the computational cost. Our nonlinear results reveal that positive currents restrain finger ramification and promote the overall stabilization of patterns. On the other hand, negative currents make the interface more unstable and lead to the formation of thin tail structures connecting the fingers and a small inner region. When no fluid is injected, and a negative current is utilized, the interface tends to approach the origin and break up into several drops. We investigate the temporal evolution of the smallest distance between the interface and the origin and find that it obeys an algebraic law $(t_∗−t)^b,$ where $t_∗$ is the estimated pinch-off time.     

Numerical investigation of distributed‐order fractional optimal control problems via Bernstein wavelets

The aim of this article is to investigate an efficient computational method for solving distributed‐order fractional optimal control problems. In the proposed method, a new Riemann‐Liouville fractional integral operator for the Bernstein wavelet is given. This approach is based on a combination of the Bernstein wavelets basis, fractional integral operator, Gauss‐Legendre numerical integration, and Newton's method for solving obtained system. Easy implementation, simple operations, and accurate solutions are the essential features of the proposed method. The error analysis of the proposed method is carried out. Examples reveal the applicability of the proposed technique.     

Hybrid Diffuse and Sharp Interface Immersed Boundary Methods for Particulate Flows in the Presence of Complex Boundaries

A coupling framework that leverages the advantages of the diffuse and sharp interface immersed boundary (IB) methods is presented for handling the interaction among particles and particles with the static complex geometries of the environment. In the proposed coupling approach, the curvilinear IB method is employed to represent the static complex geometries, a variant of the direct forcing IB method is proposed for simulating particles, and the discrete element method is employed for particle-particle and particle-wall collisions. The proposed approach is validated using several classical benchmark problems, which include flow around a sphere, sedimentation of a sphere, collision of two sedimenting spheres, and collision between a particle and a flat wall, with the present predictions showing an overall good agreement with the results reported in the literature. The capability of the proposed framework is further demonstrated by simulating the interaction between multiple particles and a wall-mounted cylinder, and the particle-laden turbulent flow over periodic hills. The proposed method provides an efficient way to simulate particle-laden turbulent flows in environments with complex boundaries.     

Weak Galerkin Method for Second-Order Elliptic Equations with Newton Boundary Condition

The weak Galerkin (WG) method is a nonconforming numerical method for solving partial differential equations. In this paper, we introduce the WG method for elliptic equations with Newton boundary condition in bounded domains. The Newton boundary condition is a nonlinear boundary condition arising from science and engineering applications. We prove the well-posedness of the WG scheme by the monotone operator theory and the embedding inequality of weak finite element functions. The error estimates are derived. Numerical experiments are presented to verify the theoretical analysis.     

A Fast Accurate Boundary Integral Method for Potentials on Closely Packed Cells

Boundary integral methods are naturally suited for the computation of harmonic functions on a region having inclusions or cells with different material properties. However, accuracy deteriorates when the cell boundaries are close to each other. We present a boundary integral method in two dimensions which is specially designed to maintain second order accuracy even if boundaries are arbitrarily close. The method uses a regularization of the integral kernel which admits analytically determined corrections to maintain accuracy. For boundaries with many components we use the fast multipole method for efficient summation. We compute electric potentials on a domain with cells whose conductivity differs from that of the surrounding medium. We first solve an integral equation for a source term on the cell interfaces and then find values of the potential near the interfaces via integrals. Finally we use a Poisson solver to extend the potential to a regular grid covering the entire region. A number of examples are presented. We demonstrate that increased refinement is not needed to maintain accuracy as interfaces become very close.     

Three Discontinuous Galerkin Methods for One- and Two-Dimensional Nonlinear Dirac Equations with a Scalar Self-Interaction

This paper develops three high-order accurate discontinuous Galerkin (DG) methods for the one-dimensional (1D) and two-dimensional (2D) nonlinear Dirac (NLD) equations with a general scalar self-interaction. They are the Runge-Kutta DG (RKDG) method and the DG methods with the one-stage fourth-order Lax-Wendroff type time discretization (LWDG) and the two-stage fourth-order accurate time discretization (TSDG). The RKDG method uses the spatial DG approximation to discretize the NLD equations and then utilize the explicit multistage high-order Runge-Kutta time discretization for the first-order time derivatives, while the LWDG and TSDG methods, on the contrary, first give the one-stage fourth-order Lax-Wendroff type and the two-stage fourth-order time discretizations of the NLD equations, respectively, and then discretize the first- and higher-order spatial derivatives by using the spatial DG approximation. The $L^2$ stability of the 2D semi-discrete DG approximation is proved in the RKDG methods for a general triangulation, and the computational complexities of three 1D DG methods are estimated. Numerical experiments are conducted to validate the accuracy and the conservation properties of the proposed methods. The interactions of the solitary waves, the standing and travelling waves are investigated numerically and the 2D breathing pattern is observed.     

On a Fast Integral Equation Method for Diffraction Gratings

The integral equation method for the simulation of the diffraction by optical gratings is an efficient numerical tool if profile gratings determined by simple cross-section curves are considered. This method in its recent version is capable to tackle profile curves with corners, gratings with thin coated layers, and diffraction scenarios with unfavorably large ratio period over wavelength. We discuss special implementational issues including the efficient evaluation of the quasi-periodic Green kernels, the quadrature algorithm, and the iterative solution of the arising systems of linear equations. Finally, as an example we present the simulation of echelle gratings which demonstrates the efficiency of our approach.     

A High-Accurate Fast Poisson Solver Based on Harmonic Surface Mapping Algorithm

Poisson's equations in a cuboid are frequently solved in many scientific and engineering applications such as electric structure calculations, molecular dynamics simulations and computational astrophysics. In this paper, a fast and highly accurate algorithm is presented for the solution of the Poisson's equation in a cuboidal domain with boundary conditions of mixed type. This so-called harmonic surface mapping algorithm is a meshless algorithm which can achieve a desired order of accuracy by evaluating a body convolution of the source and the free-space Green's function within a sphere containing the cuboid, and another surface integration over the spherical surface. Numerical quadratures are introduced to approximate the integrals, resulting in the solution represented by a summation of point sources in free space, which can be accelerated by means of the fast multipole algorithm. The complexity of the algorithm is linear to the number of quadrature points, and the convergence rate can be arbitrarily high even when the source term is a piecewise continuous function.     

A modified pseudospectral scheme for accurate solution of Bang‐Bang optimal control problems

In the present contribution, a modified Legendre pseudospectral scheme for accurate and efficient solution of bang‐bang optimal control problems is investigated. In this scheme control and state functions are considered as piecewise constant and piecewise continuous polynomials, respectively, and the switching points are also taken as decision variables. Furthermore, for simplicity in discretization, the integral formulation of the dynamical equations is considered. Thereby, the problem is converted into a mathematical programming problem which can be solved by the well‐developed parameter optimization algorithms. The main advantages of the present method are that: (i) it obtains good results even by using a small number of collocation points and the rate of convergence is high; (ii) the switching times can be captured accurately; and (iii) the wrongly chosen number of switching points can be detected by the results of the method. These are illustrated through a numerical implementation of the method on three examples and the efficiency of the method is reported. Copyright © 2010 John Wiley & Sons, Ltd.