Stationary Iterative Methods Research Articles

A new block preconditioner for weighted Toeplitz regularized least-squares problems

We introduce a new block preconditioner for the solution of weighted Toeplitz regularized least-squares problems written in augmented system form. The proposed preconditioner is obtained based on the new splitting of coefficient matrix which results in an unconditionally convergent stationary iterative method. Spectral analysis of the preconditioned matrix is investigated. In particular, we show that the preconditioned matrix has a very nice eigenvalue distribution which can lead to fast convergence of the preconditioned Krylov subspace methods such as GMRES. Numerical experiments are reported to demonstrate the performance of preconditioner used with GMRES method in the solution of augmented system form of weighted Toeplitz regularized least-squares problems.

Depth-dependent attenuation and backscattering characterization of optical coherence tomography by stationary iterative method.

Extracting optical properties of tissue [e.g., the attenuation coefficient () and the backscattering fraction] from the optical coherence tomography (OCT) images is a valuable tool for parametric imaging and related diagnostic applications. Previous attenuation estimation models depend on the assumption of the uniformity of the backscattering fraction () within layers or whole samples, which does not accurately represent real-world conditions. Our aim is to develop a robust and accurate model that calculates depth-wise values of attenuation and backscattering fractions simultaneously from OCT signals. Furthermore, we aim to develop an attenuation compensation model for OCT images that utilizes the optical properties we obtained to improve the visual representation of tissues. Using the stationary iteration method under suitable constraint conditions, we derived the approximated solutions of and on a single scattering model. During the iteration, the estimated value of can be rectified by introducing the large variations of , whereas the small ones were automatically ignored. Based on the calculation of the structure information, the OCT intensity with attenuation compensation was deduced and compared with the original OCT profiles. The preliminary validation was performed in the OCT A-line simulation and Monte Carlo modeling, and the subsequent experiment was conducted on multi-layer silicone-dye- phantoms and ex vivo cow eyes. Our method achieved robust and precise estimation of and for both simulated and experimental data. Moreover, corresponding OCT images with attenuation compensation provided an improved resolution over the entire imaging range. Our proposed method was able to correct the estimation bias induced by the variations of and provided accurate depth-resolved measurements of both and simultaneously. The method does not require prior knowledge of the morphological information of tissue and represents more real-life tissues. Thus, it has the potential to help OCT imaging based disease diagnosis of complex and multi-layer biological tissue.

A New Quasi-Newton Method with PCG Method for Nonlinear Optimization Problems

The major stationary iterative method used to solve nonlinear optimization problems is the quasi-Newton (QN) method. Symmetric Rank-One (SR1) is a method in the quasi-Newton family. This algorithm converges towards the true Hessian fast and has computational advantages for sparse or partially separable problems [1]. Thus, investigating the efficiency of the SR1 algorithm is significant. It's possible that the matrix generated by SR1 update won't always be positive. The denominator may also vanish or become zero. To overcome the drawbacks of the SR1 method, resulting in better performance than the standard SR1 method, in this work, we derive a new vector <img src=image/13429423_01.gif> depending on the Barzilai-Borwein step size to obtain a new SR1 method. Then using this updating formula with preconditioning conjugate gradient (PCG) method is presented. With the aid of inexact line search procedure by strong Wolfe conditions, the new SR1 method is proposed and its performance is evaluated in comparison to the conventional SR1 method. It is proven that the updated matrix of the new SR1 method, <img src=image/13429423_02.gif>, is symmetric matrix and positive definite matrix, given <img src=image/13429423_03.gif> is initialized to identity matrix. In this study, the proposed method solved 13 problems effectively in terms of the number of iterations (NI) and the number of function evaluations (NF). Regarding NF, the new SR1 method also outperformed the classic SR1 method. The proposed method is shown to be more efficient in solving relatively large-scale problems (5,000 variables) compared to the original method. From the numerical results, the proposed method turned out to be significantly faster, effective and suitable for solving large dimension nonlinear equations.

Fourier Neural Solver for Large Sparse Linear Algebraic Systems

Large sparse linear algebraic systems can be found in a variety of scientific and engineering fields and many scientists strive to solve them in an efficient and robust manner. In this paper, we propose an interpretable neural solver, the Fourier neural solver (FNS), to address them. FNS is based on deep learning and a fast Fourier transform. Because the error between the iterative solution and the ground truth involves a wide range of frequency modes, the FNS combines a stationary iterative method and frequency space correction to eliminate different components of the error. Local Fourier analysis shows that the FNS can pick up on the error components in frequency space that are challenging to eliminate with stationary methods. Numerical experiments on the anisotropic diffusion equation, convection–diffusion equation, and Helmholtz equation show that the FNS is more efficient and more robust than the state-of-the-art neural solver.

Relaxation parameters and composite refinement techniques

A composite refinement technique for two stationary iterative methods, one of them contains a relaxation parameter, is introduced. Four new techniques, Jacobi successive over relaxation (SOR) composite refinement (RJSOR), SOR Jacobi composite refinement (RSORJ), Gauss–Seidel (GS) SOR composite refinement (RGSSOR) and SOR with GS composite refinement (RSORGS) are compared with their classical forms. The efficient performance of the new forms is well established and confirmed through numerical example. The computational costs and the speed of convergence are considered. The decrease in the required number of iteration is established through the calculation of the spectral radius of the iteration matrices. It is illustrated that the convergence of Jacobi and Gauss–Seidel methods engage the divergence and extend the domain of convergence in the SOR method in the refinement technique. The calculations and graphs are performed by computer algebra system, Mathematica.

Numerical investigation on evaluating nonlinear waves due to an air cushion vehicle in steady motion by a higher order desingularized boundary integral equation method

The problem of nonlinear ship waves generated by an air cushion vehicle is studied and modelled in the framework of potential flow theory. A higher order desingularized boundary integral equation method (HODBIEM) using Rankine source is developed. A stationary iterative method is adopted to solve the nonlinear free surface which is discretized into the inner area disturbed by cushion pressure and the outer undisturbed area. Techniques for numerical implementation are described in detail. To save computing time, grid density and domain size are refined.Both straight-ahead and yawed cases are considered for validations. Wave-making resistance versus various speeds in deep and shallow water are calculated for straight-ahead vehicles. The resistance curves obtained by HODBIEM agree well with results of nonlinear CFD simulations. The linear theory predicts higher peak value of resistance than nonlinear methods. The experimental data fit at low Froude number but deviate from numerical predictions when exceeding primary hump speed. Effects of speed and water depth are analyzed by wave patterns and wave profiles. Resistance, side force and yaw moment are calculated at Fn = 0.533 for drift angle from 0° to 90°. Satisfactory results are given by present method. Details about wave patterns and wave angles are discussed.

Composite refinement techniques for solving linear systems

A composite refinement approach for stationary iterative methods is introduced. Two new formulas (RJGS and RGSJ) are compared with the classical forms. Rates of convergence of the introduced composite formulas (RJGS and RGSJ) are well established. The efficient performance of the new forms is established theoretically and confirmed through numerical examples. The decrease in the required number of iterations for convergence is established through the calculation of the spectral radius of the iteration matrices. The algorithmic structure of the new formulas is announced. Three numerical examples with different convergent properties are considered. The calculations are performed with the help of computer algebra software Mathematica.

Genetic programming for iterative numerical methods

We introduce GPLS (Genetic Programming for Linear Systems) as a GP system that finds mathematical expressions defining an iteration matrix. Stationary iterative methods use this iteration matrix to solve a system of linear equations numerically. GPLS aims at finding iteration matrices with a low spectral radius and a high sparsity, since these properties ensure a fast error reduction of the numerical solution method and enable the efficient implementation of the methods on parallel computer architectures. We study GPLS for various types of system matrices and find that it easily outperforms classical approaches like the Gauss–Seidel and Jacobi methods. GPLS not only finds iteration matrices for linear systems with a much lower spectral radius, but also iteration matrices for problems where classical approaches fail. Additionally, solutions found by GPLS for small problem instances show also good performance for larger instances of the same problem.

Numerical solutions to linear transfer problems of polarized radiation

Context. Numerical solutions to transfer problems of polarized radiation in solar and stellar atmospheres commonly rely on stationary iterative methods, which often perform poorly when applied to large problems. In recent times, stationary iterative methods have been replaced by state-of-the-art preconditioned Krylov iterative methods for many applications. However, a general description and a convergence analysis of Krylov methods in the polarized radiative transfer context are still lacking. Aims. We describe the practical application of preconditioned Krylov methods to linear transfer problems of polarized radiation, possibly in a matrix-free context. The main aim is to clarify the advantages and drawbacks of various Krylov accelerators with respect to stationary iterative methods and direct solution strategies. Methods. After a brief introduction to the concept of Krylov methods, we report the convergence rate and the run time of various Krylov-accelerated techniques combined with different formal solvers when applied to a 1D benchmark transfer problem of polarized radiation. In particular, we analyze the GMRES, BICGSTAB, and CGS Krylov methods, preconditioned with Jacobi, (S)SOR, or an incomplete LU factorization. Furthermore, specific numerical tests were performed to study the robustness of the various methods as the problem size grew. Results. Krylov methods accelerate the convergence, reduce the run time, and improve the robustness (with respect to the problem size) of standard stationary iterative methods. Jacobi-preconditioned Krylov methods outperform SOR-preconditioned stationary iterations in all respects. In particular, the Jacobi-GMRES method offers the best overall performance for the problem setting in use. Conclusions. Krylov methods can be more challenging to implement than stationary iterative methods. However, an algebraic formulation of the radiative transfer problem allows one to apply and study Krylov acceleration strategies with little effort. Furthermore, many available numerical libraries implement matrix-free Krylov routines, enabling an almost effortless transition to Krylov methods.

Numerical solutions to linear transfer problems of polarized radiation

Context.The numerical modeling of the generation and transfer of polarized radiation is a key task in solar and stellar physics research and has led to a relevant class of discrete problems that can be reframed as linear systems. In order to solve such problems, it is common to rely on efficient stationary iterative methods. However, the convergence properties of these methods are problem-dependent, and a rigorous investigation of their convergence conditions, when applied to transfer problems of polarized radiation, is still lacking.Aims.After summarizing the most widely employed iterative methods used in the numerical transfer of polarized radiation, this article aims to clarify how the convergence of these methods depends on different design elements, such as the choice of the formal solver, the discretization of the problem, or the use of damping factors. The main goal is to highlight advantages and disadvantages of the different iterative methods in terms of stability and rate of convergence.Methods.We first introduce an algebraic formulation of the radiative transfer problem. This formulation allows us to explicitly assemble the iteration matrices arising from different stationary iterative methods, compute their spectral radii and derive their convergence rates, and test the impact of different discretization settings, problem parameters, and damping factors.Results.Numerical analysis shows that the choice of the formal solver significantly affects, and can even prevent, the convergence of an iterative method. Moreover, the use of a suitable damping factor can both enforce stability and increase the convergence rate.Conclusions.The general methodology used in this article, based on a fully algebraic formulation of linear transfer problems of polarized radiation, provides useful estimates of the convergence rates of various iterative schemes. Additionally, it can lead to novel solution approaches as well as analyses for a wider range of settings, including the unpolarized case.

A new iterative method for solving a class of two-by-two block complex linear systems

We present a stationary iteration method, namely Alternating Symmetric positive definite and Scaled symmetric positive semidefinite Splitting (ASSS), for solving the system of linear equations obtained by using finite element discretization of a distributed optimal control problem together with time-periodic parabolic equations. An upper bound for the spectral radius of the iteration method is given which is always less than 1. So convergence of the ASSS iteration method is guaranteed. The induced ASSS preconditioner is applied to accelerate the convergence speed of the GMRES method for solving the system. Numerical results are presented to demonstrate the effectiveness of both the ASSS iteration method and the ASSS preconditioner.

Refinement of Extended Accelerated Over-Relaxation Method for Solution of Linear Systems

Given any linear stationary iterative methods in the form z^(i+1)=Jz^(i)+f, where J is the iteration matrix, a significant improvements of the iteration matrix will decrease the spectral radius and enhances the rate of convergence of the particular method while solving system of linear equations in the form Az=b. This motivates us to refine the Extended Accelerated Over-Relaxation (EAOR) method called Refinement of Extended Accelerated Over-Relaxation (REAOR) so as to accelerate the convergence rate of the method. In this paper, a refinement of Extended Accelerated Over-Relaxation method that would minimize the spectral radius, when compared to EAOR method, is proposed. The method is a 3-parameter generalization of the refinement of Accelerated Over-Relaxation (RAOR) method, refinement of Successive Over-Relaxation (RSOR) method, refinement of Gauss-Seidel (RGS) method and refinement of Jacobi (RJ) method. We investigated the convergence of the method for weak irreducible diagonally dominant matrix, matrix or matrix and presented some numerical examples to check the performance of the method. The results indicate the superiority of the method over some existing methods.

Digit Stability Inference for Iterative Methods Using Redundant Number Representation

In our recent work on iterative computation in hardware, we showed that arbitrary-precision solvers can perform more favorably than their traditional arithmetic equivalents when the latter's precisions are either under- or over-budgeted for the solution of the problem at hand. Significant proportions of these performance improvements stem from the ability to infer the existence of identical most-significant digits between iterations. This technique uses properties of algorithms operating on redundantly represented numbers to allow the generation of those digits to be skipped, increasing efficiency. It is unable, however, to guarantee that digits will stabilize, i.e., never change in any future iteration. In this article, we address this shortcoming, using interval and forward error analyses to prove that digits of high significance will become stable when computing the approximants of systems of linear equations using stationary iterative methods. We formalize the relationship between matrix conditioning and the rate of growth in most-significant digit stability, using this information to converge to our desired results more quickly. Versus our previous work, an exemplary hardware realization of this new technique achieves an up-to 2.2× speedup in the solution of a set of variously conditioned systems using the Jacobi method.

Alternating stationary iterative methods based on double splittings

Matrix double splitting iterations are simple in implementation while solving real non-singular (rectangular) linear systems. In this paper, we present two Alternating Double Splitting (ADS) schemes formulated by two double splittings and then alternating the respective iterations. The convergence conditions are then discussed along with comparative analysis. The set of double splittings used in each ADS scheme induces a preconditioned system which helps in showing the convergence of the ADS schemes. We also show that the classes of matrices for which one ADS scheme is better than the other, are mutually exclusive. Numerical experiments confirm that the proposed ADS schemes have several computational advantages over the existing methods. Though the problems are considered in the rectangular matrix settings, the same problems are even new in non-singular matrix settings.

Solvability for Two Forms of Nonlinear Matrix Equations

In this paper, we study nonlinear matrix equations $$\begin{aligned} X^p=A+\sum \limits _{i=1}^m M_i^T(X\#B)M_i \end{aligned}$$ and $$\begin{aligned} X^p=A+\sum \limits _{i=1}^j M_i^T(X\#B)M_i+\sum \limits _{i=j+1}^m M_i^T(X^{-1}\#B)M_i, \end{aligned}$$ where p, m, j are positive integers, $$1\le j\le m$$ , A, B are $$n\times n$$ positive definite matrices and $$M_i(i=1,2,3,\ldots ,m)$$ are $$n\times n$$ nonsingular real matrices. Based on some fixed point theorems for monotone and mixed monotone operators in ordered Banach spaces and some properties of cone, we prove that these equations always have a unique positive definite solution. In addition, an iterative sequence can be given to approximate the unique positive definite solution by employing a multi-step stationary iterative method.

Incomplete block‐matrix factorization of M‐matrices using two‐step iterative method for matrix inversion and preconditioning

Using the general method of Owe Axelsson given in 1986 for incomplete factorization of M‐matrices in block‐matrix form, we give a recursive approach to construct incomplete block‐matrix factorization of M‐matrices by proposing a two‐step iterative method for the approximation of the inverse of diagonal pivoting block matrices at each stage of the recursion. For various predescribed tolerances in the accuracy of the approximation of the inverses, the obtained incomplete block‐matrix factorizations are used to precondition the iterative methods as one‐step stationary iterative (OSSI) method and biconjugate gradient stabilized method (BI‐CGSTAB). Certain applications are conducted on M‐matrices occurring from the discretization of two Dirichlet boundary value problems of Laplace's equation on a rectangle using finite difference method. Numerical results justify that the given incomplete block‐matrix factorization of M‐matrices using the two‐step iterative method to approximate the inverse of diagonal pivoting block matrices at each stage of the recursion give robust preconditioners. The obtained results are presented through tables and figures.

Regularized DPSS preconditioners for generalized saddle point linear systems

By introducing a regularization matrix and an additional iteration parameter, a new class of regularized deteriorated positive-definite and skew-Hermitian splitting (RDPSS) preconditioners are proposed for generalized saddle point linear systems. Compared with the well-known Hermitian and skew-Hermitian splitting (HSS) preconditioner and the regularized HSS preconditioner (Bai, 2019) studied recently, the new RDPSS preconditioners have much better computing efficiency especially when the (1,1) block matrix is non-Hermitian. It is proved that the corresponding RDPSS stationary iteration method is unconditionally convergent. In addition, clustering property of the eigenvalues of the RDPSS preconditioned matrix is studied in detail. Two numerical experiments arising from the meshfree discretization of a static piezoelectric equation and the finite element discretization of the Navier–Stokes equation show the effectiveness of the new proposed preconditioners.

Regularized preconditioned GMRES and the regularized iteration method

The present paper proposes a new preconditioner for solving saddle point problems. The preconditioner is obtained by replacing the (1,2) and (2,2) blocks in the original saddle point matrix A by another well chosen block. The proposed preconditioner can be used as a preconditioner corresponding to the stationary iterative method or to accelerate the convergence of the generalized minimal residual method (GMRES). The convergent and semi-convergent analyses of the regularized iteration method are presented. The eigenvalue distribution and the forms of the eigenvectors of the preconditioned matrix are analyzed. Finally, numerical results show the effectiveness of the proposed preconditioner as compared to other preconditioners.

A Fast Block $\alpha$-Circulant Preconditoner for All-at-Once Systems From Wave Equations

In this paper, we propose a fast block $\alpha$-circulant preconditioner for solving the nonsymmetric linear system arising from an all-at-once implicit discretization scheme in time for the wave equation. As a generalization of the well-known block circulant preconditioning technique, the proposed block $\alpha$-circulant preconditioner can also be efficiently inverted in a parallel-in-time manner. The complex eigenvalues of the preconditioned matrix are fully derived in explicit expression and its diagonalizability is also shown. Furthermore, a mesh-independent convergence rate of the preconditioned GMRES method is derived under certain conditions. Building upon the proposed preconditioner, two stationary iterative methods with uniform asymptotic convergence rates were also presented. The extension of our preconditioner within a simplified Newton iteration to the nonlinear wave equation is also discussed. Both linear and nonlinear numerical examples are given to illustrate the promising performance of our proposed block $\alpha$-circulant preconditioner and also validate our theoretical results on convergence analysis.

Steady flow analysis of a slender wing by lifting surface method

In the paper hereby, steady flow around a thin-walled wing is analysed by means of the Lifting Surface Method. In order to carry out tests, the wing has been divided into a finite number of quadrilateral panels. All panel edges in turn are replaced by discrete straight vortex segments which induce velocities within the flow field. The problem boils down to working out velocity circulation distribution on the wing surface. For this purpose, numerical realization has been developed in C by Minimalist GNU for Windows compiler and Code::Blocks IDE. To work out a solution to the linear non-homogeneous algebraic system, the Gauss – Seidel stationary iterative method has been applied. The obtained results for various angle of attack values are depicted by means of ParaView.