Biconjugate Gradient Research Articles

Accelerating Contaminant Transport Simulation in MT3DMS Using JASMIN-Based Parallel Computing

To overcome the large time and memory consumption problems in large-scale high-resolution contaminant transport simulations, an efficient approach was presented to parallelize the modular three-dimensional transport model for multi-species (MT3DMS) (University of Alabama, Tuscaloosa, AL, USA) program on J adaptive structured meshes applications infrastructures (JASMIN). In this approach, a domain decomposition method and a stencil-based method were used to accomplish parallel implementation, while a ghost cell strategy was used for communication. The MODFLOW-MT3DMS coupling mode was optimized to achieve the parallel coupling of flow and contaminant transport. Five types of models were used to verify the correctness and test the parallel performance of the method. The developed parallel program JMT3D (China University of Geosciences (Beijing), Beijing, China) can increase the speed by up to 31.7 times, save memory consumption by 96% with 46 processors, and ensure that the solution accuracy and convergence do not decrease as the number of domains increases. The BiCGSTAB (Bi-conjugate gradient variant algorithm) method required the least amount of time and achieved high speedup in most cases. Coupling the flow and contaminant transport further improved the efficiency of the simulations, with a 33.45 times higher speedup achieved on 46 processors. The AMG (algebraic multigrid) method achieved a good scalability, with an efficiency above 100% on hundreds of processors for the simulation of tens of millions of cells.

Retrieval of Composite Model Parameters for 3-D Microwave Imaging of Biaxial Objects by BCGS-FFT and PSO

This article presents the 3-D quantitative microwave imaging of homogeneous biaxial anisotropic objects with their optical axes being in arbitrary directions. Two sets of model parameters, the anisotropic dielectric parameters and the angles of optical axes with respect to the user-defined coordinate system, are retrieved iteratively. In the forward model, the total fields are solved from the discretized combined field volume integral equations (CFVIEs) by the stabilized biconjugate gradient fast Fourier transform (BCGS-FFT) method. In the inversion model, the optimized solutions for two sets of model parameters are obtained by the particle swarm optimization (PSO) method. If the geometric models of the homogeneous anisotropic objects are unknown, a full-wave inversion method, the variational Born iterative method (VBIM), is combined with the multiparametric structural consistency constraint (SCC) to find the spatial boundaries of the homogeneous objects in advance. Two numerical examples are used to test the proposed method, and the results show that BCGS-FFT-PSO not only has strong adaptability for optical axes angles varying in wide ranges but also can retrieve composite model parameters of multiple homogeneous biaxial objects simultaneously even when the measured data are contaminated by noise.

A meshfree collocation method based on moving Taylor polynomial approximation for high order partial differential equations

This paper presents a meshfree collocation method for solving high order partial differential equations (PDEs). The leading numerical difficulty is the approximation of high order derivatives. To make the approximation simple and efficient, a moving Taylor polynomial (MTP) approximation is presented by using movable expansion point for each sub-domain. Derivatives can be derived straightforward from the corresponding Taylor coefficients, which are determined by solving a weighted least squares problem. A distinct feature of the method is its ability to give the derivatives along with the shape function itself without further cost. To ensure the accuracy of high order approximation, stability of the weighted least squares problems for determining the Taylor coefficients is another issue should be addressed. For this purpose, the basis functions are rescaled by the size of window functions, and QR decomposition is adopted to solve the weighted least squares problems. The collocation method based on this MTP approximation does not require any grid or background cell, so it is a truly meshfree method. When solving the linear algebraic system generated by the MTP collocation method, a preconditioned sparse biconjugate gradients stabilized (BICGSTAB) solver is used to accelerate the computation speed. Numerical tests show that the proposed method is much accurate and efficient for high order PDEs.

On Efficient Numerical Solution of Linear Algebraic Systems Arising in Goal-Oriented Error Estimates

We deal with the numerical solution of linear partial differential equations (PDEs) with focus on the goal-oriented error estimates including algebraic errors arising by an inaccurate solution of the corresponding algebraic systems. The goal-oriented error estimates require the solution of the primal as well as dual algebraic systems. We solve both systems simultaneously using the bi-conjugate gradient method which allows to control the algebraic errors of both systems. We develop a stopping criterion which is cheap to evaluate and guarantees that the estimation of the algebraic error is smaller than the estimation of the discretization error. Using this criterion and an adaptive mesh refinement technique, we obtain an efficient and robust method for the numerical solution of PDEs, which is demonstrated by several numerical experiments.

An efficient quadratic finite volume method for variable coefficient Riesz space‐fractional diffusion equations

A quadratic finite volume (FV) method for steady‐state Riesz space‐fractional diffusion equations (sFDEs) with variable diffusivity coefficient is developed using piecewise quadratic basis functions, and a resulting linear algebra system with two‐by‐two block‐type Toeplitz‐like coefficient matrix is formulated. It is proved that the method requires a minimum memory of order , where is the number of spatial partitions. Moreover, as two of the produced Toeplitz‐like submatrices are not square, a new fast nonsquare Toeplitz‐like matrix‐vector product is specially designed, which requires an almost linear computational complexity of order . Then, a fast version of biconjugate gradient stabilized (BiCGSTAB) solution algorithm, named FBiCGSTAB, is proposed for the FV scheme. The FV method combined with Crank‐Nicolson (CN) time discretization is applied to solve time‐dependent sFDEs, and an efficient CN‐FV scheme is developed and analyzed. Finally, numerical results are presented to show the utility of the fast FV and fast CN‐FV methods.

A matrix-free Chimera approach based on Dirichlet–Dirichlet coupling for domain composition purposes

A novel approach using an algebraic-based Chimera type domain composition method in the context of the finite element method for non-matching overlapping unstructured grids is proposed in this work. The scheme is based on both the transfer of information across each grid interface via Dirichlet boundary conditions and a high-order interpolation algorithm to obtain one global solution of the system. The solution can be obtained iteratively, with a convergence rate that is similar to that obtained with an analogous conformal mesh, and the matrix–vector operator can be computed with completely decoupled operations on both meshes. Furthermore, the scheme can be set as a linear operator that can be fed to a matrix-free efficient iterative solver, such as the BiConjugate Gradient Stabilized method. Several numerical examples using non-matching unstructured grids that are partially and completely overlapped with different element sizes are presented, assessing the precision and convergence rate of the method.

A meshless technique based on generalized moving least squares combined with the second-order semi-implicit backward differential formula for numerically solving time-dependent phase field models on the spheres

In the current research paper, the generalized moving least squares technique is considered to approximate the spatial variables of two time-dependent phase field partial differential equations on the spheres in Cartesian coordinate. This is known as a direct approximation (it is the standard technique for generalized finite difference scheme [69,77]), and it can be applied for scattered points on each local sub-domain. The main advantage of this approach is to approximate the Laplace-Beltrami operator on the spheres using different types of distribution points simply, in which the studied mathematical models are involved. In fact, this scheme permits us to solve a given partial differential equation on the sphere directly without changing the original problem to a problem on a narrow band domain with pseudo–Neumann boundary conditions. A second-order semi-implicit backward differential formula (by adding a stabilized term to the chemical potential that is the second-order Douglas-Dupont-type regularization) is applied to approximate the temporal variable. We show that the time discretization considered here guarantees the mass conservation and energy stability. Besides, the convergence analysis of the proposed time discretization is given. The resulting fully discrete scheme of each partial differential equation is a linear system of algebraic equations per time step that is solved via an iterative method, namely biconjugate gradient stabilized algorithm. Some numerical experiments are presented to simulate the phase field Cahn-Hilliard, nonlocal Cahn-Hilliard (for diblock copolymers as microphase separation patterns) and crystal equations on the two-dimensional spheres.

REEF3D::FNPF—A Flexible Fully Nonlinear Potential Flow Solver

Abstract In situations where the calculation of ocean wave propagation and impact on structures are required, fast numerical solvers are desired in order to find relevant wave events. Computational fluid dynamics (CFD)-based numerical wave tanks (NWTs) emphasize on the hydrodynamic details such as fluid–structure interaction, which make them less ideal for the event identification due to the large computational resources involved. Therefore, a computationally efficient numerical wave model is needed to identify the events both for offshore deep-water wave fields and coastal wave fields where the bathymetry and coastline variations have strong impact on wave propagation. In the current paper, a new numerical wave model is represented that solves the Laplace equation for the flow potential and the nonlinear kinematic and dynamics free surface boundary conditions. This approach requires reduced computational resources compared to CFD-based NWTs. The resulting fully nonlinear potential flow solver REEF3D::FNPF uses a σ-coordinate grid for the computations. This allows the grid to follow the irregular bottom variation with great flexibility. The free surface boundary conditions are discretized using fifth-order weighted essentially non-oscillatory (WENO) finite difference methods and the third-order total variation diminishing (TVD) Runge–Kutta scheme for time stepping. The Laplace equation for the potential is solved with Hypre’s stabilized bi-conjugated gradient solver preconditioned with geometric multi-grid. REEF3D::FNPF is fully parallelized following the domain decomposition strategy and the message passing interface (MPI) communication protocol. The numerical results agree well with the experimental measurements in all tested cases and the model proves to be efficient and accurate for both offshore and coastal conditions.

Hybrid multi-projection method using sparse approximate inverses on GPU clusters

The state-of-the-art supercomputing infrastructures are equipped with accelerators, such as graphics processing units (GPUs), that operate as coprocessors for each workstation of the distributed memory system. The multi-projection type methods are a class of algebraic domain decomposition methods based on semi-aggregation techniques. The multi-projection type methods have improved convergence behavior, as the number of subdomains increases, due to the corresponding augmentation of the semi-aggregated local linear systems with more coarse components, while the number of fine components is reduced. Moreover, limited amount of communications among the workstations is required by the proposed method. The utilization of the available GPUs allows an increase in the number of subdomains along with finer-grained parallelism, leading to improved performance. A load-balancing algorithm that ensures the concurrency of the computations on multicore processors and GPUs is proposed. Flexible parallel preconditioned Krylov subspace iterative methods enhanced with multi-projection type methods have been designed appropriately in order to have improved performance, compared to CPU-only or GPU-only executions, by exploiting the available CPUs and GPUs of the distributed memory system concurrently. The unsymmetric local linear systems are solved by the preconditioned Bi-Conjugate Gradient STABilized (BiCGSTAB) method enhanced with the modified generic factored approximate sparse inverse preconditioner, whereas the preconditioned conjugate gradient (CG) method along with the symmetric factored approximate sparse inverse preconditioner is used for the symmetric positive definite local coefficient matrices. Numerical results regarding the convergence behavior, the performance, and the scalability of the proposed method for several problems are given.

A Symmetric Successive Overrelaxation (SSOR) based Gauss-Seidel Massive MIMO Detection Algorithm

With an increasing number of antennas and users, the complexity increases dramatically, in order to simplify the process. In this paper, We propose existed algorithms discussing their principles and their advantages as well as the comparison between each method and then we provide the possible improvement in theory changing steps in iteration process to avoid the matrix inversion, by introducing a pre-conditioned gradient (PCG) method and further improved Jacobi method and using Neumann-series terms to estimate an approximate matrix inversion and the Gauss-Seidel (GS) method with a diagonal-approximate initial solution to the method and biconjugate gradient stabilized (BICGSTAB) method and above all of the methods are based on the minimum mean square error (MMSE) detector including the application in multi-user MIMO system. Another approach is from a signal detector based on symmetric successive over relaxation (SSOR) method without matrix inversion and Neumann series based on the zero forcing signal detector. The third algorithm uses a modified version of the conjugate gradient least square (CGLS). These methods are mainly replacing the inversion process into iteration process to reduce the complexity of the algorithm not only replacing the stages but optimizing the speed from setting the initial solution as well as the improvement using the Neumann-series to replace the iteration part. Finally, we proposed a symmetric successive overrelaxation (SSOR) based Gauss-Seidel for massive MIMO detection.

On regularized Hermitian splitting iteration methods for solving discretized almost‐isotropic spatial fractional diffusion equations

SummaryThe shifted finite‐difference discretization of the one‐dimensional almost‐isotropic spatial fractional diffusion equation results in a discrete linear system whose coefficient matrix is a sum of two diagonal‐times‐Toeplitz matrices. For this kind of linear systems, we propose a class of regularized Hermitian splitting iteration methods and prove its asymptotic convergence under mild conditions. For appropriate circulant‐based approximation to the corresponding regularized Hermitian splitting preconditioner, we demonstrate that the induced fast regularized Hermitian splitting preconditioner possesses a favorable preconditioning property. Numerical results show that, when used to precondition Krylov subspace iteration methods such as generalized minimal residual and biconjugate gradient stabilized methods, the fast preconditioner significantly outperforms several existing ones.

The Scattering of Electromagnetic Fields from Anisotropic Objects Embedded in Anisotropic Multilayers

This article presents an efficient method for calculating the scattered electromagnetic (EM) fields from anisotropic objects embedded in general anisotropic multilayers, both with arbitrary anisotropy. The electric field integral equation (EFIE) is derived based on the electric dyadic Green’s functions for general anisotropic multilayers. The volume integral equation (VIE) involves a convolution form between the dyadic Green’s function and equivalent sources. To accelerate the integral computation, the fast Fourier transform technique (FFT) is adopted which reduces the computational cost from $O(N^{2})$ to $O (N \log N)$ . The biconjugate gradient stabilized method (BiCGSTAB) is then applied to solve the EFIE efficiently to calculate the total electric fields in the computational domain. In order to validate the algorithm, several examples are modeled for anisotropic objects embedded in both isotropic and anisotropic multilayers. The modeled scattered fields are compared with the finite-element method (FEM) results. The good agreement between the two results shows the algorithm works properly.

Numerical investigation on the transport equation in spherical coordinates via generalized moving least squares and moving kriging least squares approximations

The main aim of this paper is to present new and simple numerical methods for solving the time-dependent transport equation on the sphere in spherical coordinates. We use two techniques, namely generalized moving least squares and moving kriging least squares to find the new formulations for approximating the advection operator in spherical coordinates such that any singularities have been omitted. Another feature of the methods considered here is that they do not depend on the background mesh or triangulation for approximation, which they could be applied on the transport equation in spherical coordinates easily with different distribution points. Furthermore, due to the eigenvalue stability of the dicretized advection operator via two proposed approximations, an implicit-explicit linear multistep method has been applied to discretize the time variable. The fully discrete scheme obtained here yields a linear system of algebraic equations at each time step, which is solved via the biconjugate gradient stabilized algorithm with zero-fill incomplete lower-upper (ILU) preconditioner. Three well-known test problems, namely “solid body rotation”, “vortex roll-up” and “deformational flow” are solved to demonstrate our developments. Also for the first test problem, we apply a simple positivity-preserving filter at the end of each time step, which keeps the transported variable positive.

Development of hybrid low salinity water polymer flooding numerical reservoir simulator and smart proxy model for chemical enhanced oil recovery (CEOR)

A MATLAB based finite volume CEOR simulator was developed for modeling the recovery of oil by low salinity water polymer flooding (LSWPF) in a proposed geostatistical method of generating porous media. Firstly, this paper coupled ion exchange, wettability and fluid flow transport models to the proposed heterogeneous porous media to investigate LSWPF oil recovery performance. Secondly, for the first time, a new hybrid smart proxy model of particle swarm optimization artificial neural network (PSO-ANN) recovery factor (RF) correlation was derived for LSWPF oil recovery performance using the extracted weights and biases. The proposed porous media as well as stochastic multi-fractal random field and circulant embedding method (CEM) permeability field realizations did not exhibit artificial banding unlike turning bands method (TBM) realizations when compared. A fast convergence pre-conditioner solution technique of algebraic multigrid (AMG) linear solver combined with an acceleration method of biconjugate gradients stabilized method (BI-CGSTAB) was used to solve the fully implicit system of coupled linearized equations. The simulator was validated with PHREEQC geochemical package and a secondary LSWPF coreflood test with outstanding history match. Numerical 2D examples study in this paper revealed the synergetic effects of LSWPF and high salinity water polymer flooding (HSWPF) with higher oil RFs of 39% and 31%, respectively than the standalone LSWF of 23% and normal waterflood of 11% for 50 days of production. The optimal architecture of the PSO-ANN RF model is 7 inputs, 3 hidden neurons and 1 output. Also, the PSO-ANN RF model predictions were compared to the LSWPF numerical simulation results with the former remarkably having reduced computational time in seconds with training and blind testing errors <1% average absolute percentage error (AAPE) for various permeability field realizations investigated. Finally, the PSO-ANN RF model will serve as a fast-alternative reservoir management and planning tool to supplement the high-fidelity LSWPF CEOR simulator developed to minimize computational overheads, risks and uncertainties in oil recovery performance modeling.

Resolved particle simulations using the Physalis method on many GPUs

We present a distributed memory many-GPU implementation of the Physalis method for resolving spherical particles in disperse multiphase flow simulations. The current work extends a previous single-GPU computational procedure by implementing a distributed memory Poisson solver and distributed finite-size particle methods using MPI. We document the changes required to move to a distributed memory model for both the fluid and solid phases. We benchmark the code with up to one million resolved particles in a domain size of 19203 on 216 GPUs at the Maryland Advanced Research Computing Center and present strong and weak scaling results. Finally, by taking advantage of the realization that the solution procedure for the pressure Poisson equation can be implemented using a symmetric matrix, we are able to replace the biconjugate gradient stabilized algorithm used in previous work with the conjugate gradient algorithm.

Cost‐efficient numerical algorithm for solving the linear inverse problem of finding a variable magnetization

The paper is devoted to developing an original cost‐efficient algorithm for solving the inverse problem of finding a variable magnetization in a rectangular parallelepiped. The problem is ill‐posed and is described by the integral Fredholm equation.It is shown that after discretization of the area and approximation of the integral operator, this problem is reduced to solving a system of linear algebraic equations with the Toeplitz‐block‐Toeplitz matrix. We have constructed the memory efficient variant of the stabilized biconjugate gradient method BiCGSTABmem. This optimized algorithm exploits the special structure of the matrix to reduce the memory requirements and computing time.The efficient implementation is developed for multicore CPU and GPU. A series of the model problems with synthetic and real magnetic data are solved. Investigation of efficiency and speedup of parallel algorithm is performed.

Modified Conjugate Gradient Algorithms for Gram Matrix Inversion of Massive MIMO Downlink Linear Precoding

This paper deals with various low complexity algorithms for higher order matrix inversion involved in massive MIMO system precoder design. The performance of massive MIMO systems is optimized by the process of precoding which is divided into linear and nonlinear. Nonlinear precoding techniques are most complex precoding techniques irrespective of its performance. Hence, linear precoding is generally preferred in which the complexity is mainly contributed by matrix inversion algorithm. To solve this issue, Krylov subspace algorithm such as Conjugate Gradient (CG) was considered to be the best choice of replacement for exact matrix inversions. But CG enforces a condition that the matrix needs to be Symmetric Positive Definite (SPD). If the matrix to be inverted is asymmetric then CG fails to converge. Hence in this paper, a novel approach for the low complexity inversion of asymmetric matrices is proposed by applying two different versions of CG algorithms- Conjugate Gradient Squared (CGS) and Bi-conjugate Gradient (Bi-CG). The convergence behavior and BER performance of these two algorithms are compared with the existing CG algorithm. The results show that these two algorithms outperform CG in terms of convergence speed and relative residue.

On Convergence Speed of Parallel Variants of GPBiCG Method for Solving Linear Equations

The hybrid Bi-Conjugate Gradient (Bi-CG) methods such as Bi-CG stabilized (Bi-CGSTAB), Generalized Product-type based Bi-CG (GPBiCG), and BiCGstab(ℓ) are well-known for efficiently solving linear equations. GPBiCG and BiCGstab(ℓ) are more effective and robust than Bi-CGSTAB on problems with strongly nonsymmetric matrices. On present petascale high-performance computing hardware, the scalability of Krylov subspace methods has recently become increasingly prominent. The main bottleneck for efficient parallelization is the inner products which require a global reduction. The parallel variants of Bi-CGSTAB reducing the number of global communication phases and hiding the communication latency have been proposed. However, it has been reported that the convergence of the parallel variants of Bi-CGSTAB is affected by rounding errors than that of the standard Bi-CGSTAB, and is not as robust as the standard. In this paper, therefore, following [1], we design parallel variants of GPBiCG, which converges faster and is more robust than Bi-CGSTAB. Then we compare the convergence speed between the standard GPBiCG and the parallel variants by numerical experiments.

AFODSS: Heart Rate Estimation Method From Photoplethysmographic Signals With Motion Artifacts Using Fourier-Sparse Dual Optimization

Precise heart rate (HR) estimation of individuals during physical exercise from PPG sensors using direct methods often produces erroneous output due to the occurrence of exceedingly high motion artifacts (MAs). The effect of MAs can be neutralized using non-direct algorithms such as optimization techniques using some a priori about MA. In this paper, we propose a dual-optimization technique using adaptive Fourier decomposition in conjunction with sparse spectral reconstruction (AFODSS) using a multiple measurement vectors (MMVs) model and a regularized bi-conjugate gradient focal underdetermined system solver (BCG-FOCUSS) algorithm. Using AFODSS, the average absolute error (AAE) of HR obtained on the dataset recorded from 23 subjects during their intense exercise is 1.35 beats/min. The BCG-FOCUSS algorithm used in this paper helps to speed up the convergence rate and reduce time complexity when compared with other FOCUSS algorithms.

3-D Electromagnetic Scattering and Inverse Scattering by Magnetodielectric Objects With Arbitrary Anisotropy in Layered Uniaxial Media

3-D electromagnetic (EM) scattering and inverse scattering of magnetodielectric objects with arbitrary anisotropy embedded in layered uniaxial media are studied. The stabilized biconjugate gradient fast Fourier transform (BCGS-FFT) method is employed to solve the forward scattering problem formulated by the combined field volume integral equation (CFVIE). In the inversion, we use the variational Born iterative method (VBIM) enhanced by the structural consistency constraint (SCC) to reconstruct eighteen unknown dielectric parameters of the scatterers simultaneously. In order to further improve the accuracy of the reconstructed permittivity, permeability, and conductivity of the scatterers, we propose the structural continuity scanning (SCS) technique. Comparisons between the BCGS-FFT solution and simulated results from finite element method (FEM) are performed to verify the reliability and accuracy of the forward solver. Meanwhile, we reconstruct the multiple anisotropic parameters of several typical structures embedded in layered uniaxial media to show the feasibility and anti-noise ability of the inversion algorithm.