Preconditioned Conjugate Gradient Algorithm Research Articles

Overlapping domain decomposition methods for finite volume discretizations

Two-level additive overlapping domain decomposition methods are applied to solve the linear system arising from the cell-centered finite volume discretization methods (FVMs) for the elliptic problems. The conjugate gradient (CG) methods are used to accelerate the convergence. To analyze the preconditioned CG algorithm, a discrete L2 norm, an H1 norm, and an H1 semi-norm are introduced to connect the matrices resulting from the FVMs and related bilinear forms. It has been proved that, with a small overlap, the condition number of the preconditioned systems does not depend on the number of the subdomains. The result is similar to that for the conforming finite element. Numerical experiments confirm the theory.

Direct Unconstrained Optimization of Molecular Orbital Coefficients in Density Functional Theory.

One-electron orbitals in Kohn-Sham density functional theory (DFT) are typically constrained to be orthogonal during their variational optimization, leading to elaborate parameterization of the orbitals and complicated optimization algorithms. This work shows that orbital optimization can be performed with nonorthogonal orbitals if the DFT energy functional is augmented with a term that penalizes linearly dependent states. This approach, called variable-metric self-consistent field (VM SCF) optimization, allows us to use molecular orbital coefficients, natural descriptors of one-electron orbitals, as independent variables in a direct, unconstrained minimization, leading to very simple closed-form expressions for the electronic gradient and Hessian. It is demonstrated that efficient convergence of the VM SCF procedure can be achieved with a basic preconditioned conjugate gradient algorithm for a variety of systems, including challenging narrow-gap systems and spin-pure two-determinant states of singlet diradicals. This simple reformulation of the variational procedure can be readily extended to electron correlation methods with multiconfiguration states and to the optimization of excited-state orbitals.

Lightning-fast Method of Fundamental Solutions

The method of fundamental solutions (MFS) and its associated boundary element method (BEM) have gained popularity in computer graphics due to the reduced dimensionality they offer: for three-dimensional linear problems, they only require variables on the domain boundary to solve and evaluate the solution throughout space, making them a valuable tool in a wide variety of applications. However, MFS and BEM have poor computational scalability and huge memory requirements for large-scale problems, limiting their applicability and efficiency in practice. By leveraging connections with Gaussian Processes and exploiting the sparse structure of the inverses of boundary integral matrices, we introduce a variational preconditioner that can be computed via a sparse inverse-Cholesky factorization in a massively parallel manner. We show that applying our preconditioner to the Preconditioned Conjugate Gradient algorithm greatly improves the efficiency of MFS or BEM solves, up to four orders of magnitude in our series of tests.

A small-scale active anechoic chamber

This paper presents a solution method to suppress reflections from the walls in an anechoic room using active control, in which the control coefficients are obtained with an efficient algorithm. One of the main challenges in suppressing the reflections, is the estimation of the reflected sound field, which is a non-measurable quantity. The reflected sound field is computed using the Kirchhoff-Helmholtz integral with a set of microphones on a circle. The effectiveness of the method is shown on an experimental setup. The experimental setup is a small-scale version of an acoustic anechoic chamber. An anechoic chamber, with larger dimensions and higher frequencies, will result in a system with a relatively large number of sources and sensors. It is computationally expensive to find a set of control filters for larger-scale systems. The conjugate gradient algorithm using block-circulant preconditioning is an efficient method to solve for the optimal set of control filters. The preconditioned conjugate gradient algorithm is used together with the Kirchhoff-Helmholtz integral. Real-time measurements show that the reflections are effectively suppressed, expressed in terms of the reverberation time. This results in a small-scale active anechoic chamber.

High-efficiency computation for electromagnetic forming process: An explicit-implicit GPU approach

Low computational efficiency is one of the crucial challenges in electromagnetic forming process simulation. However, the huge gap in switching from the electromagnetic field to mechanics hinders the development of high-efficiency parallel coupling. A graphics processing unit (GPU)-based explicit-implicit coupling method is proposed to address the high-computational cost in simulations of the electromagnetic forming process, in which the large deformation (explicit) and electromagnetic field (implicit) are intimately coupled. The electromagnetic-mechanical coupling is achieved by constructing five parallel algorithms, which include the parallel operations of vector-vector, matrix-matrix, matrix-vector, vector assembly, and matrix-free assembly parallel preconditioned conjugate gradient algorithm. The high-efficiency GPU parallel computing is realized through the special data storage structure and parallel strategies. The explicit-implicit GPU method is validated against the experimental results, and several examples are presented for simulations of the electromagnetic forming process to illustrate both the accuracy and high-efficiency parallel acceleration performance of the present program.

Multi-dimensional visual data completion via weighted hybrid graph-Laplacian

Graph Laplacian has been extensively studied as a priori to assist in visual completion of partially observed entries. However, a simple treatment using first-order operation may be suboptimal since the inherent smooth property within the visual data is ill-considered. To relieve this issue, in this paper, we propose an improved priori, namely weighted hybrid graph Laplacian (WHGL), to better characterize the spatially hierarchical knowledge hidden in the multi-dimensional data. On the whole, two individual layers are presented to encode the structured properties of the given data. In the first layer, the underlying adjacency matrix of the local and nonlocal relationships within the data is composed to capture the pixel-level structure information. Global smoothness is considered as the second layer to reveal the piecewise correlations of different visual elements. Mathematically, we generalize the graph Laplacian in the first- and second-order to construct the respective layers. Finally, a re-weighting strategy is used to mitigate the discontinuity of the reconstructed data under extreme cases, e.g. 95% voxels missing. With an initial step borrowed from some off-the-shelf methods, our model can be solved under limited iterations by the Preconditioned Conjugate Gradients algorithm. The experimental results demonstrate the superiority of the proposed method.

Three-Dimensional Direct Current Resistivity Finite Element Method Forward Modelling Based on Jacobi Preconditioned Conjugate Gradient Algorithm

The calculation speed and accuracy of resistivity forward model are important for direct current exploration. This study proposes a calculation method for the forward model. First, we structured the unstructured grid to enable dissection of the three-dimensional finite element space. We then solved the large sparse linear equations, using the preconditioned conjugate gradient method. The diagonal matrix in the Jacobi iteration was used as the pre-processed matrix, which had the characteristics of a convenient inversion solution and consumed less memory; therefore, it substantially increased the speed to solve large linear equations. The forward modelling performed using synthetic data shows that the preconditioned conjugate gradient algorithm exhibits fast calculation speed and high accuracy. These results are in accordance with the potential distribution rule of the underground electric field. The results have value for three-dimensional direct current exploration in the field of engineering. Furthermore, it provides a reference for solving large linear equations of the finite element method in other fields.

Two New Preconditioned Conjugate Gradient Methods for Minimization Problems

In application to general function, each of the conjugate gradient and Quasi-Newton methods has particular advantages and disadvantages. Conjugate gradient (CG) techniques are a class of unconstrained optimization algorithms with strong local and global convergence qualities and minimal memory needs. Quasi-Newton methods are reliable and efficient on a wide range of problems and they converge faster than the conjugate gradient method and require fewer function evaluations but they have the disadvantage of requiring substantially more storage and if the problem is ill-conditioned, they may take several iterations. A new class has been developed, termed preconditioned conjugate gradient (PCG) method. It is a method that combines two methods, conjugate gradient and Quasi-Newton. In this work, two new preconditioned conjugate gradient algorithms are proposed namely New PCG1 and New PCG2 to solve nonlinear unconstrained optimization problems. A new PCG1 combines conjugate gradient method Hestenes-Stiefel (HS) with new self-scaling symmetric Rank one (SR1), and a new PCG2 combines conjugate gradient method Hestenes-Stiefel (HS) with new self-scaling Davidon, Flecher and Powell (DFP). The algorithm uses the strong Wolfe line search condition. Numerical comparisons with standard preconditioned conjugate gradient algorithms show that for these new algorithms, computational scheme outperforms the preconditioned conjugate gradient.

Algorithm Research Based on the PDE Sensitivity Filter

The PDE sensitivity filter can eliminate the checkerboard patterns and numerical instability existing in the topology optimization results of continuum structures, and the essence of the PDE sensitivity filter is the Helmholtz partial differential equation with Neumann boundary conditions. To solve the large-scale PDE sensitivity filter problem, the conjugate gradient algorithm, the multigrid algorithm and the multigrid preconditioned conjugate gradient algorithm were used to solve the algebraic equations obtained by finite element analysis, and the effects of accuracy, filter radius and grid numbers on the efficiency of topology optimization were studied. The results show that, compared with the conjugate gradient algorithm and the multi-grid algorithm, the multi-grid preconditioned conjugate gradient algorithm has the least number of iterations and the shortest running time, which greatly improves the efficiency of topology optimization.

A cascade of preconditioned conjugate gradient networks for accelerated magnetic resonance imaging

Background and objective: Recent unfolding based compressed sensing magnetic resonance imaging (CS-MRI) methods only reinterpret conventional CS-MRI optimization algorithms and, consequently, inherit the weaknesses of the alternating optimization strategy. In order to avoid the structural complexity of the alternating optimization strategy and achieve better reconstruction performance, we propose to directly optimize the ℓ1 regularized convex optimization problem using a deep learning approach. Method: In order to achieve direct optimization, a system of equations solving the ℓ1 regularized optimization problem is constructed from the optimality conditions of a novel primal-dual form proposed for the effective training of the sparsifying transform. The optimal solution is obtained by a cascade of unfolding networks of the preconditioned conjugate gradient (PCG) algorithm trained to minimize the mean element-wise absolute difference (ℓ1 loss) between the terminal output and ground truth image in an end-to-end manner. The performance of the proposed method was compared with that of U-Net, PD-Net, ISTA-Net+, and the recently proposed projection-based cascaded U-Net, using single-coil knee MR images of the fastMRI dataset. Results: In our experiment, the proposed network outperformed existing unfolding-based networks and the complex version of U-Net in several subsampling scenarios. In particular, when using the random Cartesian subsampling mask with 25 % sampling rate, the proposed model outperformed PD-Net by 0.76 dB, ISTA-Net+ by 0.43 dB, and U-Net by 1.21 dB on the positron density without suppression (PD) dataset in term of peak signal to noise ratio. In comparison with the projection-based cascade U-Net, the proposed algorithm achieved approximately the same performance when the sampling rate was 25% with only 1.62% number of network parameters at the cost of a longer reconstruction time (approximately twice). Conclusion: A cascade of unfolding networks of the PCG algorithm was proposed to directly optimize the ℓ1 regularized CS-MRI optimization problem. The proposed network achieved improved reconstruction performance compared with U-Net, PD-Net, and ISTA-Net+, and achieved approximately the same performance as the projection-based cascaded U-Net while using significantly fewer network parameters.

Pipelined Preconditioned Conjugate Gradient Methods for real and complex linear systems for distributed memory architectures

Preconditioned Conjugate Gradient (PCG) is a popular method for solving large and sparse linear systems of equations. The performance of PCG at scale is affected due to the costly global synchronization steps that arise in dot-products on distributed memory systems. Pipelined PCG (PIPECG) removes the costly global synchronization steps from PCG by only executing a single non-blocking allreduce per iteration and overlapping it with independent computations.In our previous work, we have developed a novel pipelined PCG algorithm called PIPECG-OATI (One Allreduce per Two Iterations) for real linear systems which executes a single non-blocking allreduce per two iterations and provides a large overlap of global communication with independent computations at higher number of cores. Our method achieves this overlap by using iteration combination and by introducing new recurrence and non-recurrence computations. We implement optimizations in the PIPECG-OATI method to use cache memory efficiently.In this work, we present PIPECG-OATI-c method for linear systems with complex Hermitian positive definite and complex symmetric matrices. We compare our method with various pipelined CG methods on a variety of problems and demonstrate that our method always gives the least run times. We performed experiments with our method using 20M and 30M unknowns on up to 16K cores and obtained up to 2.48X performance improvement over PCG and 2.14X performance improvement over PIPECG methods. We also experimented with up to 1-billion unknowns on 16K cores, the largest problem size explored for the CG problem, to our knowledge, and obtained about 25% improvement over PCG.

Least-Squares Reverse Time Migration for Reflection-Angle-Dependent Reflectivity

Least-squares reverse time migration (LSRTM) can estimate high-quality reflectivity of subsurface medium from seismic data. However, the subsurface reflectivity depends on reflection angles, and its variations over reflection angles are extremely important because they can be used to estimate sub-surface physical properties for seismic interpretation. We present a new formulation of the LSRTM method that can estimate reflection-angle-dependent reflectivity from seismic data. We derive a forward modeling operator which predicts the reflection data without calculating the reflection angles, and verify that it approximately equals to the reflection-angle-dependent wave-equation-based Kirchhoff modeling operator under the assumption that the velocity perturbation is small and the reflection angle is smaller than the critical angle. Based on the proposed modeling operator associated with the adjoint of the angle-dependent wave-equation-based Kirchhoff modeling operator, we reformulate LSRTM as an inverse problem to invert for reflection-angle-dependent reflectivity using a preconditioned conjugate gradient algorithm. The algorithm uses a low-rank filter as the preconditioner to attenuate migration artifacts. Imaging tests on synthetic and field seismic data are used to verify validity and superiority of the proposed method. The tests illustrate that the proposed method can produce the reflection-angle-dependent reflectivity with much higher signal-to-noise ratio, resolution and amplitude fidelity than reverse time migration. Compared with conventional LSRTM, it can produce more focused stacked image when the migration velocity contains errors. Moreover, conventional LSRTM only produces the angle-independent reflectivity, whereas the proposed method has the feasibility to produce the reflection-angle-dependent reflectivity.

On numerical solution of nonlinear parabolic multicomponent diffusion-reaction problems

Abstract This work continues our previous analysis concerning the numerical solution of the multi-component mass transfer equations. The present test problems are two-dimensional, parabolic, non-linear, diffusion- reaction equations. An implicit finite difference method was used to discretize the mathematical model equations. The algorithm used to solve the non-linear system resulted for each time step is the modified Picard iteration. The numerical performances of the preconditioned conjugate gradient algorithms (BICGSTAB and GMRES) in solving the linear systems of the modified Picard iteration were analysed in detail. The numerical results obtained show good numerical performances.

Gravity Gradient Inversion of Gravity Field and Steady‐State Ocean Circulation Explorer Satellite Data for the Lithospheric Density Structure in the Qinghai‐Tibet Plateau Region and the Surrounding Regions

AbstractThe continuous uplift and evolution of the Qinghai‐Tibet Plateau region (TPR) represents one of the most important Cenozoic geological events on Earth. The deformation mechanism of the lithosphere in the TPR is still under debate. Relevant viewpoints can be informed by investigations of the lithospheric density structure in the Qinghai‐Tibet Plateau and surrounding regions (TP‐SRs). Considering investigating the fine‐scale lithospheric density structure of the large‐scale areas, the GOCE satellite data covering entire TP‐SRs are used in this study. GOCE gravity gradient data are first subjected to several corrections. Then, based on gravity gradient curvature (SI) calculation complemented with the preconditioned conjugate gradient algorithm, the distribution of SI and the lithospheric density structure at depths of 0–80 km in TP‐SRs are obtained. The effects of temperature differences and material composition on the lithospheric density distribution are further distinguished. The results show that in the crust, the density anomalies in the TPR are blocked by the surrounding hard blocks. In the upper mantle, the distribution of the density anomalies in the TPR is predominantly affected by the northward subduction of the Indian Plate and the high heat flux environment. The Longmen Shan Thrust Belt is found to have penetrated the lower crust but has not reached the depth of the upper mantle. The density anomalies in the basins on the northern side of the TPR exhibit features indicative of independent stable blocks.

Elastic Least-Squares Reverse-Time Migration Based on a Modified Acoustic-Elastic Coupled Equation for OBS Four-Component Data

Elastic least-squares reverse-time migration (ELSRTM) yields subsurface high-resolution elastic images from multicomponent seismic data. For ocean bottom seismic (OBS) acquisitions, both pressure and displacement components are measured using four-component (4C) modern receivers. However, these 4C data cannot be effectively simulated and migrated using conventional ELSRTM approaches due to the limitation of the standard elastic wave equation. For this reason, we propose a new ELSRTM method based on a modified acoustic-elastic coupled (AEC) equation. Owing to this equation, either the prediction or the migration of OBS 4C data can be implemented using a one-time wavefield simulation. Incorporation with a perturbation imaging condition, the proposed method can invert for P- and S-wave velocity and density and generate amplitude-preserving images. Besides, it can provide subsurface impedance reflectivity models by stacking the obtained multi-parameter images. In the least-squares migration, a preconditioned conjugate gradient algorithm is implemented based on a multi-parameter diagonal Hessian. Numerical experiments on an uncorrelated layer model and a portion of the Southern Yellow China sea model can demonstrate the effectiveness of this method on improving imaging resolution and accelerating convergence rate.

Three-dimensional forward modelling of gravity field vector and its gradient tensor using the compact difference schemes

SUMMARY The traditional gravity forward modelling methods for solving partial differential equations (PDEs) only can yield second-order accuracy. When computing the gravity field vector and gradient tensor from the obtained potential, those numerical differentiation approaches will inevitably lose accuracy. To mitigate this issue, we propose an efficient and accurate 3-D forward modelling algorithm based on a fourth-order compact difference scheme. First, a 19-point fourth-order compact difference scheme with general meshsizes in x-, y- and z-directions is adopted to discretize the governing 3-D Poisson’s equation. The resulting symmetric positive-definite linear systems are solved by the pre-conditioned conjugate gradient algorithm. To obtain the first-order (i.e. the gravity field vector) and second-order derivatives (i.e. the gravity gradient tensor) with fourth-order accuracy, we seek to solve a sequence of tridiagonal linear systems resulting from the above mentioned finite difference discretization by using fast Thomas algorithm. Finally, two synthetic models and a real topography relief are used to verify the accuracy of our method. Numerical results show that our method can yield a nearly fourth-order accurate approximation not only to the gravitational potential, but also to the gravity field vector and its gradient tensor, which clearly demonstrates its superiority over the traditional PDE-based methods.

Inversion of high-amplitude magnetic total field anomaly: an application to the Mengku iron-ore deposit, northwest China

In magnetic prospecting, the total field anomaly formula that represents the projection of the magnetic anomaly vector on the geomagnetic field is widely used because it simplifies the calculation of forward modelling and inversion of magnetic data. However, the projection anomaly yields errors relative to the true observed magnetic anomaly, especially for high-amplitude magnetic anomalies such as in iron orebody and unexploded ordnance prospecting. In this study, we analyse the difference between the projection anomaly and observed modulus difference anomaly with physical parameters, and propose to directly invert for the modulus difference anomaly by constructing a nonlinear matrix equation between the model corrections and data corrections. The inversion is then implemented using a preconditioned conjugate gradient algorithm. Synthetic and field magnetic data were used to test the inversion method. Comparison of the two types of total field anomalies shows that the error of the projection anomaly increased with increasing total-field magnetic anomaly. When the total-field magnetic anomaly was < 5,000 nT, the difference between the projection anomaly and modulus difference anomaly results can be ignored. For high-amplitude magnetic anomalies, the modulus difference anomaly inversion produced more accurate representations of both the shape and location of the magnetic sources.

Susceptibility and remanent magnetization inversion of magnetic data with a priori information of the Köenigsberger ratio

SUMMARY Magnetic susceptibility and natural remanent magnetization of rocks are useful parameters to study geological structures and geodynamic processes. Traditional widely used algorithms for the inversion of magnetic data can recover the distribution of the apparent susceptibility or total magnetization intensity, but do not provide information on the remanent magnetization. In this paper, we propose a framework to directly invert for the magnetic susceptibility and the natural remanent magnetization vector using surface or airborne magnetic data, assuming that the Köenigsberger ratio of the rock is known or approximately deducible. The susceptibility and remanence are computed using two different approaches: (1) the susceptibility, intensity, and direction of the remanent magnetization are continuously recovered for each discretized cell and (2) the remanence direction is assumed to be uniform in each subzone and is iteratively computed as discrete values. Both processes are implemented using the preconditioned conjugate gradient algorithm. The method is tested on three synthetic models and one field data set from the Zaohuohexi iron-ore deposit, Qinghai Province, northwest China. The results of the continuous inversion show the trend of the remanent magnetization directions, while the discrete inversion yields more specific values. This inversion framework can determine the source bodies’ geometry and position, and also provide superposed and comprehensive information on the natural remanent magnetization, which may be useful to investigate geological bodies bearing stable primary remanent magnetization.

Convolution operators for visual tracking based on spatial–temporal regularization

In recent years, the method based on discriminative correlation filter has been shown excellent performance in short-term visual tracking. However, discriminative correlation filter-based method heavily suffers from the problem of the multiple peaks and model drift in responds maps incurred by occlusion and rotation. To solve the above problem, we proposed convolution operators for visual tracking based on spatial–temporal regularization. Firstly, we add spatial–temporal regularization in loss function, which will guarantee continuity of the model in time. And we use preconditioned conjugate gradient algorithm to obtain filter coefficients. Secondly, we proposed channel reliability to estimate quality of the learned filter and fuse the different reliability coefficients to weight response map in location. We set a threshold to reduce the number of iteration in location and accelerate the compute speed of algorithm. Finally, we use two different correlation filters to estimate location and scale of target, respectively. Extensively experiment in five video sequences show that our tracker has been significantly improved performance in case of occlusion and rotation. The AUC in success plot improves 33.2% than ECO-HC and 41.5% than STRCF, respectively.

Application of Tikhonov Regularization for 1D Geothermal Heat Flux Ill-Posed Inverse Problem: A Case Study on Chad Sedimentary Basin

Subsurface heat flux information is important in geothermal exploration. With the information, geophysicists can map exactly the thermal potential in a particular area. Based on the surface heat flux, inverse modeling produces the 1D subsurface heat flux distribution. However, inverse problems in the geothermal system are generally ill-posed. Small changes in the data can cause large changes in the solution and the solution may not be unique. To solve the mentioned non-linear and ill-posed equation above, Tikhonov regularization is a choice for stabilizing the inverse calculation. This paper demonstrates how Tikhonov regularization is useful to solve subsurface heat flux distribution both in the synthetic model and real model. Based on surface heat flux distribution from the direct problem, the preconditioned conjugate gradient algorithm calculates the subsurface heat flux. With the correct choice of the regularization parameter, the inverse model fits the initial model. For the testing purposes in real-world conditions, Chad sedimentary basin located in Chad and Nigeria is used as a model. A high geothermal gradient is found in this area. Therefore, geothermal explorations are on the rise recently. Its thermal conductivity, heat production, and stratigraphy data from previous researches provide information about the initial model. The heat flux curve generated from inversion matches the initial noisy model with the error of around 10−9 mW/m2. Therefore, to answer the increasing energy demand, this method can be highly applicable to future geothermal prospecting.