Preconditioned Bi-CGSTAB Method Research Articles

General framework for re-assuring numerical reliability in parallel Krylov solvers: A case of bi-conjugate gradient stabilized methods

Parallel implementations of Krylov subspace methods often help to accelerate the procedure of finding an approximate solution of a linear system. However, such parallelization coupled with asynchronous and out-of-order execution often makes more visible the non-associativity impact in floating-point operations. These problems are even amplified when communication-hiding pipelined algorithms are used to improve the parallelization of Krylov subspace methods. Introducing reproducibility in the implementations avoids these problems by getting more robust and correct solutions. This paper proposes a general framework for deriving reproducible and accurate variants of Krylov subspace methods. The proposed algorithmic strategies are reinforced by programmability suggestions to assure deterministic and accurate executions. The framework is illustrated on the preconditioned BiCGStab method and its pipelined modification, which in fact is a distinctive method from the Krylov subspace family, for the solution of non-symmetric linear systems with message-passing. Finally, we verify the numerical behavior of the two reproducible variants of BiCGStab on a set of matrices from the SuiteSparse Matrix Collection and a 3D Poisson’s equation.

Multiphysics Computing of Challenging Antenna Arrays Under a Supercomputer Framework

A parallel multiphysics simulation solver is developed to solve electromagnetic-thermal-mechanical coupling for some challenging large-scale antenna arrays. To achieve high scalability of supercomputer architectures, we reconstruct the preconditioned BiCGSTAB method and the non-overlapping domain decomposition method, so that the most resource-intensive matrix factorization steps can be performed in parallel independently within subdomains. The electromagnetic and thermal fields are solved separately, while coupled through the dissipated power and the temperature-dependent material parameters; after thermal steady state is reached, the mechanical simulation is stimulated subject to the temperature rise. The accuracy of electromagnetic-thermal coupling and thermal stress solution are first validated, and then the strong/weak parallel scalability experiments of the developed multiphysics solver are performed on supercomputer. Finally, an extremely challenging antenna array is simulated using the proposed solver, where to our best knowledge we bring the scale of multiphysics simulations excited by frequency-domain electromagnetic fields to the order of billion unknowns for the first time.

A high-throughput hybrid task and data parallel Poisson solver for large-scale simulations of incompressible turbulent flows on distributed GPUs

The solution of the pressure Poisson equation arising in the numerical solution of incompressible Navier–Stokes equations (INSE) is by far the most expensive part of the computational procedure, and often the major restricting factor for parallel implementations. Improvements in iterative linear solvers, e.g. deploying Krylov-based techniques and multigrid preconditioners, have been successfully applied for solving the INSE on CPU-based parallel computers. These numerical schemes, however, do not necessarily perform well on GPUs, mainly due to differences in the hardware architecture. Our previous work using many P100 GPUs of a flagship supercomputer showed that porting a highly optimized MPI-parallel CPU-based INSE solver to GPUs, accelerates significantly the underlying numerical algorithms, while the overall acceleration remains limited (Zolfaghari et al. [3]). The performance loss was mainly due to the Poisson solver, particularly the V-cycle geometric multigrid preconditioner. We also observed that the pure compute time for the GPU kernels remained nearly constant as grid size was increased. Motivated by these observations, we present herein an algebraically simpler, yet more advanced parallel implementation for the solution of the Poisson problem on large numbers of distributed GPUs. Data parallelism is achieved by using the classical Jacobi method with successive over-relaxation and an optimized iterative driver routine. Task parallelism is enhanced via minimizing GPU-GPU data exchanges as iterations proceed to reduce the communication overhead. The hybrid parallelism results in nearly 300 times less time-to-solution and thus computational cost (measured in node-hours) for the Poisson problem, compared to our best-case scenario CPU-based parallel implementation which uses a preconditioned BiCGstab method. The Poisson solver is then embedded in a flow solver with explicit third-order Runge-Kutta scheme for time-integration, which has been previously ported to GPUs. The flow solver is validated and computationally benchmarked for the transition and decay of the Taylor-Green Vortex at Re=1600 and the flow around a solid sphere at ReD=3700. Good strong scaling is demonstrated for both benchmarks. Further, nearly 70% lower electrical energy consumption than the CPU implementation is reported for Taylor-Green vortex case. We finally deploy the solver for DNS of systolic flow in a bileaflet mechanical heart valve, and present new insight into the complex laminar-turbulent transition process in this prosthesis.

A Parallel Preconditioned Solver for the 3D Helmholtz Equation

Discretization of the 3D Helmholtz equation leads to a linear system of a large size. The resulting linear system is difficult to be solved by a sequential method. In this paper, we propose a parallel preconditioned Bi-CGSTAB method for solving the system based on MPI(Message Passing Interface). We precondition the 3D Helmholtz equation by the complex shifted-Laplacian preconditioner, and employ the Krylov subspace method Bi-CGSTAB combined with the multi-grid to solve the resulting system. Numerical experiments are presented to illustrate the efficiency of the parallel Preconditioned solver.

A preconditioned fast finite difference scheme for space-fractional diffusion equations in convex domains

A fast finite difference method is developed for solving space-fractional diffusion equations with variable coefficient in convex domains using a volume penalization approach. The resulting coefficient matrix can be written as the discretized matrix from the extended rectangular domain plus a diagonal matrix with jumping entries due to the penalization parameter. An efficient preconditioner is constructed based on the combination of two approximate inverse circulant matrices. The preconditioned BiCGSTAB method, with the proposed preconditioner, is implemented for solving the resulting linear system. Numerical results are carried out to demonstrate the utility of the proposed algorithm.

Numerical Simulation on the Adhesive Contact Between Rough Surfaces with Bi-Conjugate Gradient Stabilized Method

A numerical model with bi-conjugate gradient stabilized method (Bi-CGSTAB) is proposed to simulate the adhesive contact between rough surfaces. This model provides a self-consistent solution of surface separation and contact pressure throughout rough surface with the adhesive interactions governed by the Lennard-Jones potential. This model uses the inexact Newton method and the preconditioned Bi-CGSTAB method for each Newton iteration. Both left and right preconditioners are employed. This model is proved to be accurate and feasible and is suitable for the moderate adhesive contact between rough surfaces.

Parallel multilevel recursive approximate inverse techniques for solving general sparse linear systems

In this article, a new parallel multilevel algebraic recursive generic approximate inverse solver (PMARGAIS) is proposed. PMARGAIS utilizes the parallel modified generic factored approximate sparse inverse (PMGenFAspI) matrix technique designed for shared memory parallel systems. PMARGAIS requires a block independent set reordering scheme, to create a hierarchy of levels. A modified block breadth first search (MBBFS) is proposed for reducing memory requirements and retaining load balancing. The SVD method is used to compute the inverse of the independent blocks that are formed from the reordering scheme, and computes accurately the Schur complement that is used as a coefficient matrix on the next level, resulting in a hybrid direct-iterative method for large linear systems. The solution of the linear system at the last level is performed with the parallel explicit preconditioned BiCGSTAB method in conjunction with the PMGenFAspI matrix. The parallelization of the proposed methods uses the vector units of modern CPUs. Implementation details are provided and numerical results are given demonstrating the applicability and effectiveness of the proposed schemes.

Efficient spectral–Galerkin methods for fractional partial differential equations with variable coefficients

Efficient spectral-Galerkin algorithms are developed to solve multi-dimensional fractional elliptic equations with variable coefficients in conserved form as well as non-conserved form. These algorithms are extensions of the spectral-Galerkin algorithms for usual elliptic PDEs developed in 24. More precisely, for separable FPDEs, we construct a direct method by using a matrix diagonalization approach, while for non-separable FPDEs, we employ a preconditioned BICGSTAB method with a suitable separable FPDE with constant-coefficients as preconditioner. The cost of these algorithms is of O ( N d + 1 ) flops where d is the space dimension. We derive rigorous weighted error estimates which provide more precise convergence rate for problems with singularities at boundaries. We also present ample numerical results to validate the algorithms and error estimates.

A multigrid-based preconditioned solver for the Helmholtz equation with a discretization by 25-point difference scheme

In this paper, a preconditioned iterative method is developed to solve the Helmholtz equation with perfectly matched layer (Helmholtz-PML equation). The complex shifted-Laplacian is generalized to precondition the Helmholtz-PML equation, which is discretized by an optimal 25-point finite difference scheme that we presented in Chen et al. (2011). A spectral analysis is given for the discrete preconditioned system from the perspective of linear fractal mapping, and Bi-CGSTAB is used to solve it. The multigrid method is employed to invert the preconditioner approximately, and a new matrix-based prolongation operator is constructed in the multigrid cycle. Numerical experiments are presented to illustrate the efficiency of the multigrid-based preconditioned Bi-CGSTAB method with the new prolongation operator. Numerical results are also given to compare the performance of the new prolongation operator with that of the prolongation operator based on the algebraic multigrid (AMG) principle.

Multiple Solution of Systems of Linear Algebraic Equations by an Iterative Method with Recomputed Preconditioners

An algorithm for multiple solution of systems of linear algebraic equations by the preconditioned BiCGStab method is proposed. When the number of iterations for solving a current system becomes too large, the precontioner is recomputed. Results of numerical experiments on computing the capacitance matrices of a microstrip line in a wide range of its sizes are presented. Bibliography: 7 titles.

A dispersion minimizing finite difference scheme and preconditioned solver for the 3D Helmholtz equation

In this paper, a new 27-point finite difference method is presented for solving the 3D Helmholtz equation with perfectly matched layer (PML), which is a second order scheme and pointwise consistent with the equation. An error analysis is made between the numerical wavenumber and the exact wavenumber, and a refined choice strategy based on minimizing the numerical dispersion is proposed for choosing weight parameters. A full-coarsening multigrid-based preconditioned Bi-CGSTAB method is developed for solving the linear system stemming from the Helmholtz equation with PML by the finite difference scheme. The shifted-Laplacian is extended to precondition the 3D Helmholtz equation, and a spectral analysis is given. The discrete preconditioned system is solved by the Bi-CGSTAB method, with a multigrid method used to invert the preconditioner approximately. Full-coarsening multigrid is employed, and a new matrix-based prolongation operator is constructed accordingly. Numerical results are presented to demonstrate the efficiency of both the new 27-point finite difference scheme with refined parameters, and the preconditioned Bi-CGSTAB method with the 3D full-coarsening multigrid.

A multigrid-based preconditioned Krylov subspace method for the Helmholtz equation with PML

In this paper, we generalize the complex shifted Laplacian preconditioner to the complex shifted Laplacian-PML preconditioner for the Helmholtz equation with perfectly matched layer (Helmholtz-PML equation). The Helmholtz-PML equation is discretized by an optimal 9-point difference scheme, and the preconditioned linear system is solved by the Krylov subspace method, especially by the biconjugate gradient stabilized method (Bi-CGSTAB). The spectral analysis of the linear system is given, and a new matrix-based interpolation operator is proposed for the multigrid method, which is used to approximately invert the preconditioner. The numerical experiments are presented to illustrate the efficiency of the preconditioned Bi-CGSTAB method with the multigrid based on the new interpolation operator, also, numerical results are given for comparing the performance of the new interpolation operator with that of classic bilinear interpolation operator and the one suggested in Erlangga et al. (2006) [10].

ICONE19-43939 Phenomenological evaluation of laser-irradiated welding processes with a combined use of higher-accuracy experiments and computational science methodologies : (5) Numerical simulation of the welding processes with a multi-dimensional multi-physics analysis code SPLICE

It is quite important that the establishing the welding and solidification processes, generation of the residual stress, in a laser welding for various reactor components. In order to standardize laser welding repair processes and controlling a residual stress which is induced by laser welding, we constructed the fully parallelized laser welding simulation code using one-fluid model (solid, liquid and gas phases are simultaneously calculated by one set of governing equations) and some advanced numerical models, e.g., VSIA (Volume and Surface Integrated Average) based CIP finite volume method for the discretization, the THINC/WLIC scheme for an accurate interface capturing and the robust and stable pressure Poisson equation solver, AMG preconditioned BiCGSTAB method. In the simulation, the base material is a pure aluminum which was included to the code as a physical parameter and we considered the surface tension force and its effect of a temperature gradient named Marangoni effect. As a result, reasonable results were obtained that is welding bead which is one of the representative behavior of a low power density laser welding and the appropriate shape of the vertical cross section of the molten pool. Therefore, the model can be applied to practical laser welding problems and contribute the standardization of a laser welding repair technology.

Preconditioned Nullspace Method for the Two-Dimensional Oseen Problem

The Oseen problem, which arises in the simulation of the time-dependent Navier-Stokes equations for incompressible fluid flow, leads to indefinite, nonsymmetric and possibly ill-conditioned linear systems of equations. This paper presents a method for obtaining a reduced linear system from the original system, which is then solved by the preconditioned BiCGStab method. The system reduction is obtained through an efficient implicit representation of a basis of discretely divergence-free functions, also known as the nullspace method. We will show a close relationship between this preconditioned nullspace method and the BFBt-preconditioned pressure Schur complement problem arising in typical block preconditioners. We will use this relationship to analyze the spectrum of the preconditioned reduced system. We will also present numerical tests for the two-dimensional Oseen problem to illustrate the performance of this method, which is similar to the performance of the BFBt-preconditioned pressure Schur complement solver. In particular, it is robust w.r.t. the strength of the viscosity for simple constant wind but not for more complicated flows, and it shows only a slight dependence on the mesh size.

Optimization of the numerical algorithms of the ADREA-I mesoscale prognostic meteorological model for real-time applications

The real-time applicability of the ADREA-I prognostic mesoscale meteorological model was enhanced by applying the preconditioned BiCGSTAB method for the numerical solution of the pressure equation in combination with increasing the magnitude of the time steps up to the values allowed by the Courant number. The ILU, MILU, ILUT and ILUM preconditioning methods with different ordering strategies were used. The implementation was developed for arbitrarily complex geometries. The application of MILU(1) preconditioning and ILUT preconditioning with red–black ordering of the unknowns (RB + ILUT) has resulted in up to six times shorter overall computational time in comparison to the previously implemented line relaxation (LR) method. The feasibility of increasing the time steps has been proved by comparing the results of the 24-h ADREA-I forecasts with the observations during a real sea breeze case in Attiki, Greece: decreasing the time steps by a factor of 10 in comparison with the values allowed by the Courant number leads to decrease of the statistical error indicators by only 1–5%.

Robust numerical valuation of European and American options under the CGMY process

We develop an implicit discretization method for pricing European and American options when the underlying asset is driven by an infinite activity Lévy process. For processes of finite variation, quadratic convergence is obtained as the mesh and timestep are refined. For infinite variation processes, better than first-order accuracy is achieved. The jump component in the neighborhood of log jump size zero is specially treated by using a Taylor expansion approximation and the drift term is dealt with using a semi-Lagrangian scheme. The resulting partial integro-differential equation is then solved using a preconditioned BiCGSTAB method coupled with a fast Fourier transform. Proofs of fully implicit timestepping stability and monotonicity are provided. The convergence properties of the BiCGSTAB scheme are discussed and compared with a fixed point iteration. Numerical tests showing the convergence and performance of this method for European and American options under processes of finite and infinite variation are presented.

Partitioning Sparse Matrices for Parallel Preconditioned Iterative Methods

This paper addresses the parallelization of the preconditioned iterative methods that use explicit preconditioners such as approximate inverses. Parallelizing a full step of these methods requires the coefficient and preconditioner matrices to be well partitioned. We first show that different methods impose different partitioning requirements for the matrices. Then we develop hypergraph models to meet those requirements. In particular, we develop models that enable us to obtain partitionings on the coefficient and preconditioner matrices simultaneously. Experiments on a set of unsymmetric sparse matrices show that the proposed models yield effective partitioning results. A parallel implementation of the right preconditioned BiCGStab method on a PC cluster verifies that the theoretical gains obtained by the models hold in practice.

Preconditioned techniques for solving large sparse linear systems arising from the discretization of the elliptic partial differential equations

In this paper, we use the BiCG, BiCGSTAB methods as preconditioned techniques. Also we compare the preconditioned Krylov subspace methods such as GMRES, GMRES(m), QMR, BiCG, CGS, BiCGSTAB for solving linear systems arising from a class of fourth-order approximations for solving the elliptic partial differential equation Au xx + Bu yy = f( x, y, u, u x , u y ), where A and B are constants. Numerical results are given to show the efficiency of the proposed preconditioned BiCGSTAB method.

Ｌａｎｄａｕ‐Ｌｉｆｓｈｉｔｚ‐Ｇｉｌｂｅｒｔ方程式のＣｒａｎｋ‐Ｎｉｃｏｌｓｏｎ法による陰解法

The efficiency of implicit methods in comparison with the classical Runge-Kutta method is investigated when the former are applied to the Landau-Lifshitz-Gilbert equation. The Crank-Nicolson and the two-stage fourth-order implicit Runge-Kutta methods were chosen and implemented by means of Newton's method. Implicit solutions using Newton's method call for the solution of large-scale simultaneous algebraic equations. To decrease the computing time, only the nearest-neighbor interaction was taken into consideration in the coefficient matrix of the simultaneous algebraic equations, which were solved by using the preconditioned Bi-CGSTAB method. Through calculations of three-dimensional magnetization reversal in a fine magnetic particle, it was confirmed that the examined methods yield a much larger maximum usable time step than the classical Runge-Kutta method. It was found that the Crank-Nicolson method takes about a quarter of the computing time required by the classical Runge-Kutta method when a dissipation constant α of 0.01 and a cell size of 8 nm are used. The maximum usable time step yielded by the two-stage fourth-order Runge-Kutta method was found to be almost four times as large as that yielded by the Crank-Nicolson method, though the former method required a longer computing time than the latter in the present calculation.

A NUMERICAL STUDY OF SOLIDIFICATION IN THE PRESENCE OF A FREE SURFACE UNDER MICROGRAVITY CONDITIONS

A numerical study of the relative importance of Marangoni effects under microgravity conditions is presented. The mathematical formulation adopted is based on the enthalpy porosity method. One of the advantages of the fixed grid method is that a unique set of equations and boundary conditions is used for the whole domain, including both solid and liquid phases. The governing equations written in a vorticity-velocity formulation are discretized using a finite volume technique on a staggered grid. A fully implicit method has been adopted for the mass and momentum equations, while the temperature field is solved separately in order to evaluate the variation in the local liquid mass fraction. The resulting algebraic system of equations is solved using a preconditioned BI-CGStab method. Numerical results modelling the free surface, including the effects on it of Marangoni convection, are presented. The influence of the presence of argon in the gap above the free surface is investigated. During the numerical simulations presented in this paper 161 2 41 and 641 2 161 uniform meshes on the whole computational domain for values of Marangoni number ( Ma ) up to 16,120 and Rayleigh number ( Ra ) of 5 have been used.