Higham - Accuracy and Stability of Numerical Algorithms (523152), страница 30
Текст из файла (страница 30)
A matrix is irreducible if it cannot be symmetrically permuted toblock triangular form. A Perron vector of B > 0 is a nonnegative eigenvectorcorresponding to the eigenvalue p(B), where p denotes the spectral radius.Theorem 7.8 (Bauer). Letand |A - 1 ||A| are irreducible. Thenbe nonsingular and suppose that |A||A- |1(7.21)The minimum is attained for D1 = diag(x)-1 and D2 = diag(|A- 1 |x), wherex > 0 is a right Perron vector of |A||A - 1 | (so that |A||A- 1 | = p(|A||A - 1 |)x).Proof. See Problem 7.9.For the Kahan example (7.14),and, in fact,3 for D =two-sided scaling is nearly optimal in this case.so a symmetric140P ERTURBATION T HEORYFORLINEAR S YSTEMS7.4. The Matrix InverseWe briefly discuss component wise perturbation theory for the matrix inverse.With X := A-1 and X + ∆X := (A + ∆X)-1, a componentwise conditionnumber isIn general, the inequality- is strict.
but there is equality when |A -1| = D1 A - 1D 2for D i of the form diag(±1), [407, 1992, Thm. 1.10], [439, 198 2]. Anothercomponentwise condition number is evaluated in Problem 7.10. We saw inTheorem 6.5 that the reciprocal of the normwise condition number for matrixinversion is the normwise relative distance to singularity. Is the same truefor an appropriate componentwise condition number? The componentwisedistance to singularity,has been characterized by Rohn [877,1989],[878,1990]aswhere the maximum is taken over all signature mat rices Si = diag(±1) andwhereThis formula involves 4n eigenproblems and thus is computationally intractable(in fact it has been shown to be NP-hard by Poljak and Rohn [836, 1993]).Demmel [285, 1992] shows by complexity arguments that there can be nosimple relationship between dE (A) and the quantitywhich is anupper bound for µ E (A).
He also presents evidence for the conjecture thatfor a constant γn . The lower bound always holds and Demmel identifiesseveral classes of matrices for which the upper bound holds. This conjectureis both plausible and aesthetically pleasing because d|A| (A) is invariant undertwo-sided diagonal scalings of A and p(|A- 1 ||A|) is the minimum-normcondition number achievable by such scalings, as shown by Theorem 7.8.7.5. ExtensionsThe componentwise analyses can be extended in three main ways.7.6 N UMERICAL S TABILITY141(1) We can use more general measures of size for the data and the solution.Higham and Higham [528, 1992] measure ∆A , ∆b, and ∆x bywhereand the e ij , fi , and gi aretolerances. They show that the corresponding backward error is given by theexplicit formulawhere r = b - Ay, Dj = diag(e j1, .
. . , ejn, fj), and p -1 + q -1 = 1; boundsfor the corresponding condition number are also obtained. Theorem 7.3, andTheorem 7.4 with the-norm, correspond to p =and giIfwe take p =and g = |x|, we are measuring the change in the solution ina componentwise relative sense, asand thecondition number is [528, 1992]This latter case has also been considered by Rohn [876, 1989] and Gohbergand Koltracht [455, 1993].
It is also possible to obtain individual boundsforby refraining from taking norms in the analysis: seeChandrasekaran and Ipsen [197, 1995] and Problem 7.1.(2) The backward error results and the perturbation theory can be extended to systems with multiple right-hand sides. For the general vp measuredescribed in (1), the backward error can be computed by finding the minimum p-norm solutions to n underdetermined linear systems. For details, seeHigham and Higham [528, 1992].(3) Structure in A and b can be preserved in the analysis.
For example, if Ais symmetric or Toeplitz then its perturbation can be forced to be symmetricor Toeplitz too, while still using componentwise measures. References includeHigham and Higham [527, 1992] and Gohberg and Koltracht [455, 1993] forlinear structure, and Bartels and D. J. Higham [76, 1992] for Vandermondestructure. A symmetry-preserving normwise backward error is explored byBunch, Demmel, and Van Loan [163, 1989], while Smoktunowicz [930, 1995]considers the componentwise case (see Problem 7.11).
Symmetry-preservingnormwise condition numbers are treated by D. J. Higham [526, 1995].7.6. Numerical StabilityThe backward errors examined in this chapter lead to definitions of numericalstability of algorithms for solving linear systems. Precise and formal definitions of stability can be given, but there are so many possibilities, across142PERTURBATION THEORYFORLI N E A R S Y S T E M Sdifferent problems, that to define and name each one tends to cloud the issuesof interest. We therefore adopt an informal approach.A numerical method for solving a square, nonsingular linear system Ax = bis normwise backward stable if it produces a computed solution such thatis of order the unit roundoff.
How large we allowto be.while still declaring the method backward stable, depends on the context,. Itis usually implicit in this definition that= O(u) for all A and b, anda method that yields= O(u) or a particular A and b is said to haveperformed in a normwise backward stable manner.The significance of normwise backward stability is that the computed solution solves a slightly perturbed problem, and if the data A and b containuncertainties bounded only normwisewithand similarly for b), thenmay be the exact solution to the problem wewanted to solve, for all we know.Componentwise backward stability is defined in a similar way: we now reto be of order u.
This is aquire the componentwise backward errormore stringent requirement than normwise backward stability. The roundingerrors incurred by a met hod that is componentwise backward stable are insize and effect equivalent to the errors incurred in simply converting the dataA and b to floating point numbers before the solution process begins.If a method is normwise backward stable then, by Theorem 7.2, the forward erroris bounded by a multiple of κ(A )u.
However, a met hodcan produce a solution whose forward error is bounded in this way without thenormwise backward errorbeing of order u. Hence it is useful to definea method for whichas normwise forward stable.By similar reasoning involvingwe say a method is componentwiseforward stable ifTable 7.1 summarizes thedefinitions and the relations between them. There are several examples inthis book of linear-equation-solving algorithms that are forward stable butnot backward stable: Cramer’s rule for n = 2 (§1.10.1).
Gauss-Jordan elimination (§13.4), and the seminormal equations method for underdeterminedsystems (§20.3).Other definitions of numerical stability can be useful (for example, rowwise backward stability means thatand they will beintroduced when needed.7.7. Practical Error BoundsSuppose we have a computed approximationto the solution of a linearsystem Ax = b, whereWhat error bounds should we compute?1437.7 P RACTICAL E RROR B OUNDSTable 7.1.
Backward and forward stability.Componentwise backward stabilityComponentwise forward stabilityNormwise backward stabilityNormwise forward stabilityThe backward error can be computed exactly, from the formulae(7.23)at the cost of one or two matrix-vector products, for r = b andThe only question is what to do if the denominator is so small as to causeoverflow or division by zero in the expression forThis could happen, for example, when E = |A| and f = |b| and, for some i , aijxj = 0 forall j, as is most likely in a sparse problem. LAPACK’s xyyRFS (“refine solution”) routines apply iterative refinement in fixed precision, in an attemptIf the ith component of the denominator in (7.23)to satisfyis less than safe_min/u, where safe_min is the smallest number such that1/safe_min does not overflow, then they add (n + 1) safe_min to the ith components of the numerator and denominator.
A more sophisticated strategy isadvocated for sparse problems by Arioli, Demmel, and Duff [24, 1989]. Theysuggest modifying the formula (7.23) by replacing |bi | in the denominator bywhen the ith denominator is very small. See [24, 1989] fordetails and justification of this strategy.Turning to the forward error, one approach is to evaluate the forward errorbound from Theorem 7.2 or Theorem 7.4, withequal to the correspondingbackward error. Because x in (7.9) is unknown, we should use the modifiedbound(7.24)If we have a particular E and f in mind for backward error reasons, then it isnatural to use them in (7.24). However, the size of the forward error boundvaries with E and f, so it is natural to ask which choice minimizes the bound.P ERTURBATION T HEORY144FORLINEAR S YSTEMSLemma 7.9.
The upper bound in (7.23) is at least as large as the upper boundin(7.25)and is equal to it whenProof. First note that r = b (7.25). Now, for z > 0,is a multiple of |r|.implieswith equality if z is a multiple of r. Taking z =with equality whenis preserved whenwhich impliesgivesis a multiple of |r|. The truth of this statement-norms are taken, so the result follows.Since the bound (7.25) is obtained by taking absolute values in the equationit is clearly the smallest possible such bound subjectto ignoring signs in A-1 and r. It is reasonable to ask why we do not takeas our error bound.