Higham - Accuracy and Stability of Numerical Algorithms (523152), страница 32
Текст из файла (страница 32)
Suppose thathas wi nonzeros in its ith row, i = 1:n. Showthat the inequality (7.27) can be replaced bywhereThis bound is potentially much smaller than (7.27)for large, sparse matrices.7.13. (D. J. Higham, after Fletcher [376, 198 5]) Suppose the nonsingular,square matrix A is perturbed to A + ∆A and b to b + ∆b. Then, to first order,the solution of Ax = b is perturbed to x + ∆x, whereP ERTURBATION T HEORY150FORLINEAR S YSTEMSSuppose that the perturbations have the formwhere theand δi are independent random variables, each having zeromean and variance σ2.
(As usual, the matrix E and vector f represent fixedtolerances.) Let ε denote the expected value.(a) Show thatwhere square brackets denote the operation of elementwise squaring: [B] ij=(b) Hence explain whymay be regarded as an “expected condition number” for the linear systemAx = b.(c) For the case where eij||A|| 2 and fj||b|| 2.
compare condexp(A, x)with the “worst-case” condition number κA , b (A, x) for the 2-norm.7.14. (Horn and Johnson [581,1991,p. 331]) Prove that for any nonsingularwhere is the Hadamard product (A B = (aijbij)) andis defined as inB||2 < ||A||2 ||B||2.) Discuss theProblem 7.9. (Hint: use the inequality ||Aattainability of this bound.Previous Home NextChapter 8Triangular SystemsIn the end there is left the coefficient of one unknown and the constant term.An elimination between this equation andone from the previous set that contains two unknownsyields an equation with the coefficient ofanother unknown and another constant term, etc.The quotient of the constant term by the unknownyields the value of the unknown in each case.-JOHN V.
ATANASOFF, Computing Machine for the Solution ofLarge Systems of Linear Algebraic Equations (1940)The solutions of triangular systems are usually computed to high accuracy.This fact . . . cannot be proved in genera/, for counter examples exist.However, it is true of many special kinds of triangular matrices andthe phenomenon has been observed in many others.The practical consequences of this fact cannot be over-emphasized.-G. W. STEWART, Introduction to Matrix Computations (1973)In practice one almost invariably finds thatif L is ill-conditioned, so thatthen the computed solution of Lx = b (or the computed inverse)is far more accurate than [standard norm bounds] would suggest.-J. H. WILKINSON, Rounding Errors in Algebraic Processes (1963)151152TRIANGULAR SYSTEMTriangular systems play a fundamental role in matrix computations. Manymethods are built on the idea of reducing a problem to the solution of oneor more triangular systems, including virtually all direct methods for solvinglinear systems.
On serial computers triangular systems are universally solvedby the standard back and forward substitution algorithms. For parallel computation there are several alternative methods, one of which we analyse in§8.4.Backward error analysis for the substitution algorithms is straightforwardand the conclusion is well known: the algorithms are extremely stable. Thebehaviour of the forward error, however, is intriguing, because the forwarderror is often surprisingly small --much smaller than we would predict fromthe normwise condition number κ, or, sometirnes, even the componentwisecondition number cond. The quotes from Stewart and Wilkinson at the startof this chapter emphasize the high accuracy that is frequently observed inpractice. The analysis we give in this chapter provides a partial explanationfor the observed accuracy of the substitution algorithms.
In particular, itreveals three important but nonobvious properties:• the accuracy of the computed solution from substitution depends stronglyon the right-hand side:• a triangular matrix may be much more or less ill conditioned than itstranspose; and• the use of pivoting in LU, QR, and Cholesky factorizations can greatlyimprove the conditioning of a resulting triangular system.As well as deriving backward and forward error bounds, we show how tocompute upper and lower bounds for the inverse of a triangular matrix.8.1.
Backward Error AnalysisRecall that for an upper triangular matrixbe solved using the formula xi =components of x in order from last to first.the system Ux = b canwhich yields theAlgorithm 8.1 (back substitution). Given a nonsingular upper triangularmatrixthis algorithm solves the system Ux = b.xn = bn /u n nfor i = n - 1:-1:1s = bifor j = i + 1:ns = s - uijxij8.1 B ACKWARD E RROR A NALYSIS153endx i = s/u iiendWe will not state the analogous algorithm for solving a lower triangular system, forward substitution. All the results below for back substitutionhave obvious analogues for forward substitution. Throughout this chapter Tdenotes a matrix that can be upper or lower triangular.To analyse the errors in substitution we need the following lemma.Lemma 8.2.
Let y =metic according tobe evaluated in floating point arith-s =cfor i = 1:1k-1s = s - ai biendy = s/bkThen the computedsatisfies(8.1)where |θi | < γi = iu/(1 - iu).Proof. Analysis very similar to that leading to (3.2) shows thatsatisfieswhere:=The final division yields, using (2.5),so that, after dividing through bywe haveThe result is obtained on invoking Lemma 3.1.Two remarks are in order. First, we chose the particular form of (8.1),in which c is not perturbed, in order to obtain a backward error result forUx = b in which b is not perturbed. Second, we carefully kept track of theterms 1 + δi in the proof, so as to obtain the best possible constants. Directapplication of the lemma to Algorithm 8.1 yields a backward error result.TRIANGULAR SYSTEMS154Theorem 8.3. The computed solutionfrom Algorithm 8.1 satisfiesTheorem 8.3 holds only for the particular ordering of arithmetic operationsused in Algorithm 8.1.
A result that holds for any ordering is a consequenceof the next lemma.Lemma 8.4. If y =is evaluated in floating point arithmetic, then. no matter what the order of evaluation,wherefor all i. If bk = 1, so that there is no division, thenfor all i.Proof. The result is not hard to see after a little thought , but a formalproof is tedious to write down. Note that the ordering used in Lemma 8.2 isthe one for which this lemma is least obvious! The last part of the lemmais useful when analysing unit lower triangular systems, and in various othercontexts.Theorem 8.5.
Let the triangular system TX = b, whereis nonsingular, be solved by substitution, with any ordering. Then the computedsolutionsatisfiesIn technical terms, this result says thathas a tiny componentwise relativebackward error. In other words, the backward error is about as small as wecould possibly hope.In most of the remaining error analyses in this book, we will derive results that, like the one in Theorem 8.5, do not depend on the ordering ofthe arithmetic operations.
Results of this type are more general, usually noless informative. and easier to derive, than ones that depend on the ordering. However, it is important to realise that the actual error does depend onthe ordering, possibly strongly so for certain data. This point is clear fromChapter 4 on summation.1558.2 F ORWARD E RROR A NALYSIS8.2. Forward Error AnalysisFrom Theorems 8.5 and 7.4 there follows the forward error boundwhereThis bound can, of course, be arbitrarily smaller than the correspondingbound involvingfor the reasons explained in Chapter 7. For further insight, note that, in terms of the traditional conditionnumber, κ(T), ill conditioning of a triangular matrix stems from two possible sources: variation in the size of the diagonal elements and rows withoff-diagonal elements which are large relative to the diagonal elements. Significantly, because of its row scaling invariance, cond(T , x) is susceptible onlyto the second source.Despite its pleasing properties, cond(T , x) can be arbitrarily large.
Thisis illustrated by the upper triangular matrix(8.2)for which(8.3)We have cond(U(a),e) = cond(U( a))2an-1 as aTherefore WCcannot assert that all triangular systems are solved to high accuracy. Nevertheless, for any T there is always at least one system for which high accuracyis obtained: the system TX = e1 if T is upper triangular, or TX = en if Tis lower triangular.
In both cases cond(T, x) = 1, and the solution comprisesthe computation of just a single scalar reciprocal.To gain further insight we consider special classes of triangular matrices,beginning with one produced by certain standard factorizations with pivoting.In all the results below, the triangular matrices are assumed to be n × n andnonsingular, andis the computed solution from substitution.Lemma 8.6.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
