Special Issue Dedicated to Professor Jie Shen's 60th Birthday

In this paper, an efficient spectral Petrov-Galerkin time-stepping method for solving nonlinear Hamiltonian systems is presented and studied. Conservation properties of the proposed method (including symplectic structure preserving and energy conservation) are discussed. Iterative algorithm on how to discretize the nonlinear term is introduced and the uniqueness, stability and convergence properties of the iterative algorithm are also established. Finally, numerical experiments are presented to verify the efficiency of our algorithm.

This paper is concerned with the acoustic scattering of a point incident wave by a sound hard obstacle embedded in a two-layered lossy background medium which is separated by an infinite rough surface. Given the point incident wave, the direct scattering problem is to determine the acoustic wave field for the given obstacle and infinite rough surface; the inverse scattering problem is to determine both the obstacle and the infinite rough surface from the reflected and transmitted wave fields measured on two plane surfaces enclosing the structure. For the direct scattering problem, the well-posedness is studied by using the method of boundary integral equations. For the inverse scattering problem, we prove that the obstacle and the infinite rough surface can be uniquely determined by the measured wave fields corresponding to a single point incident wave. To prove the local stability, the domain derivative of the wave field with respect to the change of the shapes of the obstacle and the infinite rough surface is examined. The local stability indicates that the Hausdorff distance of two domains is bounded above by the distance of corresponding wave fields if the two domains are close enough.

We establish error bounds of the finite difference time domain (FDTD) methods for the long time dynamics of the nonlinear Klein-Gordon equation (NKGE) with a cubic nonlinearity, while the nonlinearity strength is characterized by ε² with 0<ε≤1 a dimensionless parameter. When 0<ε≪1, it is in the weak nonlinearity regime and the problem is equivalent to the NKGE with small initial data, while the amplitude of the initial data (and the solution) is at O(ε). Four different FDTD methods are adapted to discretize the problem and rigorous error bounds of the FDTD methods are established for the long time dynamics, i.e. error bounds are valid up to the time at O(1/ε^β) with 0 ≤ β ≤ 2, by using the energy method and the techniques of either the cut-off of the nonlinearity or the mathematical induction to bound the numerical approximate solutions. In the error bounds, we pay particular attention to how error bounds depend explicitly on the mesh size h and time step τ as well as the small parameter ε∈(0,1], especially in the weak nonlinearity regime when 0<ε≪1. Our error bounds indicate that, in order to get "correct" numerical solutions up to the time at O(1/ε^β), the ε-scalability (or meshing strategy) of the FDTD methods should be taken as: h=O(ε^β/2) and τ=O(ε^β/2). As a by-product, our results can indicate error bounds and ε-scalability of the FDTD methods for the discretization of an oscillatory NKGE which is obtained from the case of weak nonlinearity by a rescaling in time, while its solution propagates waves with wavelength at O(1) in space and O(ε^β) in time. Extensive numerical results are reported to confirm our error bounds and to demonstrate that they are sharp.

In this paper we propose and analyze an energy stable numerical scheme for the square phase field crystal (SPFC) equation, a gradient flow modeling crystal dynamics at the atomic scale in space but on diffusive scales in time. In particular, a modification of the free energy potential to the standard phase field crystal model leads to a composition of the 4-Laplacian and the regular Laplacian operators. To overcome the difficulties associated with this highly nonlinear operator, we design numerical algorithms based on the structures of the individual energy terms. A Fourier pseudo-spectral approximation is taken in space, in such a way that the energy structure is respected, and summation-by-parts formulae enable us to study the discrete energy stability for such a high-order spatial discretization. In the temporal approximation, a second order BDF stencil is applied, combined with an appropriate extrapolation for the concave diffusion term(s). A second order artificial Douglas-Dupont-type regularization term is added to ensure energy stability, and a careful analysis leads to the artificial linear diffusion coming at an order lower than that of surface diffusion term. Such a choice leads to reduced numerical dissipation. At a theoretical level, the unique solvability, energy stability are established, and an optimal rate convergence analysis is derived in the $ℓ$^∞(0,$T$;$ℓ$²)∩$ℓ$²(0,$T$;$H_N^3$) norm. In the numerical implementation, the preconditioned steepest descent (PSD) iteration is applied to solve for the composition of the highly nonlinear 4-Laplacian term and the standard Laplacian term, and a geometric convergence is assured for such an iteration. Finally, a few numerical experiments are presented, which confirm the robustness and accuracy of the proposed scheme.

n this paper, we propose a family of time-stepping schemes for approximating general nonlinear Schrödinger equations. The proposed schemes all satisfy both mass and energy conservation (in a modified form for the latter). Truncation and dispersion error analyses are provided for four proposed schemes. 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 time of the proposed schemes.

We consider a midpoint scheme to approximate analytical solutions to a white noise driven BBM equation that reads du−du_xx+u_x ◦dW+u^pu_xdt=0. We prove the well-posedness of the time-discrete approximation scheme and we provide the strong error order, which is 1.

In this paper, we propose and analyze a Müntz spectral method for a class of two-dimensional space-fractional convection-diffusion equations. The proposed methods make new use of the fractional polynomials, also known as Müntz polynomials, which can be regarded as continuation of our previous work. The extension is twofold. Firstly, the existing Müntz spectral method for fractional differential equation with fractional derivative order 0<µ<1 is generalized to 0<µ≤2, which is nontrivial since the classical Müntz polynomials only have no more than H^1/2+λ regularity, where 0 < λ ≤ 1 is the characteristic parameter of the Müntz polynomial. Secondly, 1D Müntz spectral method is extended to the 2D space-fractional convection-diffusion equation. Compared to the time-fractional diffusion equation, some new operators such as suitable H¹-projectors are needed to analyze the error of the numerical solution. The main contribution of the present paper consists of an efficient method combining the Crank-Nicolson scheme for the temporal discretization and a new spectral method using the Müntz Jacobi polynomials for the spatial discretization of the 2D space-fractional convection-diffusion equation. A detailed convergence analysis is carried out, and several error estimates are established. Finally, a series of numerical experiments is performed to verify the theoretical claims.

In this paper, we propose highly efficient, unconditionally energy-stable numerical schemes to approximate the isotropic phase field model of solid-state dewetting problems by using the invariant energy quadratization (IEQ) method. The phase field model is governed by the isotropic Cahn-Hilliard equation with degenerate mobilities and dynamic contact line boundary conditions. By using the backward differential formula to discretize temporal derivatives, we construct linearly first- and second-order IEQ schemes for solving the model. It can be rigorously proved that these numerical schemes are unconditionally energy-stable and satisfy the total mass conservation during the evolution. By performing numerical simulations, we demonstrate that these IEQ-based schemes (including the first-order IEQ/BDF1, second-order IEQ/BDF2) are highly efficient, accurate and energy-stable. Furthermore, many interesting dewetting phenomena (such as the hole dynamics, pinch-off), are investigated by using the proposed IEQ schemes.

This paper is devoted to finite element analysis for the magneto-heat coupling model which governs the electromagnetic fields in large power transformers. The model, which couples Maxwell's equations and Heat equation through Ohmic heat source, is nonlinear. First we derive an equivalent weak formulation for the nonlinear magneto-heat model. We propose a linearized and temporally discrete scheme to approximate the continuous problem. The well-posedness and error estimates are proven for the semi-discrete scheme. Based on the results, we propose a fully discrete finite element problem and prove the error estimates for the approximate solutions. To validate the magneto-heat model and verify the efficiency of the finite element method, we compute an engineering benchmark problem of the International Compumag Society, P21^b-MN. The numerical results agree well with experimental data.

In this paper, we develop the energy law preserving method for a phase-field model of Cahn-Hilliard type describing binary mixtures. A new class of dynamic boundary conditions in a rather general setting proposed in [1] is adopted here. The model equations are discretized by a continuous finite element method in space and a midpoint scheme in time. The discrete energy law of the numerical method for the model with the dynamic boundary conditions is derived. By a few two-phase examples, we demonstrate the performance of the energy law preserving method for the computation of the phase-field model with the new class of dynamic boundary conditions, even in the case of relatively coarse mesh.

In this paper, we rigorously prove the convergence of fully discrete first- and second-order exponential time differencing schemes for solving the Cahn-Hilliard equation. Our analyses mainly follow the standard procedure with the consistency and stability estimates for numerical error functions, while the technique of higher-order consistency analysis is adopted in order to obtain the uniform L^∞ boundedness of the numerical solutions under some moderate constraints on the time step and spatial mesh sizes. This paper provides a theoretical support for numerical analysis of exponential time differencing and other related numerical methods for phase field models, in which an assumption on the uniform L^∞ boundedness is usually needed.

In this paper, a high order arbitrary Lagrangian-Eulerian (ALE) finite difference weighted essentially non-oscillatory (WENO) method for Hamilton-Jacobi equations is developed. This method is based on moving quadrilateral meshes, which are often used in Lagrangian type methods. The algorithm is formed in two parts: spatial discretization and temporal discretization. In the spatial discretization, we choose a new type of multi-resolution WENO schemes on a nonuniform moving mesh. In the temporal discretization, we use a strong stability preserving (SSP) Runge-Kutta method on a moving mesh for which each grid point moves independently, with guaranteed high order accuracy under very mild smoothness requirement (Lipschitz continuity) for the mesh movements. Extensive numerical tests in one and two dimensions are given to demonstrate the flexibility and efficiency of our moving mesh scheme in solving both smooth problems and problems with corner singularities.

The plasmonics of graphene and other two-dimensional materials has attracted enormous amounts of attention in the scientific literature over the past decade. Both the possibility of exciting plasmons in the terahertz to midinfrared regime, and the active tunability of graphene via electrical gating or chemical doping has generated a great deal of excitement among engineers seeking sensing devices which operate in this regime. Consequently, there is significant demand for robust and highly accurate computational capabilities which can incorporate such materials. Standard volumetric approaches can answer this demand, but require vast computational resources in exchange. Here we describe an algorithm which addresses this issue in two ways, first, we model the graphene layer with a surface current which is applicable to a wide class of two-dimensional materials. In addition, we reformulate the governing volumetric equations in terms of surface quantities using Dirichlet-Neumann Operators. These surface equations can be numerically simulated in an efficient, stable, and accurate fashion using a novel High-Order Perturbation of Envelopes methodology. We utilize an implementation of this algorithm to study absorbance spectra of TM polarized plane-waves scattered by a periodic grid of graphene ribbons.

A new multi-component diffuse interface model with the Peng-Robinson equation of state is developed. Initial values of mixtures are given through the NVT flash calculation. This model is physically consistent with constant diffusion parameters, which allows us to use fast solvers in the numerical simulation. In this paper, we employ the scalar auxiliary variable (SAV) approach to design numerical schemes. It reformulates the proposed model into a decoupled linear system with constant coefficients that can be solved fast by using fast Fourier transform. Energy stability is obtained in the sense that the modified discrete energy is non-increasing in time. The calculated interface tension agrees well with laboratory experimental data.

In this paper, the DBSCAN (Density-Based Spatial Clustering of Applications with Noise) method is proposed to detect particle clusters in particle-fluid systems. The particles are grouped in one cluster when they are connected by a dense environment. The parameters that define the dense environment are determined by analyzing the structure of the system, therefore, our approach needs little human intervention. The method is illustrated by identifying the clusters in two kinds of simulation trajectories of different particle-fluid systems. The robustness of cluster identification in terms of statistical properties of clusters in the steady state is demonstrated.