Linear Systems: Direct Methods

In this chapter we consider numerical methods for solving a system of linear equations Ax = b. We assume that the given matrix A is real, n × n, and nonsingular and that b is a given real vector in Rn , and we seek a solution x that is necessarily also a vector in Rn . Such problems arise frequently in virtually any branch of science, engineering, economics, or finance. There is really no single technique that is best for all cases. Nonetheless, the many available numerical methods can generally be divided into two classes: direct methods and iterative methods. The present chapter is devoted to methods of the first type. In the absence of roundoff error, such methods would yield the exact solution within a finite number of steps. The basic direct method for solving linear systems of equations is Gaussian elimination, and its various aspects and variants occupy the first seven sections of this chapter. Section 5.1 presents the method in simple terms. The bulk of the algorithm involves only the matrix A and amounts to its decomposition into a product of two matrices that have a simpler form. This is called an LU decomposition, developed in Section 5.2. Such an alternative view is useful, for instance, when there are several right-hand-side vectors b each requiring a solution, as the LU decomposition can then be shared. The simple algorithm of Section 5.1 is not guaranteed to be stable or even well-defined in general. In Section 5.3 we modify it using pivoting strategies to make the algorithm of Gaussian elimination with pivoting practically stable in general. Then in Section 5.4 we discuss efficient methods for implementing these algorithms. The following three sections are devoted to variants of the same basic algorithm for special classes of matrices. Important simplifications can be made for symmetric positive definite matrices, as discussed in Section 5.5. If the matrix A contains only a few nonzero elements, then it is called sparse. Techniques for storing such matrices and for solving linear systems for the special case where all elements are zero outside a narrow band along the main diagonal are presented in Section 5.6. Direct methods for handling more general sparse matrices are briefly considered in Section 5.7. In Section 5.8 we consider not the numerical solution of Ax = b but rather its assessment. Since there are roundoff errors, the result of any of our algorithms is approximate, and so the question of how close the approximate solution is to the exact one arises. This in turn depends on the condition number of the matrix, a concept that is defined and further discussed here.