In this paper, a discrete multiresolution analysis of the flow field with moving interface is built up based on the interpolating errors of the sign function of level set function on a nested grid structure. By using the multiresolution coefficients to establish the local grid size and computational scheme, we constructe a class of multiscale level set methods. For flow region near the moving interface where the magnitudes of the multiresolution coefficients are large, the high-order WENO scheme is used for time evolution. While in the rest computational region, solutions are obtained directly by polynomial interpolation. Compared with the single-scale level set method, this method can capture more sophisticated and sharper motion interface with less CPU expenses.

From the viewpoint of distributed memory and parallel computing, the overlapping additive methods, such as AS (classical additive Shcwarz), RAS (restricted additive Schwarz) and PAS (parallel additive Schwarz) were inspected again. And a new kind of additive method, called "Distributed Additive Schwarz Method" (DAS), was proposed and analyzed. Numerical experiments showed that DAS was more efficient than other additive methods in terms of iteration number and iteration time.

H-tensors have wide applications in the scientific computation and the applications in engineering, but it is not easy to determine whether a given tensor is an H-tensor or not in practice. In this paper, we give some practical criteria for H-tensors by constructing different positive diagonal matrices and applying some techniques of inequalities. As an application, some sufficient conditions of the positive definiteness for an even-order real symmetric tensor are given. Advantages of results obtained are illustrated by numerical examples.

Sparse matrix algorithms are one of the most challenging topics in the supercomputing area. This paper from three aspects, i.e., high scalability, high performance and high practicability, provides an overview of the most important research work published in the last 30 years. Then a large amount of experimental results from evaluating over ten sparse BLAS algorithms on three GPUs are analyzed, and the difficulties of getting the three achievements are discussed. Finally, a serial of challenges of the next generation sparse matrix computations are proposed.

Algebraic multigrid(AMG) is one of the most efficient methods for solving large-scale linear systems arising from the discretized PDEs, it is widely used in numerical simulations for scientific and engineering applications. As the supercomputers become more and more powerful, both application features and machine architectures become more and more complicated, which lead to great challenges in scalability and robustness for AMG. This paper reviews the state-of-the art and challenges of AMG algorithms, including the parallel scalability, algorithmic scalability and float-point performance optimization for extreme-scale applications, especially for the upcoming exascale era.

This paper gives a method for removing periodical constraint conditions under the condition of maintaining features of the system of finite element algebraic equations. The feature of this method is that the number of variables and order of one-dimensioal arrangement remains unchanged, symmetry and positive definite properties are maintained and the cholesky method can still be applied. It is proved that this method is extremely convenient and effect in practice.

In this paper, a class of compact difference schemes are constructed to solve the nonlinear delay hyperbolic partial differential equations. The unique solvability, convergence and unconditional stability of the scheme are obtained. The convergence order is O(τ²+h⁴). Furthermore, the Richardson extrapolation is applied to improve the temporal accuracy of the scheme, and a solution of order four in both temporal and spatial dimensions is obtained. Numerical example shows the accuracy and efficiency of the algorithms.

As the nonlocal model can describe the singularity and discontinuity mechanism of some important physical phenomena, it has been widely used in many fields in recent years and has played a strong role in promoting the development of related disciplines. As an important nonlocal model, anomalous diffusion model is used to describe anomalous diffusion phenomena. The nonlocality and multi-scale characteristics of nonlocal models not only promote the discovery of new mathematical theories, but also provide a new perspective for the existing discrete and local continuous models and their connections. Although there have been many achievements, there is a broad space for the development of multi-scale nonlocal and anomalous diffusion models, whether from the perspective of mathematical methods, basic theories or numerical methods. It is a research focus to further develop and improve the basic mathematical theories and methods, and develop new efficient numerical schemes under the condition of real solution regularity, especially the numerical schemes with stability, convergence and asymptotic compatibility. In the past few years, the author of this paper has been committed to the research of mathematical theory and numerical methods of nonlocal models, and has made some meaningful research results in the design of artificial boundary conditions, nonlocal maximum principle and asymptotically compatible numerical schemes. In the numerical analysis of anomalous diffusion equations, the fast algorithm of Caputo derivative and the discrete fractional Grönwall inequality are developed, and the idea of error convolution structure is proposed to represent the local consistent error, which provides some basic analysis frameworks for the optimal error estimation of a class of widely-used variable-step-size numerical schemes. There is still a long way to completely solve various mathematical problems in nonlocal and anomalous diffusion models, which needs further study. It is hoped that this paper can play a role in promoting the in-depth development of basic theories and algorithms of multi-scale nonlocal and anomalous diffusion models.

The split Perfectly Matched Layer (PML) Absorbing Boundary Condition (ABC) accomplished the absorption of nonphysical reflections, preliminarily realized the numerical modeling of open infinite space to the finite area, however, on the edge of the ABC area, the field component in calculation need to be split, which will increasing the number of independent equations and the calculation capacity in Maxwell's equations. The Uniaxial-anisotropic Perfectly Matched Layer (UPML) boundary conditions do not need to split the field component, and the iterative formula is easy convenient to programming, except to the decay parameters in the conventional UPML, another two real variable factors are introduced, which is exclusively for the reflections of low frequency. The TMZ wave equations of improved UPML is deduced, the dielectric parameter distribution mode is given, and the program implementation steps of the algorithm was introduced in detail,then the merits were verified by a numerical modeling example of electromagnetic waves, with the split PML, the conventional UPML and the improved UPML boundary conditions, respectively, by comparing the absorption effect from the aspects of the snapshots, and the time domain reflection errors and the frequency domain reflection errors, the results show that the improved UPML boundary condition possess an apparent advantage over absorbing the low frequency reflections in the latter period of the wave propagation, which is more realistically simulates the infinite and open space with little truncation.

We designed encoding and decoding algorithms for high dimensional Hilbert order. Hilbert order has good locality, and it has wide applications in various fields in computer science, such as memory management, database, and dynamic load balancing. We analyzed existing algorithms for computing 2D and 3D Hilbert order, and designed improved algorithms for computing Hilbert order in arbitrary space dimensions. We also proposed an alternate form of Hilbert space filling curve which has the advantage of preserving the ordering between different levels.

Abstract It is shown that the simple iterative methods are available and efficient for solving implicit symplectic schemes.The numerical test supports this conclusioa.

The computation of the Moore-Penrose inverse and Drazin inverse of real matrix can be transformed into solving the problem of linear matrix equations. Then the modified conjugate gradient method can be used to get the general solution of linear matrix equations. Finally, the Moore-Penrose inverse and Drazin inverse of real matrix can be obtained directly or through matrix multiplication. The modified conjugate gradient method is different from the usual conjugate gradient method. It which does not require the positive definite, reversible or full column rank of the coefficient matrix of the involved linear algebraic equations. Thus this method is always feasible. The numerical experiments show that the algorithm is effective.

In this paper, we present a preconditioning algorithm for sparse, symmetric, diagonally-dominant(SDD) linear systems by using combinatorial preconditioning techniques. Comparing to incomplete LU factorization preconditioners and AMG preconditioners, combinatorial preconditioners have good parallel scalability. Moreover they meet the flux conservation and quivalent resistance principle. As a SDD matrix can be converted to a Laplacian matrix which is isomorphic to an undireceted graph. Based on this fact, we first construct a low stretch spanning tree of the graph corresponding to the SDD coefficient matrix by using Ofer's et al algorithm. Then we partition the tree into subtrees based on a tree-decomposition algorithm and add appropriate edges to get an stretch optimized subgraph. Finally, convert the subgraph to a SDD matrix and take it as a preconditioner. We show experimentally that these combinatorial preconditioners are more efficient than incomplete Cholesky factorization preconditioners with modification and without modification. Moreover, these combinatorial precontioners are insensitive to the matrix ordering and problem boundary.

As one of the most important data assimilation methods in numerical weather predication, 3DVAR (three-dimensional variational data assimilation) can improve the quality of data predication significantly. With gradual development of scientific research and the improvement of detecting instruments and computer technology, the traditional sequential 3DVAR system can no longer meet the demand for the high resolution and high accuracy numerical forecasting due to computational capacity and memory limit. So the parallel design and implementation of 3DVAR is very essential. In the paper, a mixed domain decomposition parallel method and message communication library are applied in the 3DVAR system of the state meteorological agency. Numerical experiments show that the 3DVAR system could run with high speedup ratio and parallel efficiency.

In this paper, we study the parallel algorithm based on CUDA and MPI for the Fast Fourier Transform on the hexagon (FFTH) and its implementation. Firstly, we design a CUDA FFTH algorithm by utilizing the hierachical parallelization mechanism and the build-in CUFFT library for classic rectangular FFTs. With respect to the serial cpu program, our CUDA program achieves 12x speedup for 3×2048² double-precision complex-to-complex FFTH. If we ignore the PCI between main memory and GPU device memory, around 30x-40x speedup can be even achieved. Although the non-tensorial FFTH is much more complicated than the rectangular FFT, our CUDA FFTH program gains the same efficiency as the rectangular CUFFT. Next, efforts are mainly contributed to optimization techniques for parallel array transposition and data sorting, which significantly improve the efficiency of the CUDA-MPI FFTH algorithm. On a 10-node cluster with 60 GPUs, our CUDA-MPI program achieves about 55x speedup with respect to the the serial cpu program for 3×8192² complex-to-complex double-precision FFTH, and it is more efficient than the MPI parallel FFTW. Our research on the CUDA-MPI algorithm for FFTH is beneficial to the exploration and development of new parallel algorithms on large-scale CPU-GPU heterogeneous computer systems.

Some bounds for the Perron root ρ of nonnegative matrices are established. Let A be any nonnegative matrix, f(A) a polynomial of A satisfied f(A)≥0 and all the row sums of f(A) be nonzero, then min(r_i(A·f(A))/r_i(f(A)))≤ρ≤r_i(A·f(A))/r_i(f(A)).
This result is a generalized form of bounds in paper[4-7] and can improve the estimation of ρ by choosing appropriate polynomials.

The oil- and gas-bearing formations are widely distributed in China, but the geological structure is complex and the natural energy is usually not sufficient. Numerical simulation methods and software are important tools for reservoir engineers to analyze and manage reservoirs and are among the main tools for determining the remaining oil distribution and enhancing oil/gas recovery ratios in the late stage of reservoir development. The fine geological models increase the spatial resolution, leading to a huge number of grid cells and costly simulation runs, which present many new challenges to numerical algorithms. In this paper, we briefly introduce a simplified composition model, its mathematical model, discretization methods, solution methods, parallel computation, and implementation, with emphasis on several preconditioning and decoupling methods now adopted by commercial and industrial software.

This paper proposes an algorithm for the generation and meshing of multi-scale composites 2D 3D unit cell model. The generation of unit cell model uses compact and reject algorithm, we use ellipse ellipsoid to represent its internal reinforcing particles. This algorithm could generate a large number of ellipses and high volume ratio of the unit cell with good randomness. When meshing, first we use adaptive line approximation algorithm that uses polygon to approximate the ellipse and polyhedron to approximate the ellipsoid. The Numerical experiments show that this algorithm generates high quality meshes in a very fast speed and the validity of the method.

Based on the method of the modified conjugate gradient to the linear matrix equation over constrained matrices, and by modifying the construction of some matrices, an iterative algorithm is presented to find the generalized reflexive solution of the matrix equation which is a special type with several matrix variables. The convergence of the iterative algorithm is proved. And the problem of the optimal approximation to the given matrix is solved in the generalized reflexive solution set of this matrix equation. When this matrix equation is consistent, its generalized reflexive solution can be obtained within finite iterative steps. And its least-norm generalized reflexive solution can be got by choosing the special initial matrices. The numerical example shows that the iterative algorithm is quite efficient.

In this paper, an e ective scheme is given, which can solve eigenvalue for large sparse matrixes based on spectrum division method. In this scheme, we divide the spectrum into several intervals and extract eigenvalues from each subinterval independently. To estimate the eigenvalue distribution, the Gerschgorin Circle Theorem is used. In addition, an inter- polation method, aimed to find the subintervals, is applied to reduce calculation amount. In this paper, a multistage parallel algorithm, combined with spectral projection method and approximate spectral projection method based on contour integral, is given. The validity, balance, and e ciency of this algorithm are veried in supercomputer "shenteng7000" and "yuan". Working on 1024 processors, the performance of this parallel algorithm improves five times compared with the general algorithm.

A new kind of transform, block-pulse function transform, and some fie1ds of its applicationare introduced briefly. Two numerical examples are given.

We investigate the time fractional anomalous diffusion equation on a bounded domain. We propose an efficient method for its numerical solution. This method is based on a finite difference in time and spectral method in space. The numerical examples show the convergence rate is O(Δt^3-α+N^-m), where α,Δt,N and m are respectively the order of time fractional derivatives, time step size, the polynomial degree and the regularity of the exact solution.

We studied existing hp adaptive strategies in the literature, and proposed a new strategy for adaptive finite element method. Numerical experiments show that our new strategy achieves exponential convergence, and is superior, in both precision of the solutions and computation time, to the strategy compared. This part of the work also serves to validate the hp adaptivity module of PHG.

Delaunay triangulation has been widely used in many fields such as compu- tational fluid dynamics, statistics, meteorology solid state physics, computational geometry and so on. Bowyer-Watson algorithm is a very popular one for generating Delaunay triangulation. In generating the Delaunay triangulation of a preassigned set of n points, the complexity of Bowyer-Watson algorithm can at most be reduced to O(n log n) for the simple reason that the complexity of its tree search process is O(nlog n). In this paper we suggest a tree search technique whose complexity is O(n). Noting that the order of point insertion can affect the efficiency of Bowyer- Watson algorithm, we propose a technique to optimize the point insertion process. Based on these two techniques, we obtain a fast algorithm for generating Delaunay triangulation.

In this paper, we consider a fast and high-precision algorithm of the basic elementary functions. Firstly, we discuss the power series expansion of functions which are related to the Bernoulli number B_2n or Euler number E_2n, such as tan x, sec x, tanh x and so on, and study the corresponding fast computation. For the basic elementary functions, hyperbolic functions and inverse hyperbolic functions, we derive a fast and arbitrarily accurate algorithm based on the power series expansion in the complex domain. The algorithm proposed in this paper is suitable for all elementary functions due to that the exponential and logarithm functions can be expressed by the power series. The feature of this algorithm is that it can be easily coded and it is self-contained for the computation of the elementary functions.

A new three-stage third-order solution (NTSTO) to symplectical schemes is obtained based on partitioned Runge-Kutta form, several numerical results show that, the scheme is excellent in suppressing residual increase. Combining finite difference(FD) in spatial discretization with the new symplectical scheme in temporal discretization, the further numerical experiments are provided and the results also show that the method is effective. The form is completely consistent with the requirements of the symplectical scheme in time domain when the elastic wave equations are discretized using spectral element methods(SEM) in space domain. So it is natural to solve the elastic wave equations using the symplectical scheme combined with the spectral element methods (NTSTO-SEM). Finally, the algorithm is employed to simulate wave propagation in transversely isotropic media, the results show that the performance is good and superior to classical algorithms, such as Newmark method combined with SEM (Newmark-SEM) and Runge-Kutta method of third-order combined with SEM (RK3-SEM).

A new algorithm named MIC-GABP algorithm is proposed for solving large scale symmetric sparse linear system. This algorithm is based on GaBP(Gauss Belief Propagation) and MIC(Many Integerated Core). We take several large scale sparse matrices from The University of Florida Sparse Matrix Collection as examples to observe the performance of our algorithm. The experimental result shows that MIC-GaBP algorithm has a higher efficiency than traditional GaBP algorithm under the same accuracy.

Sparse matrix-vector multiplication (SpMV) is an important computational kernel in scientific and engineering applications. The performance of SpMV by using traditional CSR format is often far below 10% of the peak performance on modern processors with memory hierarchy. When using the CSR format for SpMV, it is often hard to directly take advantage of the SIMD acceleration technology on mordern processors, due to irregular memory access pattern. In order to use the SIMD technology, a new storage format for sparse matrices, CSRL (Compressed Sparse Row with Local information), is proposed.The CSRL format has locality characteristic, and is SIMD-friendly. The new format reduces the number of memory access and improves the SpMV performance. Experiments show that, compared with the implementation in Intel MKL library (version 10.3), the SpMV based on the CSRL format gains an average of 29.5% and maximum of 89%performance improvement.

By using the Legendre polynomials approximation theory and Gauss-Lobatto quadrature formula, a four-stage fourth-order implicit Runge-Kutta method is presented. It is showed that the new algorithm has good stability properties in theoretical analysis, A(α)-stable and α is close to ninety degrees, and stiff stable and D is close to zero. It is almost A-stable and almost L-stable. The new method can solve stiff ordinary differential equations effectively. The numerical examples illustrate its effectiveness.

Abstract In this paper, a necessary and sufficient condition for constructing a shape preserving C1 cubic interpolation spline function is obtained, and two new simple algorithms. one for directly constructing shape preserving C1 cubic interpolation spline functions, and the other for constructing shape preserving C1 piecewise cubic interpolation functions by insetting some new nodes, are given.

The numerical method of pricing up-and-out call Parision Option based on the Black-Scholes model is focused in this article. The two-point compact scheme with second-order accuracy is used. A technique to remove the singularity of the pay-off function is used to make the result more accurate,more effective and more stable. The influence of the delaying time and the barrier on the option price is discussed.

A common task in geometric modeling is to reconstruct a data set, sampled from a surface, with a parametric polynomial or spline surface. Bézier surfaces, B-spline surfaces, or Radial basis functions are often used to deal with this problem. When the parametric domain of the data set is a convex polygon, it is necessary to subdivide the parametric domain, then consider the continuity of the adjacent patches. In this paper, reconstructing the surface by toric Bézier patch is a entire scheme without subdividing the parametric domain. Some example are given to demonstrate that the method is available and efficient.

This paper presents a parallel conjugate gradient method to solve the large-scale sparse matrix eigenvalue problem. Design a matrix's row-partition method with load balance,implement a matrix rearranging method for the overlap of computation and communication, decrease the traffic between processors through preprocess. In addition, we present a MPI and OpenMP parallel methods. We take some numerical experiments under ShenTeng7800 shows that keep constant communication time under different processors. Shows that algorithm has high parallel efficiency under multi-core.

In the paper, we give in a structural sense a necessary and sufficient condition of thecomplete controllability on biological linear compartmental model. The condition is also a criticalmethod and is easily realized on computer. The critical stepes for using the method is given,which are considered to be implemented easily on computr.

Chemotaxis is the central mechanism for bacteria to find food and flee from poisons. In multicellular organisms, chemotaxis is critical for the development and health. Along with the understanding of more biological details of chemotaxis, the mathematical models become more complete and complex. We review three models for E.coli chemotaxis at different scales: the Keller-Segel model proposed in the 70’s last century, the velocity jump model given in the 90s and the pathway-based kinetic model proposed at the beginning of this century. We review the biological background of the proposed models, the modeling ideas, some main analytical results about the features of their solutions, the numerical methods and discuss about the connections and differences between models at different scales.

The Toeplitz matrix is invertible if two sets systems of Toeplitz equations is solvable. Also, the inversion of a Toeplitz matrix can be denoted as a sum of products of circulant matrices and lower triangular Toeplitz matrices. Fast algorithm for the inversion of a Toeplitz matrix withing O(n2) operations (rather than O(n3), as required by standard matrix inversion methods) are derived.

A second order wave equation of KdV type, a important nonlinear wave equation, has broad application prospect. The equation was studied based on the multi-symplectic theory in Hamilton space. The multi-symplectic Fourier pseudospectral method is reviewed, and a semi-implicit scheme with certain discrete conservation laws is constructed to solve the second order wave equation of KdV type. The multi-symplectic scheme of the second order wave equation of KdV type with several conservation laws are presented. The numerical results for the second order wave equation of KdV type are reported, showing that the multi-symplectic Fourier pseudospectral method is an efficient algorithm with excellent long-time numerical behaviors.

Gradient Projection method is a classic algorithm for solving constrained optimization. It takes the advantage of its low per-iteration cost. However, its efficiency can be affected by the choice of step size. In this paper, we propose a new Gradient Projection method with self-adaptive step size rule. On the hand, it does not need the information of objective; on the other hand, it can accept longer step size than that based on Armijo rule to accelerate the convergence. Our preliminary experimental results show that it is an efficient algorithm.

Since nonsingular H-matrices play important roles in the theory and application of Numerical Linear Algebra, it is necessary to know whether a given matrix is a nonsingular H-matrix. In this paper, a group of some sufficient conditions for nonsingular H-matrices is given according to the properties of generalized strictly α-chain diagonally dominant matrices, generalized strictly α- diagonally dominant matrices, and the method of applying the iteration factors. In addition, it improves some related results and enriches the theory of nonsingular H-matrices. Finally, the effectiveness of these sufficient conditions is illustrated with numerical examples.

In this paper a parallel finite element solver for biomolecular simulations is introduced. This solver is based on the three dimensional parallel finite element toolbox PHG and able to simulate the diffusion process influenced by electrostatic field in biological solution. The solver is developed from our previous work with new algorithms added, now offering two time-dependent algorithms and four steady-state algorithms for the Poisson-Nernst-Planck equations, while one steady-state algorithm is provided for Smoluchowski-Poisson-Boltzmann equations. The solver is able to simulate the biomolecular models, including ion channel and nanopore, and solve electrostatic field and concentration distributions for all species. Furthermore, current and reaction rate is calculated to study the functions of ion channel and enzyme.