We present a convergence analysis of a Galerkin-type two-grid method based on plain (or unsmoothed) aggregation. The analysis assumes that the system matrix is symmetric nonnegative definite (SNND) and can be decomposed into (or bounded below by) a sum of SNND matrices, each associated with an individual aggregate. Then, the convergence factor of the considered two-grid method can be bounded above (up to the smoother's scaling) as a function of a maximum over all aggregates of a quantity associated to each aggregate and named aggregate's quality. The presented analysis is a generalization of the analysis from [A. Napov and Y. Notay, Numer. Linear Algebra Appl. 18(2011), pp. 539-564], and is suitable for finite element discretizations, for which the system matrix subdivision into element matrices is a natural starting point. The analysis is further applied in the context of standard finite element discretizations of linear elasticity problems. More specifically, we derive a sharp estimate of the quality parameter of a two-node aggregate, that is, an aggregate containing all unknowns associated to two connected nodes of the discretization mesh. This estimate is given here as a function of the average stiffness tensor (which we do not assume isotropic), the vector connecting the nodes, the geometry of the elements containing both nodes and the diagonal entries of the system matrix corresponding to the nodes. Eventually, we illustrate the accuracy and usefulness of the presented bounds with some numerical experiments.
We consider the numerical solution of discrete Oseen problems. We focus on the recently proposed transform-then-solve approach, which amounts to first applying a specific algebraic transformation to the linear system of equations arising from the discretization and then solving the transformed system with an algebraic multigrid method. Promising results have been previously obtained with a two-grid variant, and here we bring two key improvements to make the approach robust in a multilevel setting. For a model problem with a constant convection field, it is shown, in a local Fourier analysis setting, that the two-grid method is convergent at any level of the hierarchy, with bounds that are independent of both the mesh size and the Reynolds number. Numerical results with a K-cycle multigrid scheme confirm the theoretical expectations for constant coefficient problems. For problems with variable convective flow, the method appears also robust with grid independence convergence, although a mild dependency with respect to the Reynolds number shows up at very low viscosity in driven cavity problems. The method is based on point coarsening and therefore, besides the system matrix, requires knowing, for each discrete unknown, which grid point it is associated with.
We consider numerical solution of finite element discretizations of the Stokes problem. We focus on the transform-then-solve approach, which amounts to first apply a specific algebraic transformation to the linear system of equations arising from the discretization, and then solve the transformed system with an algebraic multigrid method. The approach has recently been applied to finite difference discretizations of the Stokes problem with constant viscosity, and has recommended itself as a robust and competitive solution method. In this work, we examine the extension of the approach to standard finite element discretizations of the Stokes problem, including problems with variable viscosity. The extension relies, on one hand, on the use of the successive over-relaxation method as a multigrid smoother for some finite element schemes. On the other hand, we present strategies that allow us to limit the complexity increase induced by the transformation. Numerical experiments show that for stationary problems our method is competitive compared to a reference solver based on a block diagonal preconditioner and MINRES, and suggest that the transform-then-solve approach is also more robust. In particular, for problems with variable viscosity, the transform-then-solve approach demonstrates significant speed-up with respect to the block diagonal preconditioner. The method is also particularly robust for time-dependent problems whatever the time step size.
We propose a new sparse matrix format which captures the matrix structure typical for discretized partial differential equations with piecewise-constant coefficients. The format uses a stencil representation for some blocks of matrix rows, typically corresponding to regions with constant coefficients, whereas other rows are encoded in the general compressed sparse row format. The stencil representation saves memory and is suitable for SIMD-like parallelism as available on GPUs. Further, this format is well suited for the implementation of algebraic multigrid methods, and we present a proof-of-concept GPU-accelerated aggregation-based algebraic multigrid solver based on this format. This solver is compared on a few model problems (2-dimension and 3-dimension Poisson-like) with the compressed sparse row-based solver AmgX from NVIDIA and with the CUDA version of the BoomerAMG solver. For the considered tested problems with one million unknowns or more, the presented solver outperforms AmgX and BoomerAMG in terms of both run time and memory usage, and the performance gap increases with the system size.
We consider the numerical solution of discrete Oseen problems. We propose a new approach that consists of first applying a simple algebraic transformation to the linear system, which is afterwards preconditioned with an aggregation-based algebraic two-grid method. An algebraic analysis is provided, which proves uniform convergence in norm with respect to problem parameters if a few constants can be uniformly bounded. A further analysis of these constants shows that they can be bounded in the case of a constant convection field, provided that the coarsening of the pressure unknowns is also driven by the convection field. This makes the method essentially different from a similar method developed for Stokes equations which initially inspired the present work. Technically, this means that one has to either use point-based coarsening, or that an auxiliary convection-diffusion matrix has to be built on the pressure space to guide the coarsening, which makes the method only semialgebraic. Using this ingredient, promising results are obtained, showing that the number of iterations is, in practice, uniformly bounded with respect to both the mesh size and the Reynolds number even for challenging recirculating convection fields or in presence of outflow boundary conditions.
Space-time multigrid refers to the use of multigrid methods to solve discretized partial differential equations considering at once multiple time steps. A new theoretical analysis is developed for the case where one uses coarsening in space only. It proves bounds on the 2-norm of the iteration matrix that connect it to the norm of the iteration matrix when using the same multigrid method to solve the corresponding stationary problem. When using properly defined wavefront type smoothers, the bound is uniform with respect to the mesh size, the time step size, and the number of time steps, and addresses both the two-grid case and the W-cycle. On the other hand, for time-parallel smoothers, the results clearly show the condition to be satisfied by the time step size to have similar performance as with wavefront type smoothers. The analysis also leads to the definition of an effective smoothing factor that allows one to quickly check the potentialities of a given smoothing scheme. The accuracy of the theoretical estimates is illustrated on a numerical example, highlighting the relevance of the effective smoothing factor and the usefulness in following the provided guidelines to have robustness with respect to the time step size.
We consider the iterative solution of linear systems with a symmetric saddle point system matrix. We address a family of solution techniques that exploit the knowledge of a preconditioner (or approximate solution procedure) both for the top left block of the matrix on the one hand and for the Schur complement resulting from its elimination on the other hand. This includes many "segregated" or "Schur complement" iterations such as the inexact Uzawa method and pressure correction techniques, and also many "block" preconditioners, based on the approximate block factorization of the system matrix. An analysis is developed which proves convergence in norm of stationary iterations. It is more rigorous than eigenvalue analyses which ignore nonnormality effects, while being more general than previous norm analyses. The analysis also clarifies the relations that exist between the many members of this family of methods and offers practical guidelines to select the scheme most appropriate to a situation at hand.
Core results about the algebraic analysis of two-grid methods are extended in relations bounding the field of values (or numerical range) of the iteration matrix. On this basis, bounds are obtained on its norm and numerical radius, leading to rigorous convergence estimates. Numerical illustrations show that the theoretical results deliver qualitatively good predictions, allowing one to anticipate success or failure of the two-grid method. They also indicate that the field of values and the associated numerical radius are much more reliable convergence indicators than the eigenvalue distribution and the associated spectral radius. On this basis, some discussion is developed about the role of local Fourier or local mode analysis for nonsymmetric problems.
A method is investigated for solving stationary or time-dependent discrete Stokes equations. It uses one of the standard flavors of algebraic multigrid for coupled partial differential equations, which, however, is not applied directly to the linear system stemming from discretization, but to an equivalent system obtained with a simple algebraic transformation (which may be seen as a form of preconditioning in the literal sense). A two-grid analysis is provided, showing that the eigenvalues of the preconditioned matrix are within a region of the complex plane that is both bounded and away from the origin, independently of the mesh or grid size, as well as of other main problem parameters. On the other hand, whereas the approach can in principle be combined with any type of algebraic multigrid scheme, an investigation of the properties of the coarse grid matrices reveals that plain aggregation has to be preferred to maintain nice two-grid convergence at coarser levels. Eventually, numerical experiments are reported showing that the resulting method is both robust and cost effective, being significantly faster than a state-of-the-art competitor which combines MINRES with optimal block diagonal preconditioning.
Standard discretizations of Stokes problems lead to linear systems of equations in saddle point form, making difficult the application of algebraic multigrid methods. In this paper, a new approach is proposed. It consists in first transforming the system by pre- and post-multiplication with simple, algebraic, sparse block triangular matrices. This is a form of pre-conditioning in the literal sense, designed to make sure that the transformed matrix is well adapted to multigrid. In particular, after transformation, all the diagonal blocks are symmetric and positive definite, and correspond to, or resemble, a discrete Laplace operator. Then, to each of these diagonal blocks is associated a prolongation that works well for it, using any relevant algebraic or geometric multigrid method. Next, a multigrid scheme for the global system is naturally set up by combining these partial prolongations with a Galerkin coarse grid matrix. For this approach combined with damped Jacobi-smoothing, a uniform two-grid convergence bound is derived for the global system under the assumption that the two-grid schemes for the different diagonal blocks are themselves uniformly convergent. This result is illustrated by a few examples, showing further that time-dependent problems and variable viscosity can be handled in a natural way, without requiring parameter adjustment. A numerical comparison also shows that the new approach can be more effective than state-of-the-art block preconditioning techniques.
We consider the algebraic convergence theory that gives its theoretical foundations to classical algebraic multigrid methods. All the main results constitutive of the approach are properly extended to singular compatible systems, including the recent sharp convergence estimates for both symmetric and nonsymmetric systems. In fact, all results carry over smoothly to the singular case, which therefore does not require a specific treatment (except a proper handling of issues associated with singular coarse grid matrices). Regarding problems with a low-dimensional null space, the presented results thus mainly confirm what has been observed often at a more practical level in many previous works. As illustrated by the discussion of the application to the curl-curl equation, the potential impact is greater for problems with large-dimensional null space. Indeed, it turns out that the design of multilevel methods can then actually be easier for the singular case than it is for nearby regularized problems.
We consider the iterative solution of linear systems whose matrices are Laplacians of undirected graphs. Designing robust solvers for this class of problems is challenging due to the diversity of connectivity patterns encountered in practical applications. Our starting point is a recently proposed aggregation-based algebraic multigrid method that combines the recursive static elimination of the vertices of degree 1 with the degree-aware rooted aggregation (DRA) algorithm. The latter always produces aggregates big enough to ensure that the preconditioner cost per iteration is low. Here we further improve the robustness of the method by controlling the quality of the aggregates. More precisely, "bad" vertices are removed from the aggregates formed by the DRA algorithm until a quality test is passed. This ensures that the two-grid condition number is nicely bounded, whereas the cost per iteration is kept low by reforming too small aggregates when it happens that the mean aggregate size is not large enough. The effectiveness and the robustness of the resulting method are assessed on a large set of undirected graphs by comparing with the variant without quality control, as well as with another state-of-the art graph Laplacian solver.
The paper considers the parallel implementation of an algebraic multigrid method. The sequential version is well suited to solve linear systems arising from the discretization of scalar elliptic PDEs. It is scalable in the sense that the time needed to solve a system is (under known conditions) proportional to the number of unknowns. The associate software code is also robust and often significantly faster than other algebraic multigrid solvers. The present work addresses the challenge of porting it on massively parallel computers. In this view, some critical components are redesigned, in a relatively simple yet not straightforward way. Thanks to this, excellent weak scalability results are obtained on three petascale machines among the most powerful today available.
This special issue of Computing and Visualization in Science contains selected papers from the 2014 European Multigrid conference (EMG 2014), which took place in Leuven, Belgium, from 9 to 12 September, 2014.European Multigrid (EMG) is a series of conferences on the theme of multigrid methods and related fields.EMG is one of the most important conference series on this topic worldwide.Previous EMG meetings have been held in Cologne (1981 and 1985), Bonn (1990), Amsterdam (1993), Stuttgart (1996), Gent (1999), Hohenwart (2002), Scheveningen (2005), Bad Herrenhalb (2008), Ischia (2010) and Schwetzingen (2012).The EMG 2014 conference took place at the Irish College, in Leuven, an old university college founded in 1607.It attracted more than 80 participants of 15 countries, of which 20 were PhD-students; 12 participants attended the conference from outside Europe.Further details are available on the conference homepage http://metronu.ulb.ac.be/EMG2014/
Many simulation codes in physics or engineering require the repeated solution of large linear systems stemming from (or closely related to) the discretization of scalar elliptic PDEs. For large 3D simulations, it is nowadays standard to use multigrid methods. These methods are indeed scalable in the sense that the overall computational work to obtain the solution up to a prescribed tolerance is proportional to the number of unknowns. However, developing a truly efficient implementation for massively parallel computers is still considered a challenge. It is especially the case for algebraic multigrid schemes, whose setup stage requires only the system matrix, and that are therefore not limited to certain types of discretization. In this talk we consider in particular a nowadays popular aggregation-based algebraic multigrid method, as implemented in the AGMG software code (which has several hundreds of users in both academia and industry). To improve performance on massive parallel systems, some critical algorithmic components have been redesigned, in a relatively simple yet not straightforward way. Thanks to this, excellent weak scalability results have been obtained on some of the Europe’s top supercomputers.
About thirty years ago, Achi Brandt wrote a seminal paper providing a convergence theory for algebraic multigrid methods [Appl. Math. Comput., 19 (1986), pp. 23-56]. Since then, this theory has been improved and extended in a number of ways, and these results have been used in many works to analyze algebraic multigrid methods and guide their developments. This paper makes a concise exposition of the state of the art. Results for symmetric and nonsymmetric matrices are presented in a unified way, highlighting the influence of the smoothing scheme on the convergence estimates. Attention is also paid to sharp eigenvalue bounds for the case where one uses a single smoothing step, allowing straightforward application to deflation-based preconditioners and two-level domain decomposition methods. Some new results are introduced whenever needed to complete the picture, and the material is self-contained thanks to a collection of new proofs, often shorter than the original ones.
We investigate the use of algebraic multigrid (AMG) methods for thesolution of large sparse linear systems arising from thediscretization of scalar elliptic partial differential equations withLagrangian finite elements of order at most 4. The resulting systemmatrices do not have the M-matrix property that is required bystandard analyses of classical AMG and aggregation-based AMGmethods. A unified approach is presented that allows us to extend theseanalyses. It uses an intermediate M-matrix and highlights the role ofthe spectral equivalence constant that relates this matrix to theoriginal system matrix. This constant is shown to be boundedindependently of the problem size and jumps in the coefficients of thepartial differential equations, provided that jumps are located atelements' boundaries. For two-dimensional problems, it is furthershown to be uniformly bounded if the angles in the triangulation alsosatisfy a uniform bound. This analysis validates the application ofthe AMG methods to the considered problems. On the other hand,because the intermediate M-matrix can be computed automatically, analternative strategy is to define the AMG preconditioners from thismatrix, instead of defining them from the original matrix. Numericalexperiments are presented that assess both strategies using publicly availablestate-of-the-art implementations of classical AMG and aggregation-based AMG methods.
Parallel performance study of block-preconditioned iterative methods on multicore computer systems
In this work we benchmark the performance of a preconditioned iterative method, used in large scale computer simulations of a geophysical application, namely, the elastic Glacial Isostatic Adjustment model. The model is discretized using the finite element method that gives raise to algebraic systems of equations with matrices that are large, sparse, nonsymmetric, indefinite and with a saddle point structure. The efficiency of solving systems of the latter type is crucial as it is to be embedded in a time-evolution procedure, where systems with matrices of similar type have to be solved repeatedly many times. The implementation is based on available open source software packages - Deal.II, Trilinos, PARALUTION and AGMG. These packages provide toolboxes with state-of-the-art implementations of iterative solution methods and preconditioners for multicore computer platforms and GPU. We present performance results in terms of numerical and the computational efficiency, number of iterations and execution time, and compare the timing results against a sparse direct solver from a commercial finite element package, that is often used by applied scientists in their simulations.