Large tridiagonal and block tridiagonal linear systems on vector and parallel computers (Q1095581)

For solving large tridiagonal and block tridiagonal linear systems, many of the otherwise powerful, iterative algorithms are not suitable enough for execution on modern vector and parallel computers. This is because of the recurrence relations contained in the algorithms. A number of techniques have been developed which help to make the solution of tridiagonal linear systems more suitable for ``modern architectures''. These basic techniques, an overview of which is given in Section 2, may easily be generalized for use in connection with block tridiagonal linear systems. These generalized techniques are discussed in Section 3 within the bounds of a model problem of the form \(Ax=y\). For many systems, the convergence behavior is considerably improved when an incomplete decomposition \(LD^{-1}U\) of A, \(A=LD^{-1}U-R\), is used as a preconditioning matrix for, e.g., the conjugate gradient method. Here D is a diagonal matrix, L a lower and U an upper triangular matrix. This ``standard incomplete decomposition'' (briefly, s.i.d.) exists when A is an M-matrix; it is unique if some given conditions are satisfied. Some of the earlier papers of the author are associated with these questions. In order to obtain still better suitability the author proposes in Section 4, as the contribution of the present paper, a new technique that, instead of s.i.d., leads to an ``incomplete parallel decomposition'' (briefly, i.p.d.). This decomposition, \(A\approx PQ\), satisfies all the above-mentioned requirements except that the factors no longer are triangular matrices. In fact, i.p.d. is a generalization of the ``twisted factorization'', known in 1D-case and discussed in Section 2. Such a generalization exists if A is an M-matrix. In subsection 4.2 the algorithms needed for the construction of the i.p.d. as well as the back substitutions, required for the preconditioning, are given in detail within the bounds of the same model problem. It is shown that the i.p.d. of A can be computed almost entirely in two parallel parts and also the back substitution can be done largely in parallel. The effectiveness of the new method is illustrated by numerical examples in section 5.