Drummond, Tony

Duff, Iain S.

Guivarch, Ronan

Ruiz, Daniel

Zenadi, Mohamed

Scientific Publications of the University of Toulouse II Le Mirail

Partitioning Strategies for the Block Cimmino algorithmMohamed Zenadi, Tony Drummond, Iain S. Duff, Ronan Guivarch, DanielRuizTo cite this version:Mohamed Zenadi, Tony Drummond, Iain S. Duff, Ronan Guivarch, Daniel Ruiz. Partition-ing Strategies for the Block Cimmino algorithm. International Conference On PreconditioningTechniques For Scientific And Industrial Applications, Preconditioning 2011, May 2011, Bor-deaux, France. <inria-00581285>HAL Id: inria-00581285https://hal.inria.fr/inria-00581285Submitted on 30 Mar 2011HAL is a multi-disciplinary open accessarchive for the deposit and dissemination of sci-entific research documents, whether they are pub-lished or not. The documents may come fromteaching and research institutions in France orabroad, or from public or private research centers.L’archive ouverte pluridisciplinaire HAL, estdestine´e au de´poˆt et a` la diffusion de documentsscientifiques de niveau recherche, publie´s ou non,e´manant des e´tablissements d’enseignement et derecherche franc¸ais ou e´trangers, des laboratoirespublics ou prive´s.Partitioning Strategies for the Block Cimmino algorithmList of authors:M. Zenadi ∗T. Drummond † I. S. Duff ‡ R. Guivarch ∗ D. Ruiz ∗We consider the Block Cimmino projection method for the solution of the linear system ofequations Ax = b, where A is a nonsingular large sparse unsymmetric matrix of order n. Theblocks are obtained by partitioning the system into p strips of rows, with p ≤ n, as in:A1...Apx =b1...bp . (1)The Block Cimmino method projects the current iterate simultaneously onto the manifoldscorresponding to the strips, and takes a convex combination of all the resulting vectors. Ouraim is to focus on preprocessing and partitioning strategies of the original system to bring it tothe form (1).If the matrixA is ill conditioned then there are some linear combinations of rows that are almostequal to zero and, as mentioned by [2], these linear combinations may occur inside the blocks oracross the blocks after row partitionings of the form (1). Assuming that the projections in theBlock Cimmino algorithm are computed exactly on the subspaces, then the rate of convergenceof the Block Cimmino algorithm depends only on the conditioning across the blocks. If weconsider additionally Conjugate Gradient acceleration of the Block Cimmino method, as in[1, 2], the convergence behaviour of the resulting method is directly linked to the spectrum ofthe n×n matrix ERJ formed as the sum of the previously mentioned projections, and given byERJ =p∑i=1AiT(AiAiT)−1Ai , (2)assuming, for simplicity, that the block-rows Ai have full row rank. An efficient implementationof Block Cimmino requires a combination of a robust method for computing the projectionsand a partitioning strategy that minimizes the ill-conditioning across the blocks.If the block rows Ai are nearly mutually orthogonal, i.e. AAT is strongly block-diagonallydominant, we can expect that the method will converge very quickly, if the projections arecomputed accurately. Conversely, the structure of AAT tells us about the orthogonality of thesubspaces represented by block partitions of A. In all our experiments we first normalize thematrix AAT using MC77 [5] so that the 2-norms of the rows and columns are close to one. Thenthe entries in AAT correspond to the cosine of the principal angle between every pair of rows,and we may expect that if such a cosine is relatively small the corresponding pair of rows are∗Universite´ de Toulouse, INPT(ENSEEIHT)-IRIT, France†Lawrence Berkeley National Laboratory‡CERFACS and Rutherford Appleton Laboratory1almost orthogonal. That is, if the block (i, j)th entry of the scaled AAT is close to zero then thesubspaces corresponding to the blocks Ai and Aj can be considered to be orthogonal. If AATis block tridiagonal, the blocks of A are such that the even numbered blocks are orthogonalas are also the odd-numbered blocks. Thus if we solve the projected subproblems accurately(using say a direct method) then we also solve the subproblems corresponding to the odd andeven numbered blocks accurately. The partition of A is thus of the form [A] =[B1TB2T]T,where B1 ={⋃iAi / i odd}and B2 ={⋃iAi / i even}. Such a partitioning is called a two-block partitioning. We denote by m1 and m2 the number of rows in B1 and B2 respectively. Itwas shown by [4] that, with such a two-block partitioning, the spectrum of matrix ERJ isλk = 1 + cosΨkλk = 1− cosΨk−m2λk = 1k = 1, . . . ,m2k = m2 + 1, . . . , 2m2k = 2m2 + 1, . . . , nwhere {Ψk}m21are the principal angles between R(B1T ) and R(B2T ).The main idea behind the so called two-block partitioning strategy is to exploit structuralorthogonality between the subspaces R(AiT ) defined by the partitioning (1). This structuralorthogonality can be analysed on the basis of the sparsity pattern of the normal equationsmatrix AAT . The idea is to find a permutation to transform the normal equations matrixAAT into block tridiagonal form. That is we determine a permutation matrix P such thatB = PAATPT is in block tridiagonal form using an implementation of algorithm of [3]. Wethen use this permutation to solve the row-wise permuted system of equationsPAx = Pb (3)using the block Cimmino algorithm.The main problem with this first preprocessing strategy for general sparse matrices, is that theblock tridiagonal structure has diagonal blocks with very differing sizes leading to a partitioningwith a bad degree of parallelism because of unbalanced tasks in the block-row projections.In our second preprocessing strategy, we take into account the numerical values in the matrixAAT , and do not enforce strict orthogonality between every second block of rows. We drop smallentries of AAT and by using the relationship to principal angles can relax the orthogonalitycondition in a controlled way in an attempt to obtain much better parallelism compared to theprevious strategy. We thus only keep nonzero entries above a given tolerance and permute theresulting filtered matrix into block tridiagonal form using the Cuthill-McKee algorithm and usethe row partitioning defined by this block tridiagonal structure to solve (3).Of course, since numerical values have been dropped from the normal equations matrix, theresulting partition will not provide two subsets of structurally orthogonal blocks of rows, but,if the values dropped are sufficiently small, we may expect that the numerical properties of theresulting iteration matrix will be relatively close to that of the “strict” two-block partitioningcase. Additionally, since the filtered matrix has less entries than the original normal equations2matrix, the resulting block tridiagonal permuted matrix will normally have a smaller bandwidththan B and this may help to define a partitioning on A with more blocks, better balancedprojections, and a higher degree of parallelism.To illustrate this preprocessing strategy, we show preliminary results obtained with the bmw3 2matrix (of size N = 227362) from the Tim Davis collection. The original matrix, before prepro-cessing, has the symmetric pattern shown in the left subplot above, in which there is no evidentstructural orthogonality between sets of rows.We perform some runs on a 24 (Intel i7) cores parallel machine, with a generic partitioninginto 24 equal sized blocks. The Block Cimmino method needs then 915 iterations to reach ascaled residual value of the order 10−7. It requires 5.46 sec for the analysis, 2.41 sec for thefactorization of the augmented systems and 63.25 sec for the iterations.Following the second preprocessing strategy, we filter AAT matrix under the value 0.1, whichcorresponds to the cosine of an angle of 84◦, and we obtain the pattern in the right subplot afterCuthill-McKee and a reduction of nonzero entries in AAT of 88%. With this block tridiagonalstructure, we can define 18 almost balanced partitions, the biggest one of size 16129 rows andthe smallest of 11062 rows. The Block Cimmino method needs then 296 iterations to convergeto the same level, with 7.19 sec for the analysis, 3.29 sec for the factorization, 24.62 sec for theiterations, and with 2.48 sec for the extra work in the preprocessing.We intend to develop and compare the benefits of such strategies on parallel architectures,different test examples and also with variants where we may drop entries in A in combinationwith entries in AATReferences[1] M. Arioli, I. S. Duff, J. Noailles, and D. Ruiz, (1992), A Block Projection Method For Sparse Matrices,SIAM J. Sci. Stat. Comput., 47–70.[2] R. Bramley and A. Sameh, (1990), Row projection methods for large nonsymmetric linear systems, Tech.Rep. 957, Center for Supercomputing Research and Development, University of Illinois, Urbana, IL. Also,SISSC, Vol. 13, January 1992.[3] E. Cuthill and J. McKee, (1969), Reducing the bandwidth of sparse symmetric matrices, in Proceedings24th National Conference of the Association for Computing Machinery, New Jersey, Brandon Press, 157–172.[4] T. Elfving, (1980), Block-iterative methods for consistent and inconsistent linear equations, Numer. Math.,35, 1–12.[5] D. Ruiz, (2001), A scaling algorithm to equilibrate both row and column norms in matrices, Tech. Rep.RAL-TR-2001-034, Rutherford Appleton Laboratory.3

Partitioning Strategies for the Block Cimmino algorithm (Precond 2011)

Abstract

Similar works

Full text

Available Versions

Scientific Publications of the University of Toulouse II Le Mirail