Modified Sparse Approximate Inverses (MSPAI) for Parallel Preconditioning
A. Kallischko · mediaTUM – the media and publications repository of the Technical University Munich (Technical University Munich) · 2008
samkeit von MSPAI wird an einigen numerischen Beispielen wie Gebietszerlegung und Stokes-Problemen nachgewiesen.Neben der theoretischen Entwicklung des MSPAI thematisiert diese Arbeit auch eingehend dessen effiziente Implementierung.Zunächst werden dünnbesetzte Verfahren zur Lösung von Least-Squares-Problemen in Verbindung mit QR-Updates untersucht und mit den bisherigen Standardmethoden verglichen.Außerdem wird eine Caching-Strategie eingeführt, die dabei hilft, redundante QR-Zerlegungen im Fall hochstrukturierter Matrizen einzusparen.Schließlich rundet die Unterstützung maximaler Oberpattern für die Dünnbesetztheitsstrukturen die Implementierung ab.Laufzeittests beweisen die Effektivität dieser Maßnahmen in Form von deutlich reduzierten Laufzeiten im Vergleich mit der aktuellen Standardimplementierung von SPAI.Without computers or calculators, the solution of large systems still was impossible.In the 1880s, it took about 7 weeks to solve a linear system Ax = b of dimension n = 29.Hereby, people had to use logarithm tables.In 1951, it still took over 3 hours to solve a system of dimension 17 [56].On an Intel Centrino with 1600 MHz, MATLAB nowadays needs 0.6 seconds to compute the solution for n = 1000.Despite this enormous raise in computation power, Gaussian elimination is only suitable for systems up to n ≈ 10000, since it is an O n 3 algorithm.In order to solve much larger systems, iterative methods were developed.They can also take into account sparsity or certain band or block structures, which makes them usually highly efficient and the methods of choice today.These iterative approaches enable us to solve the large linear systems which arise in many numerical applications such as the discretization of partial differential equations or image restoration problems.Even these highly elaborate iterative techniques can be prohibitively slow if the coefficient matrix A of the system is ill-conditioned.Then, either convergence can not be achieved at all or a huge number of steps becomes necessary.In order to accelerate the solution process, preconditioning techniques have been investigated since the early 1970s.Hereby, the linear system is transformed into an equivalent one with much lower condition number leading to considerably lower iteration counts and thus reduced solution time.In Scientific Computing, it is common knowledge that the choice of a good preconditioner is even more important than the actual choice of the iterative solver.Many preconditioning techniques were developed so far [25].Among the most robust and flexible ones is the class based on Frobenius norm minimization.Extending this approach, Grote and Huckle [48] developed the sparse approximate inverse (SPAI) preconditioner which is able to update a given start sparsity pattern automatically.Furthermore, SPAI is inherently parallel, since its computation splits up into independent least squares problems, one for each column.There is also the factorized variant FSPAI [62] for symmetric positive definite matrices.Another important group of preconditioners are the modified incomplete factorizations [4,49].They yield preconditioners which preserve the row sums, i.e. the action on the vector of all ones.In many applications, this leads to a significantly improved convergence behavior.Another type of preconditioners, which are constructed such that they are optimal with respect to some given subspace, is classical probing [24].However, Finally, we state a theorem considering SPAI and FSPAI for M-matrices, a class of matrices frequently arising from the discretization of partial differential equations.Chapter 3 starts off with the description of the classical probing technique and the group of modified preconditioners.Both have in common that they do not only yield explicit approximations to the system matrix, but also include additional information in the preconditioner, hereafter referred to as probing information.Typically, probing and modified-type preconditioners are restricted to the conservation of row sums in the preconditioner.Moreover, classical probing can only compute very sparse approximations with simple structures, since people have to invert them in every iteration of the iterative solver.Additionally, these methods are hard to parallelize.To overcome these restrictions and drawbacks, we generalize, in the third section, the class of Frobenius norm minimization based preconditioners to target form [57]. Furthermore, we extend the idea from n × n to rectangular m × n matrices.Thus, we gain inherent algorithmic parallelism due to the Frobenius norm minimization.Moreover, we can add several rows to the input matrices which contain the probing information in form of probing vectors.Here is the next gain in generality: we can add as many probing vectors of any type we consider as meaningful for the actual application.We also benefit from the other advantages, which we have already seen for SPAI.The computation is not restricted to a special sparsity pattern, and we can perform automatic pattern updates.Since it is based on modified preconditioners, we call this formulation modified sparse approximate inverse (MSPAI ).With this, we can compute inverse approximations like SPAI, adding probing information.However, we can also get direct approximations of a matrix, for instance a sparser one which acts on certain subspaces as the original matrix does.For both inverse or direct approximations, we also derive formulations of MSPAI for factorized approximations.This allows improving any given factorized preconditioner by adding probing information.Another application is the MSPAI probing of Schur complements, which is the application domain classical probing was designed for.Which concrete probing vectors can or should be chosen for MSPAI is discussed in Section 3.4.For this purpose, we provide a theoretical basis, i.e. a purely structural condition considering probing vectors' sparsity structure and the sparsity structure of the preconditioner.The final section in this chapter derives several symmetrization techniques, since matrices resulting from Frobenius norm minimization processes are usually unsymmetric.However, for symmetric and even symmetric positive definite input matrices, we want to have a preconditioner reflecting these properties.We present one method based on Frobenius norm minimization and another approach involving two consecutive iteration steps.A symmetrization technique for factorized preconditioners concludes the chapter.Chapter 4 gives insight to a topic which is crucial for every numerical algorithm: its efficient implementation.MSPAI can only be employed and tested efficiently if we provide a fast code.We review the use of sparse methods for the solution of least squares problems and combine them with the idea of QR updates [17,46].The latter arise in SPAI-type computations involving pattern updates.An experiment will also evaluate the potential of single precision QR methods for our purpose.Section 4.2 introduces the idea of a caching algorithm for MSPAI computation.It is based on the observation that in highly structured matrices, we perform many redundant factorizations.With a caching strategy holding already computed factorizations, we try to exploit these redundancies to reduce the overall computation time.After that, we investigate the benefit from a maximum sparsity pattern [61] in parallel computing environments -an idea which was already discussed theoretically in Chapter 2. A profound comparison of the current standard SPAI implementation SPAI 3.2 and our MSPAI 1.0 code rounds up the chapter.For each single improvement, we provide detailed