Preconditioned Conjugate Gradient Research Articles

Heterogeneous parallel computing method for 3D transient nonlinear thermomechanical problems on CPU-GPU platforms

A heterogeneous parallel computing method for 3D transient nonlinear thermomechanical problems is proposed based on CPU-GPU platforms. The method is developed based on the stable node-based smoothed finite element method (SNS-FEM) and element-by-element (EBE) strategy. The SNS-FEM method ensures the accuracy of temperature and displacement solution. To obtain higher computational efficiency, a series of parallel computing strategies are proposed in this paper. A one-to-one mapping relationship is established between smoothing domains calculation and threads in GPU. To improve the efficiency of constructing the smoothing domain, a semi-parallel construction method is developed. And a preindexing method is proposed to improve the computational efficiency of solving linear equations by the parallel preconditioned conjugate gradient (PCG) method. The developed heterogeneous parallel computing method is quite accurate and efficient in computing the model composed of tetrahedral elements, obtaining a speedup ratio of 52.31 when the model contains 67,204 degrees of freedom.

A kernel-free boundary integral method for the nonlinear Poisson-Boltzmann equation

This work proposes a boundary integral formulation based Cartesian grid method for the nonlinear Poisson-Boltzmann (PB) interface problem in biophysics. The method (a) does not have the limitation associated with the standard finite difference method for complex interfaces, (b) avoids the generation of any unstructured volume or surface grids needed by the finite element method, and (c) solves the nonlinear and variable coefficient PDE in the framework of boundary integral equations. The method solves the nonlinear PB equation with the Newton iterative method. It first reformulates the linearized PB equation in each Newton iteration as a Fredholm system of two boundary integral (BI) equations of the second kind, which is well-conditioned, and then solves the BI system with a Krylov subspace method. The evaluation of boundary and volume integrals involved in the solution of the BI system is done with a kernel-free method, which does not need to know Green's function associated with the variable coefficient PB equation. The kernel-free method evaluates volume and boundary integrals by solving equivalent simple interface problems on Cartesian grids with either a fast Fourier transform based Poisson solver or a geometric multigrid preconditioned conjugate gradient solver, either of which involves computational work essentially linearly proportional to the number of nodes on the Cartesian grid. Numerical examples in both two and three space dimensions are presented to demonstrate the efficiency and accuracy of the proposed numerical method.

Unleashing the Power of Graph Spectral Sparsification for Power Grid Analysis via Incomplete Cholesky Factorization

Graph spectral sparsification based preconditioning technique has shown promising results for power grid analysis. However, the conventional methods converge slowly for high accuracy requirement. In this work, we propose an efficient approach to address this issue. Instead of using Cholesky factorization, we employ incomplete Cholesky factorization to factorize the spectral sparsifier. We also propose a concept of graph spectral pattern, which can further reduce the preconditioned conjugate gradient (PCG) iterations using less number of nonzeros. Experiments show that under 10-6 relative tolerance, our proposed preconditioning technique achieves 1.17× speedup compared to AMGPCG in average; compared to the conventional spectral sparsification based preconditioning techniques, our proposed approach achieves up to 8.53× speedup of the factorization, 8.74× speedup of the PCG iteration, and 5.6× speedup of the total time. Moreover, the speedup of the total time continues to enlarge for higher accuracy requirement, e.g., 10-12. Last, but not the least, our method is compatible with existing graph spectral sparsification algorithms for power grid analysis.

Impact of correlated observation errors on the conditioning of variational data assimilation problems

SummaryAn important class of nonlinear weighted least‐squares problems arises from the assimilation of observations in atmospheric and ocean models. In variational data assimilation, inverse error covariance matrices define the weighting matrices of the least‐squares problem. For observation errors, a diagonal matrix (i.e., uncorrelated errors) is often assumed for simplicity even when observation errors are suspected to be correlated. While accounting for observation‐error correlations should improve the quality of the solution, it also affects the convergence rate of the minimization algorithms used to iterate to the solution. If the minimization process is stopped before reaching full convergence, which is usually the case in operational applications, the solution may be degraded even if the observation‐error correlations are correctly accounted for. In this article, we explore the influence of the observation‐error correlation matrix () on the convergence rate of a preconditioned conjugate gradient (PCG) algorithm applied to a one‐dimensional variational data assimilation (1D‐Var) problem. We design the idealized 1D‐Var system to include two key features used in more complex systems: we use the background error covariance matrix () as a preconditioner (B‐PCG); and we use a diffusion operator to model spatial correlations in and . Analytical and numerical results with the 1D‐Var system show a strong sensitivity of the convergence rate of B‐PCG to the parameters of the diffusion‐based correlation models. Depending on the parameter choices, correlated observation errors can either speed up or slow down the convergence. In practice, a compromise may be required in the parameter specifications of and between staying close to the best available estimates on the one hand and ensuring an adequate convergence rate of the minimization algorithm on the other.

Hybrid eigensolvers for nuclear configuration interaction calculations

We examine and compare several iterative methods for solving large-scale eigenvalue problems arising from nuclear structure calculations. In particular, we discuss the possibility of using block Lanczos method, a Chebyshev filtering based subspace iterations and the residual minimization method accelerated by direct inversion of iterative subspace (RMM-DIIS) and describe how these algorithms compare with the standard Lanczos algorithm and the locally optimal block preconditioned conjugate gradient (LOBPCG) algorithm. Although the RMM-DIIS method does not exhibit rapid convergence when the initial approximations to the desired eigenvectors are not sufficiently accurate, it can be effectively combined with either the block Lanczos or the LOBPCG method to yield a hybrid eigensolver that has several desirable properties. We will describe a few practical issues that need to be addressed to make the hybrid solver efficient and robust.

A robust, open-source implementation of the locally optimal block preconditioned conjugate gradient for large eigenvalue problems in quantum chemistry

We present two open-source implementations of the locally optimal block preconditioned conjugate gradient (lobpcg) algorithm to find a few eigenvalues and eigenvectors of large, possibly sparse matrices. We then test lobpcg for various quantum chemistry problems, encompassing medium to large, dense to sparse, well-behaved to ill-conditioned ones, where the standard method typically used is Davidson’s diagonalization. Numerical tests show that while Davidson’s method remains the best choice for most applications in quantum chemistry, LOBPCG represents a competitive alternative, especially when memory is an issue, and can even outperform Davidson for ill-conditioned, non-diagonally dominant problems.

Thermal Convection in a Heated-Block Duct with Periodic Boundary Conditions by Element-by-Element Treatment

The periodic nature of stream-wise flow occurs in a cooling channel so frequently due to the multiple heat sources in electronic equipment, demanding the creation of an effective technique to improve the heat-cooling convection. This work explores thermal convection enhancement in a heated-block duct for periodic boundary conditions using the element-by-element (EBE) treatment in a semi-implicit projection finite element method (FEM) through a preconditioned conjugate gradient (PCG) solver. The effects of changing the Reynolds numbers (100, 175, and 250) on rectangular cylinders installed in the channel under periodic boundary conditions were studied using time-mean Nusselt number enhancement, friction factor enhancement, and thermal performance coefficient. The results show that the rectangular cylinders installed stream-wise above an upstream block promote thermal convection in the heated-block duct due to modifying the flow of no cylinders. However, increasing the number of rectangular cylinders increases the friction factor enhancement. As a result, the case for periodic boundary conditions with a rectangular cylinder above every two blocks has the best thermal performance coefficient.

Matrix-oriented FEM formulation for reaction-diffusion PDEs on a large class of 2D domains

For the spatial discretization of elliptic and parabolic partial differential equations (PDEs), we provide a Matrix-Oriented formulation of the classical Finite Element Method, called MO-FEM, of arbitrary order k∈N. On a quite general class of 2D domains, namely separable domains, and even on special surfaces, the discrete problem is then reformulated as a multiterm Sylvester matrix equation where the additional terms account for the geometric contribution of the domain shape. By considering the IMEX Euler method for the PDE time discretization, we obtain a sequence of these equations. To solve each multiterm Sylvester equation, we apply the matrix-oriented form of the Preconditioned Conjugate Gradient (MO-PCG) method with a matrix-oriented preconditioner that captures the spectral properties of the Sylvester operator. Solving the Poisson problem and the heat equation on some separable domains by MO-FEM-PCG, we show a gain in computational time and memory occupation wrt the classical vector PCG with same preconditioning or wrt a LU based direct method. As an application, we show the advantages of the MO-FEM-PCG to approximate Turing patterns on some separable domains and cylindrical surfaces for a morphochemical reaction-diffusion PDE system for battery modelling.

Gravity interpretation using modified preconditioned conjugate gradient and reweighting of inversion methods: a case study from Iran

Gravity interpretation using modified preconditioned conjugate gradient and reweighting of inversion methods: a case study from Iran

Optimal Scheme to Achieve Energy Conservation in Induced Dipole Models.

Induced dipole models have proven to be effective tools for simulating electronic polarization effects in biochemical processes, yet their potential has been constrained by energy conservation issue, particularly when historical data is utilized for dipole prediction. This study identifies error outliers as the primary factor causing this failure of energy conservation and proposes a comprehensive scheme to overcome this limitation. Leveraging maximum relative errors as a convergence metric, our data demonstrates that energy conservation can be upheld even when using historical information for dipole predictions. Our study introduces the multi-order extrapolation method to quicken induction iteration and optimize the use of historical data, while also developing the preconditioned conjugate gradient with local iterations to refine the iteration process and effectively remove error outliers. This scheme further incorporates a "peek" step via Jacobi under-relaxation for optimal performance. Simulation evidence suggests that our proposed scheme can achieve energy convergence akin to that of point-charge models within a limited number of iterations, thus promising significant improvements in efficiency and accuracy.

Low-Rank Matrix Sensing-Based Channel Estimation for mmWave and THz Hybrid MIMO Systems

This paper studies the channel estimation for wideband multiple-input multiple-output (MIMO) systems equipped with hybrid analog/digital transceivers operating in the millimeter-wave (mmWave) or terahertz (THz) bands. By exploiting the low-rank property of the concatenated channel matrix of the delay taps, we formulate the channel estimation problem as a low-rank matrix sensing (LRMS) problem and solve it using a low-complexity generalized conditional gradient-alternating minimization (GCG-ALTMIN) algorithm. This LRMS-based solution can accommodate different precoder/combiner and training structures. In addition, it does not require knowledge about the array responses at the transceivers, in contrast to most existing solutions allowing low training overhead. Furthermore, a preconditioned conjugate gradient (PCG) algorithm-based implementation and a low-rank matrix completion (LRMC) formulation are proposed to further reduce the computational complexity. In order to enhance the channel estimation performance for fat and tall channel matrices, we introduce a matrix reshaping approach that can preserve the channel rank by exploiting the shift-invariance property of uniform arrays. We also introduce a spectrum denoising (SD) approach for further improving the performance when the array responses are known and the number of paths is small. These approaches can effectively enhance the performance at a given training overhead. Simulation results suggest that the proposed solutions can achieve higher channel estimation accuracy and reduce the computational complexity as compared to several representative channel estimation schemes.

Analytical and Numerical Results for the Diffusion-Reaction Equation When the Reaction Coefficient Depends on Simultaneously the Space and Time Coordinates

We utilize the travelling-wave Ansatz to obtain novel analytical solutions to the linear diffusion–reaction equation. The reaction term is a function of time and space simultaneously, firstly in a Lorentzian form and secondly in a cosine travelling-wave form. The new solutions contain the Heun functions in the first case and the Mathieu functions for the second case, and therefore are highly nontrivial. We use these solutions to test some non-conventional explicit and stable numerical methods against the standard explicit and implicit methods, where in the latter case the algebraic equation system is solved by the preconditioned conjugate gradient and the GMRES solvers. After this verification, we also calculate the transient temperature of a 2D surface subjected to the cooling effect of the wind, which is a function of space and time again. We obtain that the explicit stable methods can reach the accuracy of the implicit solvers in orders of magnitude shorter time.

Iterative solution of spatial network models by subspace decomposition

We present and analyze a preconditioned conjugate gradient method (PCG) for solving spatial network problems. Primarily, we consider diffusion and structural mechanics simulations for fiber based materials, but the methodology can be applied to a wide range of models, fulfilling a set of abstract assumptions. The proposed method builds on a classical subspace decomposition into a coarse subspace, realized as the restriction of a finite element space to the nodes of the spatial network, and localized subspaces with support on mesh stars. The main contribution of this work is the convergence analysis of the proposed method. The analysis translates results from finite element theory, including interpolation bounds, to the spatial network setting. A convergence rate of the PCG algorithm, only depending on global bounds of the operator and homogeneity, connectivity and locality constants of the network, is established. The theoretical results are confirmed by several numerical experiments.

GPU-based mesh reduction strategy utilizing active nodes for structural topology optimization

This paper presents a novel mesh reduction strategy on graphics processing unit (GPU) that aims to reduce the computational complexity of conventional structural topology optimization. In the proposed strategy, the effective number of design variables is reduced by using the concept of active nodes and active elements in the finite element mesh. A novel mesh numbering scheme is also introduced to facilitate parallel identification of active nodes using a proposed GPU-based algorithm. The preconditioned conjugate gradient (PCG) solver is further developed using the proposed strategy and the numbering scheme. The proposed strategy is tested on three structural topology optimization problems using solid isotropic material with penalization method. Results of the proposed strategy demonstrate up to 8× speedup over the standard GPU implementation considering the entire mesh.

Adaptive FEM with quasi-optimal overall cost for nonsymmetric linear elliptic PDEs

Abstract We consider a general nonsymmetric second-order linear elliptic partial differential equation in the framework of the Lax–Milgram lemma. We formulate and analyze an adaptive finite element algorithm with arbitrary polynomial degree that steers the adaptive meshrefinement and the inexact iterative solution of the arising linear systems. More precisely, the iterative solver employs, as an outer loop, the so-called Zarantonello iteration to symmetrize the system and, as an inner loop, a uniformly contractive algebraic solver, for example, an optimally preconditioned conjugate gradient method or an optimal geometric multigrid algorithm. We prove that the proposed inexact adaptive iteratively symmetrized finite element method leads to full linear convergence and, for sufficiently small adaptivity parameters, to optimal convergence rates with respect to the overall computational cost, i.e., the total computational time. Numerical experiments underline the theory.

Convergence acceleration of preconditioned conjugate gradient solver based on error vector sampling for a sequence of linear systems

AbstractIn this article, we focus on solving a sequence of linear systems that have identical (or similar) coefficient matrices. For this type of problem, we investigate subspace correction (SC) and deflation methods, which use an auxiliary matrix (subspace) to accelerate the convergence of the iterative method. In practical simulations, these acceleration methods typically work well when the range of the auxiliary matrix contains eigenspaces corresponding to small eigenvalues of the coefficient matrix. We develop a new algebraic auxiliary matrix construction method based on error vector sampling in which eigenvectors with small eigenvalues are efficiently identified in the solution process. We use the generated auxiliary matrix for convergence acceleration in the following solution step. Numerical tests confirm that both SC and deflation methods with the auxiliary matrix can accelerate the solution process of the iterative solver. Furthermore, we examine the applicability of our technique to the estimation of the condition number of the coefficient matrix. We also present the algorithm of the preconditioned conjugate gradient method with condition number estimation.

A mixed precision LOBPCG algorithm

The locally optimal block preconditioned conjugate gradient (LOBPCG) algorithm is a popular approach for computing a few smallest eigenvalues and the corresponding eigenvectors of a large Hermitian positive definite matrix A. In this work, we propose a mixed precision variant of LOBPCG that uses a (sparse) Cholesky factorization of A computed in reduced precision as the preconditioner. To further enhance performance, a mixed precision orthogonalization strategy is proposed. To analyze the impact of reducing precision in the preconditioner on performance, we carry out a rounding error and convergence analysis of PINVIT, a simplified variant of LOBPCG. Our theoretical results predict and our numerical experiments confirm that the impact on convergence remains marginal. In practice, our mixed precision LOBPCG algorithm typically reduces the computation time by a factor of 1.4--2.0 on both CPUs and GPUs.

A modified and efficient phase field model for the biological transport network

This paper aims to establish the biological transport network based on the phase field model. In order to ensure that the topological shape is formed under the guidance of electrical conductivity, we generate the biological networks with sufficient information based on the distributions of the venation of the leaf represented by the reaction-diffusion model. We modify the original energy of the network generating model by considering the auxin gradient property. By applying the gradient flow method to minimize the modified energy, we derive the Poisson type equation for pressure, the reaction-diffusion type equation for the network conductance, and the Allen-Cahn type equation for the phase field. The proposed model is significant on the investigation of phase transitions by considering the gradient properties on the boundaries. We have innovatively added conductivity and phase-field coupling terms that inhibit perpendicular transport of nutrients, making it easy to generate thin branches from the trunk through this model. In order to obtain the second-order temporal accuracy, we take the Crank–Nicolson method for the governing system. To obtain the second-order spatial accuracy, we discretize the coupling system with the central finite difference method and linearize the nonlinear terms semi-explicitly to form a linear system at each time step. The discrete energy dissipation is provably preserved and we can use a larger time step. We apply the preconditioned conjugate gradient method with the multigrid method as a preconditioner to implement a practical algorithm with only linear algebraic complexity. The proposed algorithm is easy to implement and achieves a fast convergence. Various numerical tests are demonstrated to verify the efficiency, stability, and robustness of the proposed method.

Mixed finite elements based on superconvergent patch recovery for strain gradient theory

Strain gradient theory is one of the most powerful candidates to simulate micro/nano-sized structures considering size effect. However, too many degrees of freedom (DOFs) have been placed to satisfy the C1-continuity condition for the development of finite elements based on the strain gradient theory, greatly reducing the practicality of these elements. In this work, four-node quadrilateral and eight-node hexahedral elements are developed to simulate size effect based on strain gradient theory. Since only displacement DOFs are placed on each node without using displacement gradient DOFs, there is dramatic reduction in computational cost. The displacements are approximated with shape functions of standard finite elements. The displacement gradient field is approximated with linear polynomials and the coefficients are determined using the superconvergent patch recovery. It was found that the superconvergent properties of Barlow points are still valid even though the gradient recovery is performed during the calculation of the element stiffness matrix. The preconditioned conjugate gradient method is utilized to solve the system equations. Numerical examples demonstrate the accuracy and efficiency of the proposed elements. Size effect observed in experiments is captured by using the developed finite elements.

On solving large-scale multistage stochastic optimization problems with a new specialized interior-point approach

A novel approach based on a specialized interior-point method (IPM) is presented for solving large-scale stochastic multistage continuous optimization problems, which represent the uncertainty in strategic multistage and operational two-stage scenario trees. This new solution approach considers a split-variable formulation of the strategic and operational structures. The specialized IPM solves the normal equations by combining Cholesky factorizations with preconditioned conjugate gradients, doing so for, respectively, the constraints of the stochastic formulation and those that equate the split-variables. We show that, for multistage stochastic problems, the preconditioner (i) is a block-diagonal matrix composed of as many shifted tridiagonal matrices as the number of nested strategic-operational two-stage trees, thus allowing the efficient solution of systems of equations; (ii) its complexity in a multistage stochastic problem is equivalent to that of a very large-scale two-stage problem. A broad computational experience is reported for large multistage stochastic supply network design (SND) and revenue management (RM) problems. Some of the most difficult instances of SND had 5 stages, 839 million linear variables, 13 million quadratic variables, 21 million constraints, and 3750 scenario tree nodes; while those of RM had 8 stages, 278 million linear variables, 100 million constraints, and 100,000 scenario tree nodes. For those problems, the proposed approach obtained the solution in 1.1 days using 174 gigabytes of memory for SND, and in 1.7 days using 83 gigabytes for RM; while CPLEX v20.1 required more than 53 days and 531 gigabytes for SND, and more than 19 days and 410 gigabytes for RM.