Higham - Accuracy and Stability of Numerical Algorithms (523152), страница 37
Текст из файла (страница 37)
Explainingthis fact remains one of the major unsolved problems in numerical analysis.The best attempt to date is by Trefethen and Schreiber [1019, 1990]. whodevelop a statistical model of the average growth fact or for partial pivotingand complete pivoting. Their model supports their empirical findings that forvarious distributions of random mat rices the average growth fact or (normalized by the standard deviation of the initial matrix elements) is close to n 2 / 3for partial pivoting and n ½ for complete pivoting. Extensive experiments byEdelman suggest that for random matrices from the normal N(0, 1) distribution the unnormalized growth factor for partial pivoting grows like n ½ [345,1 9 9 5 ].We turn now to complete pivoting. Wilkinson [1085, 1961, pp.
282-285]showed that(9.13)and that this bound is not attainable. The bound is a much more slowlygrowing function than 2n -1, but can still be quite large (e.g., it is 3570 forn = 100). As for partial pivoting, in practice the growth factor is usuallysmall. Wilkinson stated that “no matrix has been encountered in practice forwhich p 1 /pn was as large as 8” [1085, 1961, p.
285] and that “no matrix hasyet been discovered for which f(r) > r” [1089, 1965, p. 213] (here. pi is the(n - i + 1)st pivot and f(r)Cryer [256, 1968] defined(9.14)The following results are known:• g(2) = 2 (trivial).• g(3) = 2¼; Tornheim [1012,• g (4) = 4; Cryer [256,1968]• g(5) < 5.005; Cohen [229,1965]and Cohen [229,and Cohen [229,1974].1974].1974].9.4 S PECIAL M ATRICES181Tornheim [1012, 1965] (see also Cryer [256, 1968]) showed that> nfor any n × n Hadamard matrix Hn (a bound which, as we saw above, holdsfor any form of pivoting).
For n such that a Hadamard matrix does not exist,the best known lower bound is g(n) >= (n + 1)/2 (see (9.11)).Cryer [256, 1968] conjectured that for real matrices< n, with equality if and only if A is a Hadamard matrix. The conjecture< n becameone of the most famous open problems in numerical analysis, and has beeninvestigated by many mathematicians. The conjecture was finally shown to befalse in 1991. Using a package LANCELOT [236, 1992] designed for large-scalenonlinear optimization, Gould [474, 1991] discovered a 13 × 13 matrix for whichthe growth factor is 13.0205 in IEEE double precision floating point arithmetic. Edelman subsequently showed, using the symbolic manipulation packages Mathematica and Maple, that a growth factor 13.02 can be achieved inexact arithmetic by making a small perturbation (of relative size 10-7) to oneelement of Gould’s matrix [338, 1992], [348, 1991].
A more striking counterex= 32.986341 [338,ample to the conjecture is a matrix of order 25 for which1992]. Interesting problems remain, such as determining limng(n)/n andevaluatingfor Hadamard matrices (see Problem 9.15).For complex matrices the maximum growth factor is at least n for anyn, since> n (see (9.12)). The growth can exceed n, even for n = 3:Tornheim [1012, 1965] constructed the examplefor which= 3.079.9.4.
Special MatricesFor matrices with certain special properties, more can be said about the behaviour of GE and, in particular, the size of the growth factor.As a first example, supposeis diagonally dominant by rows,or diagonally dominant by columns, that is, A* is diagonally dominant byrows. Then GE without pivoting is perfectly stable.Theorem 9.8 (Wilkinson).
Ifis diagonally dominant by rows orcolumns then A has an LU factorization without pivoting and the growth factorpn < 2. If A is diagonally dominant by columns then |lij| < 1 for all i andj in the LU factorization without pivoting (hence GEPP does not require anyrow interchanges).182LU FACTORIZATIONANDLINEAR E QUATATIONSFigure 9.1. A banded matrix.Proof.
The result follows immediately from the more general Theorems 12.5 and 12.6 for block diagonally dominant matrices.Note that for a matrix diagonally dominant by rows the multipliers canbe arbitrarily large but, nevertheless, pn < 2, so GE is perfectly stable.A smaller bound for the growth factor also holds for an upper Hessenbergmatrix.
(H is upper Hessenberg if hij = 0 for i > j + 1.)Theorem 9.9 (Wilkinson). Ifis upper Hessenberg then< n.Proof. The structure of an upper Hessenberg H means that at each stageof GEPP we just add a multiple of the pivot row to the next row (afterpossibly swapping these two rows). That< n is a consequence of thefollowing claim, which is easily proved by induction: at the start of stage k ,row k + 1 of the reduced matrix is the same as row k + 1 of the original matrix,and the pivot row has elements of modulus at most k times the largest elemento f H.has lower bandwidth p and upper bandwidth q if aij = 0A matrixfor i > j + p and j > i + q; see Figure 9.1.
It is well known that in an LUfactorization of a banded matrix the factors inherit A's band structure: Lhas lower bandwidth p and U has upper bandwidth q. If partial pivoting isused then, in PA = LU. L has at most p + 1 nonzeros per column and U hasupper bandwidth p + q. (For proofs of these properties see Golub and VanLoan [470, 1989, §4.3].) It is not hard to see that for a banded matrix, γ n inTheorem 9.3 can be replaced by γ m a x (p + 1 , q+1) and 2γ n in Theorem 9.4 can bereplaced by γm a x (p + 1 , q+1) + γp + q + 1 .’The following result bounds the growth factor for partial pivoting on abanded matrix.9.5 T RIDIAGONAL M ATRICES183Theorem 9.10 (Bohte).
Ifhas upper and lower bandwidths p thenand this bound is almost attainable when n = 2p+1.Proof. See Bohte [131, 1975]. An example with n = 9 and p = 4 in whichequality is almost attained is the matrixwhereis an arbitrarily small positive number, which ensures that rows 1and 5 are interchanged on the first stage of the elimination, this being theonly row interchange required. Ignoring terms inthe last column of Uin PA = LU is [1, 1, 2, 4, 8, 16, 31, 60, 116]T and the growth factor is116.A special case of Theorem 9.10 is the easily verified result that for a tridiagonal matrix,Hence GEPP achieves a small normwise backwarderror for tridiagonal matrices. In the next section we show that for severaltypes of tridiagonal matrix GE without pivoting achieves a small componentwise backward error.9.5.
Tridiagonal MatricesConsider the nonsingular tridiagonal matrixand assume that A has an LU factorization A = LU, where(9.15)184LU F ACTORIZATIONANDLINEAR E QUATIONSGE for computing L and U is described by the recurrence relations(9.16)For the computed quantities, we haveHenceIn matrix terms these bounds may be written as(9.17)If the LU factorization is used to solve a system Ax = b by forward and backsubstitution then it is straightforward to show that the computed solutionsatisfies(9.18)Combining (9.17) and (9.18) we have, overall,(9.19)The backward error result (9.19) applies to arbitrary nonsingular tridiagonal A having an LU factorization. We are interested in determining classesof tridiagonal A for which a bound of the form |∆ A| < g (u) |A| holds.
Such abound will hold ifas noted in §9.2 (see (9.8)).Three classes of matrices for whichholds for the exact Land U are identified in the following theorem.Theorem 9.11. Letbe nonsingular and tridiagonal. If any of thefollowing conditions hold then A has an LU factorization and |L||U| = |LU|:(a) A is symmetric positive definite;(b) A is totally nonnegative, or equivalently, L > 0 and U > 0;(c) A is an M-matrix, or equivalently, L and U have positive diagonalelements and nonpositive off-diagonal elements;(d) A is sign equivalent to a matrix B of type (a)-(c), that is, A = D 1 BD2 ,where |D1 | = |D2 | = I.9.5 T RIDIAGONAL M ATRICES185Proof. For (a), it is well known that a symmetric positive definite A hasan LU factorization in which U = DL T , where D is diagonal with positivediagonal elements. Hence |L||U| = |L||D||LT| = |LDLT| = |LU|, where themiddle equality requires a little thought.
In (b) and (c) the equivalences,and the existence of an LU factorization, follow from known results on totallynonnegative matrices [258, 1976] and M-matrices [94, 1994]; |L||U| = |LU| isimmediate from the sign properties of L and U. (d) is trivial.For diagonally dominant matrices, |L||U| is not equal to |LU| = |A|, butit cannot be much bigger.is nonsingular, tridiagonal, and diagTheorem 9.12.
Supposeonally dominant by rows or columns, and let A have the LU factorizationA = LU. Then |L||U| < 3|A|.Proof. If |i - j| = 1 then (|L||U|) ij = |aij|, so it suffices to consider thediagonal elements and show that (using the notation of (9.15))|li ei-1 | + |ui | < 3|di |.The rest of the proof is for the case where A is diagonally dominant by rows;the proof for diagonal dominance by columns is similar.First, we claim that |ei | < |ui | or all i. The proof is by induction.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
