In the paper, we develop and analyze a new mass-preserving splitting domain decomposition method over multiple sub-domains for solving parabolic equations. The domain is divided into non-overlapping multi-bock sub-domains. On the interfaces of sub-domains, the interface fluxes are computed by the semi-implicit (explicit) flux scheme. The solutions and fluxes in the interiors of sub-domains are computed by the splitting one-dimensional implicit solution-flux coupled scheme. The important feature is that the proposed scheme is mass conservative over multiple non-overlapping sub-domains. Analyzing the mass-preserving S-DDM scheme is difficult over non-overlapping multi-block sub-domains due to the combination of the splitting technique and the domain decomposition at each time step. We prove theoretically that our scheme satisfies conservation of mass over multi-block non-overlapping sub-domains and it is unconditionally stable. We further prove the convergence and obtain the error estimate in L2-norm. Numerical experiments confirm theoretical results.

In this paper, a lattice Boltzmann BGK (LBGK) model is proposed for simulating incompressible axisymmetric flows. Unlike other existing axisymmetric lattice Boltzmann models, the present LBGK model can eliminate the compressible effects only with the small Mach number limit. Furthermore the source terms of the model are simple and contain no velocity gradients. Through the Chapman-Enskog expansion, the macroscopic equations for incompressible axisymmetric flows can be exactly recovered from the present LBGK model. Numerical simulations of the Hagen-Poiseuille flow, the pulsatile Womersley flow, the flow over a sphere, and the swirling flow in a closed cylindrical cavity are performed. The results agree well with the analytic solutions and the existing numerical or experimental data reported in some previous studies.

We study the construction of symplectic Runge-Kutta methods for stochastic Hamiltonian systems (SHS). Three types of systems, SHS with multiplicative noise, special separable Hamiltonians and multiple additive noise, respectively, are considered in this paper. Stochastic Runge-Kutta (SRK) methods for these systems are investigated, and the corresponding conditions for SRK methods to preserve the symplectic property are given. Based on the weak/strong order and symplectic conditions, some effective schemes are derived. In particular, using the algebraic computation, we obtained two classes of high weak order symplectic Runge-Kutta methods for SHS with a single multiplicative noise, and two classes of high strong order symplectic Runge-Kutta methods for SHS with multiple multiplicative and additive noise, respectively. The numerical case studies confirm that the symplectic methods are efficient computational tools for long-term simulations.

In this paper, we develop, analyze and test local discontinuous Galerkin (LDG) methods for solving the Degasperis-Procesi equation which contains nonlinear high order derivatives, and possibly discontinuous or sharp transition solutions. The LDG method has the flexibility for arbitrary h and p adaptivity. We prove the L2 stability for general solutions. The proof of the total variation stability of the schemes for the piecewise constant P0 case is also given. The numerical simulation results for different types of solutions of the nonlinear Degasperis-Procesi equation are provided to illustrate the accuracy and capability of the LDG method.

We develop the dimension-reduced hyperbolic moment method for the Boltzmann equation, to improve solution efficiency using a numerical regularized moment method for problems with low-dimensional macroscopic variables and high-dimensional microscopic variables. In the present work, we deduce the globally hyperbolic moment equations for the dimension-reduced Boltzmann equation based on the Hermite expansion and a globally hyperbolic regularization. The numbers of Maxwell boundary condition required for well-posedness are studied. The numerical scheme is then developed and an improved projection algorithm between two different Hermite expansion spaces is developed. By solving several benchmark problems, we validate the method developed and demonstrate the significant efficiency improvement by dimension-reduction.

We couple a node-centered finite volume method to a high order finite difference method to simulate dynamic earthquake ruptures along nonplanar faults in two dimensions. The finite volume method is implemented on an unstructured mesh, providing the ability to handle complex geometries. The geometric complexities are limited to a small portion of the overall domain and elsewhere the high order finite difference method is used, enhancing efficiency. Both the finite volume and finite difference methods are in summation-by-parts form. Interface conditions coupling the numerical solution across physical interfaces like faults, and computational ones between structured and unstructured meshes, are enforced weakly using the simultaneous-approximation-term technique. The fault interface condition, or friction law, provides a nonlinear relation between fields on the two sides of the fault, and allows for the particle velocity field to be discontinuous across it. Stability is proved by deriving energy estimates; stability, accuracy, and efficiency of the hybrid method are confirmed with several computational experiments. The capabilities of the method are demonstrated by simulating an earthquake rupture propagating along the margins of a volcanic plug.

This paper presents new finite difference schemes for solving the Helmholtz equation in the polar and spherical coordinates. The most important result presented in this study is that the developed difference schemes are pollution free, and their convergence orders are independent of the wave number k. Let h denote the step size, it will be demonstrated that when solving the Helmholtz equation at large wave numbers and considering kh is fixed, the errors of the proposed new schemes decrease as h decreases even when k is increasing and kh > 1.

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.

It is well known that conventional edge elements in solving vector Maxwell's eigenvalue equations by the finite element method will lead to the presence of spurious zero eigenvalues. This problem has been addressed for the first order edge element by Kikuchi by the mixed element method. Inspired by this approach, this paper describes a higher order mixed spectral element method (mixed SEM) for the computation of two-dimensional vector eigenvalue problem of Maxwell's equations. It utilizes Gauss-Lobatto-Legendre (GLL) polynomials as the basis functions in the finite-element framework with a weak divergence condition. It is shown that this method can suppress all spurious zero and nonzero modes and has spectral accuracy. A rigorous analysis of the convergence of the mixed SEM is presented, based on the higher order edge element interpolation error estimates, which fully confirms the robustness of our method. Numerical results are given for homogeneous, inhomogeneous, L-shape, coaxial and dual-inner-conductor cavities to verify the merits of the proposed method.

Detecting corrosion by electrical field can be modeled by a Cauchy problem of Laplace equation in annulus domain under the assumption that the thickness of the pipe is relatively small compared with the radius of the pipe. The interior surface of the pipe is inaccessible and the nondestructive detection is solely based on measurements from the outer layer. The Cauchy problem for an elliptic equation is a typical ill-posed problem whose solution does not depend continuously on the boundary data. In this work, we assume that the measurements are available on the whole outer boundary on an annulus domain. By imposing reasonable assumptions, the theoretical goal here is to derive the stabilities of the Cauchy solutions and an energy regularization method. Relationship between the proposed energy regularization method and the Tikhonov regularization with Morozov principle is also given. A novel numerical algorithm is proposed and numerical examples are given.

The paper is concerned with the unconditional stability and error estimates of fully discrete Galerkin-Galerkin FEMs for the equations of incompressible miscible flows in porous media. We prove that the optimal L2 error estimates hold without any time-step (convergence) conditions, while all previous works require certain time-step restrictions. Theoretical analysis is based on a splitting of the error into two parts: the error from the time discretization of the PDEs and the error from the finite element discretization of the corresponding time-discrete PDEs, which was proposed in our previous work [26, 27]. Numerical results for both two and three-dimensional flow models are presented to confirm our theoretical analysis.

The Jiles-Atherton (J-A) model is a commonly used physics-based model in describing the hysteresis characteristics of ferromagnetic materials. However, citations and interpretation of this model in literature have been non-uniform. Solution methods for solving numerically this model has not been studied adequately. In this paper, through analyzing the mathematical properties of equations and the physical mechanism of energy conservation, we point out some unreasonable descriptions of this model and develop a relatively more accurate, modified J-A model together with its numerical solution method. Our method employs a fixed point method to compute anhysteretic magnetization. We obtain the susceptibility value of the anhysteretic magnetization analytically and apply the 4th order Runge-Kutta method to the solution of total magnetization. Computational errors are estimated and then precisions of the solving method in describing various materials are verified. At last, through analyzing the effects of the accelerating method, iterative error and step size on the computational errors, we optimize the numerical method to achieve the effects of high precision and short computing time. From analysis, we determine the range of best values of some key parameters for fast and accurate computation.

Combining asymptotic analysis and weighted L2 stability estimates, the convergence of lattice Boltzmann methods for the approximation of 1D convection-diffusion-reaction equations is proved. Unlike previous approaches, the proof does not require transformations to equivalent macroscopic equations.

In this paper we present recent developments concerning a Cell-Centered Arbitrary Lagrangian Eulerian (CCALE) strategy using the Moment Of Fluid (MOF) interface reconstruction for the numerical simulation of multi-material compressible fluid flows on unstructured grids in cylindrical geometries. Especially, our attention is focused here on the following points. First, we propose a new formulation of the scheme used during the Lagrangian phase in the particular case of axisymmetric geometries. Then, the MOF method is considered for multi-interface reconstruction in cylindrical geometry. Subsequently, a method devoted to the rezoning of polar meshes is detailed. Finally, a generalization of the hybrid remapping to cylindrical geometries is presented. These explorations are validated by mean of several test cases using unstructured grid that clearly illustrate the robustness and accuracy of the new method.

We analyze in this paper the pressure splitting scheme of a partitioned semi-implicit coupling algorithm for fluid-structure interaction (FSI) simulation. The semi-implicit coupling algorithm is developed on the ground of the artificial compressibility characteristic-based split (AC-CBS) scheme that serves not only for the fluid subsystem but also for the global FSI system. As the dual-time stepping procedure recommended for quasi-incompressible flows is incorporated into the implicit coupling stage, the fluctuating pressure may be unusually susceptible to the AC coefficient. Moreover, it is not trivial to devise an optimal AC formulation for pressure estimation. Instead, we consider a stabilized second-order pressure splitting scheme in the AC-CBS-based partitioned semi-implicit coupling algorithm. Computer simulation of a benchmark FSI experiment demonstrates that good agreement is exposed between the available and present data.

In this article the instabilities appearing in a liquid layer are studied numerically by means of the linear stability method. The fluid is confined in an annular pool and is heated from below with a linear decreasing temperature profile from the inner to the outer wall. The top surface is open to the atmosphere and both lateral walls are adiabatic. Using the Rayleigh number as the only control parameter, many kind of bifurcations appear at moderately low Prandtl numbers and depending on the Biot number. Several regions on the Prandtl-Biot plane are identified, their boundaries being formed from competing solutions at codimension-two bifurcation points.

A numerical method based on a homogeneous single-phase flow model is presented to simulate the interaction between pressure wave and flow cavitation. To account for compressibility effects of liquid water, cavitating flow is assumed to be compressible and governed by time-dependent Euler equations with proper equation of state (EOS). The isentropic one-fluid formulation is employed to model the cavitation inception and evolution, while pure liquid phase is modeled by Tait equation of state. Because of large stiffness of Tait EOS and great variation of sound speed in flow field, some of conventional compressible gasdynamics solvers are unstable and even not applicable when extended to calculation of flow cavitation. To overcome the difficulties, a Godunov-type, cell-centered finite volume method is generalized to numerically integrate the governing equations on triangular mesh. The boundary is treated specially to ensure stability of the approach. The method proves to be stable, robust, accurate, time-efficient and oscillation-free.

Novel numerical experiments are designed to investigate unsteady dynamics of the cavitating flow impacted by pressure wave, which is of great interest in engineering applications but has not been studied systematically so far. Numerical simulation indicates that cavity over cylinder can be induced to collapse if the object is accelerated suddenly and extremely high pressure pulse results almost instantaneously. This, however, may be avoided by changing the traveling speed smoothly. The accompanying huge pressure increase may damage underwater devices. However, cavity formed at relatively high upstream speed may be less distorted or affected by shock wave and can recover fully from the initial deformation. It is observed that the cavitating flow starting from a higher freestream velocity is more stable and more resilient with respect to perturbation than the flow with lower background speed. These findings may shed some light on how to control cavitation development to avoid possible damage to operating devices.

In this paper, we extend the high order finite-difference method with subcell resolution (SR) in [34] for two-species stiff one-reaction models to multispecies and multireaction inviscid chemical reactive flows, which are significantly more difficult because of the multiple scales generated by different reactions. For reaction problems, when the reaction time scale is very small, the reaction zone scale is also small and the governing equations become very stiff. Wrong propagation speed of discontinuity may occur due to the underresolved numerical solution in both space and time. The present SR method for reactive Euler system is a fractional step method. In the convection step, any high order shock-capturing method can be used. In the reaction step, an ODE solver is applied but with certain computed flow variables in the shock region modified by the Harten subcell resolution idea. Several numerical examples of multispecies and multireaction reactive flows are performed in both one and two dimensions. Studies demonstrate that the SR method can capture the correct propagation speed of discontinuities in very coarse meshes.

We study the error analysis of the weak Galerkin finite element method in [24, 38] (WG-FEM) for the Helmholtz problem with large wave number in two and three dimensions. Using a modified duality argument proposed by Zhu and Wu, we obtain the pre-asymptotic error estimates of the WG-FEM. In particular, the error estimates with explicit dependence on the wave number k are derived. This shows that the pollution error in the broken H1-norm is bounded by under mesh condition k7/2h2≤C0 or (kh)2+k(kh)p+1≤C0, which coincides with the phase error of the finite element method obtained by existent dispersion analyses. Here h is the mesh size, p is the order of the approximation space and C0 is a constant independent of k and h. Furthermore, numerical tests are provided to verify the theoretical findings and to illustrate the great capability of the WG-FEM in reducing the pollution effect.

We propose a hierarchy of novel absorbing boundary conditions for the one-dimensional stationary Schrödinger equation with general (linear and nonlinear) potential. The accuracy of the new absorbing boundary conditions is investigated numerically for the computation of energies and ground-states for linear and nonlinear Schrödinger equations. It turns out that these absorbing boundary conditions and their variants lead to a higher accuracy than the usual Dirichlet boundary condition. Finally, we give the extension of these ABCs to N-dimensional stationary Schrödinger equations.