Aasen's Method for Linear Systems with Self-Adjoint Matrices

R. J. Hanson

Visual Numerics, Inc.

July, 1997

Abstract

Aasen's algorithm is effective for solving systems of linear algebraic equations with self-adjoint matrices. The method is described by Golub and Van Loan1 for real or symmetric matrices. This note gives the formulas for complex matrices. We also mention a technique for stabilizing the diagonal elements of the tri-diagonal factor. The original research appeared in a paper by Aasen2 and was evaluated by Barwell and George3. The evaluation concluded that Aasen's method and a method due to Bunch and Kaufman 4 were competitive and effective for this class of linear algebraic equations. Both techniques require the same number of complex operations, . A complex version of the Bunch and Kaufman method is used in LAPACK User's Guide5. Visual Numerics has used the Aasen method since 1993 in the IMSL Fortran 90 routine, lin_sol_self. This software solves self-adjoint problems and includes an option for positive-definite matrices.

Introduction

For complex data, the self-adjoint linear system is , in which . Extending Aasen's formulas, we develop the factorization , in which the matrices are, respectively, a permutation, a unit-lower triangular, and a self-adjoint tri-diagonal. This decomposition allows for the solution of a linear system by solving for , then the tri-diagonal system for , and finally for with the variables rearranged, . The permutation is chosen so that . We use standard elimination with partial pivoting for solving the tri-diagonal system, disregarding the symmetry.

Computing the Decomposition

Using the notation of Golub and Van Loan, we write .

The computational steps use stores , and.

  1. Initialize values . Remark: The values are implicitly noted and are not computed.
  2. Do through Step 14, for .
  3. Set
  4. Do through Step 5 for .
  5. Set
  6. Set . Remark: Each has a mathematically zero imaginary part, which must be enforced in the factorization. We observed that without this assignment, the imaginary parts magnify. Eventually this leads to instability and incorrect results for the tri-diagonal matrix, .
  7. Set .
  8. If , go to Step 14.
  9. Compute the pivot, .
  10. Apply the pivot: Exchange entries .
  11. Apply the pivot: Exchange rows and of matrix .
  12. Apply the pivot. Exchange columns and , in rows of matrix .
  13. Set . If , .
  14. End of the loop on .

References

  1. Golub, G. H. & Van Loan, C. F., Matrix Computations, 3rd Ed., John Hopkins University Press, Baltimore, Md., USA (1996).
  2. Aasen, J. O. "On the Reduction of a Symmetric Matrix to Tri-diagonal Form." BIT, 11, p. 233-242.
  3. Barwell, V. & George, A. "A Comparison of Algorithms for Solving Symmetric Indefinite Systems of Linear Equation." ACM-Trans. Math. Soft., 2, No. 3, p. 242-251.
  4. Bunch, J. R. & Kaufman, L. "Some Stable Methods for Calculating Inertia and Solving Symmetric Linear Systems." Math. Comp., 31, p. 162-179.
  5. Anderson, E., et al. LAPACK User's Guide, SIAM Publication, Philadelphia, Pa., USA (1995).